Changeset 8999 for branches/2017

Timestamp:

2017-12-13T09:49:08+01:00 (6 years ago)

Author:

timgraham

Message:

Tidying of interp and added interpavm for vertical refinement case

File:

• branches/2017/dev_METO_MERCATOR_2017_agrif/NEMOGCM/NEMO/NST_SRC/agrif_opa_interp.F90 (modified) (5 diffs)

branches/2017/dev_METO_MERCATOR_2017_agrif/NEMOGCM/NEMO/NST_SRC/agrif_opa_interp.F90

@@ -873 +873 @@
       zrhoy = Agrif_rhoy()
       IF (before) THEN 
-         !We can't use zero as the special value because we need to include zeros
-         !when interpolating the scale factors
-         IF(Agrif_UseSpecialValue) THEN 
-!             Agrif_SpecialValue = -999._wp
-             Agrif_SpecialValue = 0._wp
-         ELSE
-             Agrif_SpecialValue = 0._wp
-         ENDIF
          DO jk=1,jpk
             DO jj=j1,j2
                DO ji=i1,i2
-                  ptab(ji,jj,jk,1) = (e2u(ji,jj) * e3u_n(ji,jj,jk) * un(ji,jj,jk)*umask(ji,jj,jk)) - &
-                                   & ((umask(ji,jj,jk)-1) * Agrif_SpecialValue)
+                  ptab(ji,jj,jk,1) = (e2u(ji,jj) * e3u_n(ji,jj,jk) * un(ji,jj,jk)*umask(ji,jj,jk)) 
 # if defined key_vertical
                   ptab(ji,jj,jk,2) = (umask(ji,jj,jk) * e2u(ji,jj) * e3u_n(ji,jj,jk))
@@ -898 +889 @@
          western_side  = (nb == 1).AND.(ndir == 1)
          eastern_side  = (nb == 1).AND.(ndir == 2)
-
-         Agrif_SpecialValue = 0._wp ! reset specialvalue to zero now interpolation completed
 
          DO ji=i1,i2
@@ -976 +965 @@
       !      
       IF (before) THEN          
-         IF(Agrif_UseSpecialValue) THEN 
-!             Agrif_SpecialValue = -999._wp
-             Agrif_SpecialValue = 0._wp
-         ELSE
-             Agrif_SpecialValue = 0._wp
-         ENDIF
          DO jk=k1,k2
             DO jj=j1,j2
                DO ji=i1,i2
-                  ptab(ji,jj,jk,1) = (e1v(ji,jj) * e3v_n(ji,jj,jk) * vn(ji,jj,jk)*vmask(ji,jj,jk)) - &
-                                   & ((vmask(ji,jj,jk)-1) * Agrif_SpecialValue)
+                  ptab(ji,jj,jk,1) = (e1v(ji,jj) * e3v_n(ji,jj,jk) * vn(ji,jj,jk)*vmask(ji,jj,jk))
 # if defined key_vertical
                   ptab(ji,jj,jk,2) = vmask(ji,jj,jk) * e1v(ji,jj) * e3v_n(ji,jj,jk)
@@ -996 +978 @@
          zrhox = Agrif_rhox()
 # if defined key_vertical
-         Agrif_SpecialValue = 0._wp !Reset special value to zero now interpolation is done
 
          southern_side = (nb == 2).AND.(ndir == 1)
@@ -1403 +1384 @@
 # if defined key_zdftke || defined key_zdfgls
 
-   SUBROUTINE interpavm( ptab, i1, i2, j1, j2, k1, k2, before )
+   SUBROUTINE interpavm( ptab, i1, i2, j1, j2, k1, k2, m1, m2, before )
       !!----------------------------------------------------------------------
       !!                  ***  ROUTINE interavm  ***
       !!----------------------------------------------------------------------  
-      INTEGER                              , INTENT(in   ) ::   i1, i2, j1, j2, k1, k2
-      REAL(wp),DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) ::   ptab
+      INTEGER                              , INTENT(in   ) ::   i1, i2, j1, j2, k1, k2, m1, m2
+      REAL(wp),DIMENSION(i1:i2,j1:j2,k1:k2,m1:m2), INTENT(inout) ::   ptab
       LOGICAL                              , INTENT(in   ) ::   before
+      REAL(wp), DIMENSION(k1:k2) :: tabin
+      REAL(wp) :: h_in(k1:k2)
+      REAL(wp) :: h_out(1:jpk)
+      REAL(wp) :: zrhoxy
+      INTEGER  :: N_in, N_out, ji, jj, jk
       !!----------------------------------------------------------------------  
       !      
-      IF( before ) THEN
-         ptab (i1:i2,j1:j2,k1:k2) = avm_k(i1:i2,j1:j2,k1:k2)
-      ELSE
-         avm  (i1:i2,j1:j2,k1:k2) = ptab (i1:i2,j1:j2,k1:k2)
+      zrhoxy = Agrif_rhox()*Agrif_rhoy()
+      IF (before) THEN         
+         DO jk=k1,k2
+            DO jj=j1,j2
+              DO ji=i1,i2
+                    ptab(ji,jj,jk,1) = avm_k(ji,jj,jk)
+              END DO
+           END DO
+        END DO
+#ifdef key_vertical         
+        DO jk=k1,k2
+           DO jj=j1,j2
+              DO ji=i1,i2
+                 ptab(ji,jj,jk,2) = wmask(ji,jj,jk) * e1e2t(ji,jj) * e3w_n(ji,jj,jk) 
+              END DO
+           END DO
+        END DO
+#else
+      ptab(i1:i2,j1:j2,k1:k2,2) = 0._wp
+#endif
+      ELSE 
+#ifdef key_vertical         
+         avm_k(i1:i2,j1:j2,1:jpk) = 0.
+         DO jj=j1,j2
+            DO ji=i1,i2
+               N_in = 0
+               DO jk=k1,k2 !k2 = jpk of parent grid
+                  IF (ptab(ji,jj,jk,2) == 0) EXIT
+                  N_in = N_in + 1
+                  tabin(jk) = ptab(ji,jj,jk,1)
+                  h_in(N_in) = ptab(ji,jj,jk,2)/(e1e2t(ji,jj)*zrhoxy)
+               END DO
+               N_out = 0
+               DO jk=1,jpk ! jpk of child grid
+                  IF (wmask(ji,jj,jk) == 0) EXIT 
+                  N_out = N_out + 1
+                  h_out(jk) = e3t_n(ji,jj,jk)
+               ENDDO
+               IF (N_in > 0) THEN
+                  CALL reconstructandremap(tabin(1:N_in),h_in,avm_k(ji,jj,1:N_out),h_out,N_in,N_out)
+               ENDIF
+            ENDDO
+         ENDDO
+#else
+         avm_k(i1:i2,j1:j2,k1:k2) = ptab (i1:i2,j1:j2,k1:k2,1)
+#endif
       ENDIF
       !