Changeset 11802 for NEMO/branches/2019/dev_r11233_AGRIF-05_jchanut_vert_coord_interp/src/NST/agrif_oce_interp.F90

Timestamp:

2019-10-25T17:15:20+02:00 (5 years ago)

Author:

jchanut

Message:

#2222, 1) add linear interpolation in vremap module.
2) Switch remapping of viscosity from polynomial to linear.
3) Move to truly volume weighted averages for parent to child update.

File:

• NEMO/branches/2019/dev_r11233_AGRIF-05_jchanut_vert_coord_interp/src/NST/agrif_oce_interp.F90 (modified) (5 diffs)

NEMO/branches/2019/dev_r11233_AGRIF-05_jchanut_vert_coord_interp/src/NST/agrif_oce_interp.F90

@@ -719 +719 @@
                tsa(ji,jj,:,:) = 0._wp
                N_in = mbkt_parent(ji,jj)
+               IF ( tmask(ji,jj,1) == 0._wp) N_in = 0
                zhtot = 0._wp
                DO jk=1,N_in !k2 = jpk of parent grid
@@ -833 +834 @@
                N_in = mbku_parent(ji,jj)
                zhtot = 0._wp
+               IF ( umask(ji,jj,1) == 0._wp) N_in = 0
                DO jk=1,N_in
                   IF (jk==N_in) THEN
@@ -928 +930 @@
                va(ji,jj,:) = 0._wp
                N_in = mbkv_parent(ji,jj)
+               IF ( vmask(ji,jj,1) == 0._wp) N_in = 0
                zhtot = 0._wp
                DO jk=1,N_in
@@ -1276 +1279 @@
       REAL(wp),DIMENSION(i1:i2,j1:j2,k1:k2,m1:m2), INTENT(inout) ::   ptab
       LOGICAL                                    , INTENT(in   ) ::   before
-      REAL(wp), DIMENSION(k1:k2) :: tabin, h_in
-      REAL(wp), DIMENSION(1:jpk) :: h_out
-      INTEGER  :: N_in, N_out, ji, jj, jk
+      !
+      INTEGER  :: ji, jj, jk
+      INTEGER  :: N_in, N_out
+      REAL(wp), DIMENSION(k1:k2) :: tabin, z_in
+      REAL(wp), DIMENSION(1:jpk) :: z_out
       !!----------------------------------------------------------------------  
       !      
@@ -1289 +1294 @@
            END DO
         END DO
-#ifdef key_vertical         
+
+# if defined key_vertical
+        ! Interpolate thicknesses
+        ! Warning: these are masked, hence extrapolated prior interpolation.
         DO jk=k1,k2
            DO jj=j1,j2
               DO ji=i1,i2
-                 ptab(ji,jj,jk,2) = wmask(ji,jj,jk) * e3w_n(ji,jj,jk) 
+                  ptab(ji,jj,jk,2) = tmask(ji,jj,jk) * e3t_n(ji,jj,jk)
               END DO
            END DO
         END DO
-#endif
+
+        ! Extrapolate thicknesses in partial bottom cells:
+        ! Set them to Agrif_SpecialValue (0.). Correct bottom thicknesses are retrieved later on
+        IF (ln_zps) THEN
+           DO jj=j1,j2
+              DO ji=i1,i2
+                  jk = mbkt(ji,jj)
+                  ptab(ji,jj,jk,2) = 0._wp
+              END DO
+           END DO           
+        END IF
+     
+        ! Save ssh at last level:
+        IF (.NOT.ln_linssh) THEN
+           ptab(i1:i2,j1:j2,k2,2) = sshn(i1:i2,j1:j2)*tmask(i1:i2,j1:j2,1) 
+        ELSE
+           ptab(i1:i2,j1:j2,k2,2) = 0._wp
+        END IF      
+# 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)
+         IF (ln_linssh) ptab(i1:i2,j1:j2,k2,2) = 0._wp 
+         avm_k(i1:i2,j1:j2,k1:k2) = 0._wp
+            
+         DO jj = j1, j2
+            DO ji =i1, i2
+               N_in = mbkt_parent(ji,jj)
+               IF ( tmask(ji,jj,1) == 0._wp) N_in = 0
+               z_in(N_in+1) = ht0_parent(ji,jj) + ptab(ji,jj,k2,2)
+               DO jk = N_in, 1, -1  ! Parent vertical grid               
+                     z_in(jk) = z_in(jk+1) - ptab(ji,jj,jk,2)
+                    tabin(jk) = ptab(ji,jj,jk,1)
                END DO
                N_out = 0
-               DO jk=1,jpk ! jpk of child grid
-                  IF (wmask(ji,jj,jk) == 0) EXIT 
+          z_out(1) = 0._wp 
+               DO jk = 2, jpk       ! Child vertical grid
+                  IF (tmask(ji,jj,jk) == 0._wp) EXIT 
                   N_out = N_out + 1
-                  h_out(jk) = e3t_n(ji,jj,jk)
+                  z_out(jk) = z_out(jk-1) + e3t_n(ji,jj,jk-1)
                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,1)
+               IF (N_in*N_out > 0) THEN
+                  CALL remap_linear(tabin(1:N_in),z_in(1:N_in),avm_k(ji,jj,1:N_out),z_out(1:N_out),N_in,N_out,1)
                ENDIF
             ENDDO