Changeset 4486 for branches/2013/dev_MERGE_2013/NEMOGCM/NEMO/OPA_SRC/DYN/dynspg_ts.F90

Timestamp:

2014-02-05T12:23:56+01:00 (10 years ago)

Author:

jchanut

Message:

Finalize Time split and AGRIF (tickets #106 and #107) + ticket #1240

File:

• branches/2013/dev_MERGE_2013/NEMOGCM/NEMO/OPA_SRC/DYN/dynspg_ts.F90 (modified) (12 diffs)

branches/2013/dev_MERGE_2013/NEMOGCM/NEMO/OPA_SRC/DYN/dynspg_ts.F90

@@ -143 +143 @@
       REAL(wp) ::   zx1, zy1, zx2, zy2         !   -      -
       REAL(wp) ::   z1_12, z1_8, z1_4, z1_2    !   -      -
-      REAL(wp) ::   zu_spg, zv_spg                !   -      -
+      REAL(wp) ::   zu_spg, zv_spg             !   -      -
       REAL(wp) ::   zhura, zhvra               !   -      -
       REAL(wp) ::   za0, za1, za2, za3           !   -      -
@@ -204 +204 @@
          CALL ts_wgt(ln_bt_av, ll_fw_start, icycle, wgtbtp1, wgtbtp2)
          !
-         IF ((neuler/=0).AND.(ln_bt_fw)) CALL ts_rst( nit000, 'READ' ) 
          !
       ENDIF
@@ -449 +448 @@
          zssh_frc(:,:) = zssh_frc(:,:) + ssh_iau(:,:)
       ENDIF
+#endif
+      !                                   !* Fill boundary data arrays with AGRIF
+      !                                   ! -------------------------------------
+#if defined key_agrif
+         IF( .NOT.Agrif_Root() ) CALL agrif_dta_ts( kt )
 #endif
       !
@@ -543 +547 @@
          zwx(:,:) = e2u(:,:) * ua_e(:,:) * zhup2_e(:,:)         ! fluxes at jn+0.5
          zwy(:,:) = e1v(:,:) * va_e(:,:) * zhvp2_e(:,:)
+         !
+#if defined key_agrif
+         ! Set fluxes during predictor step to ensure 
+         ! volume conservation
+         IF( (.NOT.Agrif_Root()).AND.ln_bt_fw ) THEN
+            IF((nbondi == -1).OR.(nbondi == 2)) THEN
+               DO jj=1,jpj
+                  zwx(2,jj) = ubdy_w(jj) * e2u(2,jj)
+               END DO
+            ENDIF
+            IF((nbondi ==  1).OR.(nbondi == 2)) THEN
+               DO jj=1,jpj
+                  zwx(nlci-2,jj) = ubdy_e(jj) * e2u(nlci-2,jj)
+               END DO
+            ENDIF
+            IF((nbondj == -1).OR.(nbondj == 2)) THEN
+               DO ji=1,jpi
+                  zwy(ji,2) = vbdy_s(ji) * e1v(ji,2)
+               END DO
+            ENDIF
+            IF((nbondj ==  1).OR.(nbondj == 2)) THEN
+               DO ji=1,jpi
+                  zwy(ji,nlcj-2) = vbdy_n(ji) * e1v(ji,nlcj-2)
+               END DO
+            ENDIF
+         ENDIF
+#endif
+         !
+         ! Sum over sub-time-steps to compute advective velocities
+         za2 = wgtbtp2(jn)
+         zu_sum  (:,:) = zu_sum  (:,:) + za2 * zwx  (:,:) / e2u  (:,:)
+         zv_sum  (:,:) = zv_sum  (:,:) + za2 * zwy  (:,:) / e1v  (:,:)
+         !
+         ! Set next sea level:
          DO jj = 2, jpjm1                                 
             DO ji = fs_2, fs_jpim1   ! vector opt.
@@ -549 +587 @@
             END DO
          END DO
-         !
-         ! Sum over sub-time-steps to compute advective velocities
-         za2 = wgtbtp2(jn)
-         zu_sum  (:,:) = zu_sum  (:,:) + za2 * ua_e  (:,:) * zhup2_e  (:,:)
-         zv_sum  (:,:) = zv_sum  (:,:) + za2 * va_e  (:,:) * zhvp2_e  (:,:)
-         !
-         ! Set next sea level:
          ssha_e(:,:) = (  sshn_e(:,:) - rdtbt * ( zssh_frc(:,:) + zhdiv(:,:) )  ) * tmask(:,:,1)
          CALL lbc_lnk( ssha_e, 'T',  1._wp )
@@ -754 +785 @@
          IF( lk_bdy ) CALL bdy_dyn2d( jn, ua_e, va_e, zun_e, zvn_e, hur_e, hvr_e, ssha_e )
 #endif
-#if defined key_agrif
-         IF( .NOT.Agrif_Root() ) CALL agrif_dyn_ts( kt, jn ) ! Agrif
+#if defined key_agrif                                                           
+         IF( .NOT.Agrif_Root() )  CALL agrif_dyn_ts( jn )  ! Agrif
 #endif
          !                                             !* Swap
@@ -842 +873 @@
          !
       END DO
+      !
+#if defined key_agrif
+      ! Save time integrated fluxes during child grid integration
+      ! (used to update coarse grid transports)
+      ! Useless with 2nd order momentum schemes
+      !
+      IF ( (.NOT.Agrif_Root()).AND.(ln_bt_fw) ) THEN
+         IF ( Agrif_NbStepint().EQ.0 ) THEN
+            ub2_i_b(:,:) = 0.e0
+            vb2_i_b(:,:) = 0.e0
+         END IF
+         !
+         za1 = 1._wp / REAL(Agrif_rhot(), wp)
+         ub2_i_b(:,:) = ub2_i_b(:,:) + za1 * ub2_b(:,:)
+         vb2_i_b(:,:) = vb2_i_b(:,:) + za1 * vb2_b(:,:)
+      ENDIF
+      !
+      !
+#endif      
       !
       !                                   !* write time-spliting arrays in the restart
@@ -956 +1006 @@
             CALL iom_get( numror, jpdom_autoglo, 'vb_e'     ,    vb_e(:,:) )
          ENDIF
+#if defined key_agrif
+         ! Read time integrated fluxes
+         IF ( .NOT.Agrif_Root() ) THEN
+            CALL iom_get( numror, jpdom_autoglo, 'ub2_i_b'  , ub2_i_b(:,:) )   
+            CALL iom_get( numror, jpdom_autoglo, 'vb2_i_b'  , vb2_i_b(:,:) )
+         ENDIF
+#endif
       !
       ELSEIF( TRIM(cdrw) == 'WRITE' ) THEN
@@ -969 +1026 @@
             CALL iom_rstput( kt, nitrst, numrow, 'vb_e'     ,    vb_e(:,:) )
          ENDIF
+#if defined key_agrif
+         ! Save time integrated fluxes
+         IF ( .NOT.Agrif_Root() ) THEN
+            CALL iom_rstput( kt, nitrst, numrow, 'ub2_i_b'  , ub2_i_b(:,:) )
+            CALL iom_rstput( kt, nitrst, numrow, 'vb2_i_b'  , vb2_i_b(:,:) )
+         ENDIF
+#endif
       ENDIF
       !
@@ -1053 +1117 @@
       ENDIF
       !
+#if defined key_agrif
+      ! Restrict the use of Agrif to the forward case only
+      IF ((.NOT.ln_bt_fw ).AND.(.NOT.Agrif_Root())) CALL ctl_stop( 'AGRIF not implemented if ln_bt_fw=.FALSE.' )
+#endif
+      !
       IF(lwp) WRITE(numout,*)    '     Time filter choice, nn_bt_flt: ', nn_bt_flt
       SELECT CASE ( nn_bt_flt )
@@ -1072 +1141 @@
          CALL ctl_stop( 'dynspg_ts ERROR: Maximum Courant number is greater than 0.9: Inc. nn_baro !' )          
       ENDIF
+      !
+      IF ((neuler/=0).AND.(ln_bt_fw)) CALL ts_rst( nit000, 'READ' ) 
       !
       CALL wrk_dealloc( jpi, jpj, zcu )
@@ -1104 +1175 @@
 
 
+