From c5e1206970ac551a08591d2e986c8e602b4fc0c9 Mon Sep 17 00:00:00 2001
From: jchanut <jerome.chanut@mercator-ocean.fr>
Date: Fri, 8 Apr 2022 14:26:48 +0200
Subject: [PATCH] 2D mode management for AGRIF with RK3 and more
- Suppress the use ub2_b/vb2_b un_bf/vn_bf with AGRIF and RK3. Consider
un_adv/vn_adv in that case, i.e. the available time integrated flux from
before to after time step.
- From the change above, one can now remove the definition and the allocation
of ub2_b/vb2_b un_bf/vn_bf with RK3 (these are purely related to a
leapgrog time stepping).
- Move definition and allocation of un_adv/vn_adv from dynspg_ts into oce to prevent
from circular dependencies.
- Restored read/write of barotropic arrays in restart with RK3 in the
case where ln_bt_av=.false. (Demange's time stepping with continuous
in time barotropic time stepping). Aditionnal arrays need also to be
stored in case of refluxing with AGRIF.
---
src/NST/agrif_all_update.F90 | 6 ++-
src/NST/agrif_oce_interp.F90 | 52 +++++++++++++++-------
src/NST/agrif_oce_update.F90 | 25 +++++++++--
src/NST/agrif_user.F90 | 34 +++++----------
src/OCE/DYN/dynspg_ts.F90 | 84 +++++++++++++++++++++---------------
src/OCE/oce.F90 | 18 ++++----
src/OCE/stpmlf.F90 | 1 -
src/OCE/stprk3.F90 | 2 -
src/OCE/stprk3_stg.F90 | 1 -
9 files changed, 132 insertions(+), 91 deletions(-)
diff --git a/src/NST/agrif_all_update.F90 b/src/NST/agrif_all_update.F90
index 8fc617eba..753bee7cd 100644
--- a/src/NST/agrif_all_update.F90
+++ b/src/NST/agrif_all_update.F90
@@ -99,9 +99,13 @@ CONTAINS
& uu_b(:,:, Kbb_a), 'U',-1._wp, &
& vv_b(:,:, Kmm_a), 'V',-1._wp, &
& vv_b(:,:, Kbb_a), 'V',-1._wp, &
+# if ! defined key_RK3
& ub2_b(:,:), 'U',-1._wp, &
- & ub2_i_b(:,:), 'U',-1._wp, &
+ & un_bf(:,:), 'U',-1._wp, &
& vb2_b(:,:), 'V',-1._wp, &
+ & vn_bf(:,:), 'V',-1._wp, &
+# endif
+ & ub2_i_b(:,:), 'U',-1._wp, &
& vb2_i_b(:,:), 'V',-1._wp )
#if defined key_qco
diff --git a/src/NST/agrif_oce_interp.F90 b/src/NST/agrif_oce_interp.F90
index b7ff728cb..7d25ea554 100644
--- a/src/NST/agrif_oce_interp.F90
+++ b/src/NST/agrif_oce_interp.F90
@@ -28,7 +28,6 @@ MODULE agrif_oce_interp
USE zdf_oce
USE agrif_oce
USE phycst
-!!! USE dynspg_ts, ONLY: un_adv, vn_adv
!
USE in_out_manager
USE agrif_oce_sponge
@@ -1441,10 +1440,11 @@ CONTAINS
!!----------------------------------------------------------------------
IF( before ) THEN
! IF ( ln_bt_fw ) THEN
+# if defined key_RK3
+ ptab(i1:i2,j1:j2) = e2u(i1:i2,j1:j2) * un_adv(i1:i2,j1:j2)
+# else
ptab(i1:i2,j1:j2) = e2u(i1:i2,j1:j2) * ub2_b(i1:i2,j1:j2)
-! ELSE
-! ptab(i1:i2,j1:j2) = e2u(i1:i2,j1:j2) * un_adv(i1:i2,j1:j2)
-! ENDIF
+# endif
ELSE
zrhot = Agrif_rhot()
! Time indexes bounds for integration
@@ -1473,12 +1473,13 @@ CONTAINS
REAL(wp) :: zrhoy
!!----------------------------------------------------------------------
IF( before ) THEN
-! IF ( ln_bt_fw ) THEN
+# if defined key_RK3
+ ptab(i1:i2,j1:j2) = e2u(i1:i2,j1:j2) * un_adv(i1:i2,j1:j2) &
+ * umask(i1:i2,j1:j2,1)
+# else
ptab(i1:i2,j1:j2) = e2u(i1:i2,j1:j2) * ub2_b(i1:i2,j1:j2) &
* umask(i1:i2,j1:j2,1)
-! ELSE
-! ptab(i1:i2,j1:j2) = e2u(i1:i2,j1:j2) * un_adv(i1:i2,j1:j2)
-! ENDIF
+# endif
ELSE
zrhoy = Agrif_Rhoy()
!
@@ -1508,12 +1509,21 @@ CONTAINS
jmin = MAX(j1, 2) ; jmax = MIN(j2, jpj-1)
DO ji=imin,imax
DO jj=jmin,jmax
+# if defined key_RK3
+ ptab(ji,jj) = 0.25_wp *(vmask(ji,jj ,1) &
+ & * ( vn_adv(ji+1,jj )*e1v(ji+1,jj ) &
+ & -vn_adv(ji-1,jj )*e1v(ji-1,jj ) ) &
+ & -vmask(ji,jj-1,1) &
+ & * ( vn_adv(ji+1,jj-1)*e1v(ji+1,jj-1) &
+ & -vn_adv(ji-1,jj-1)*e1v(ji-1,jj-1) ) )
+# else
ptab(ji,jj) = 0.25_wp *(vmask(ji,jj ,1) &
& * ( vb2_b(ji+1,jj )*e1v(ji+1,jj ) &
& -vb2_b(ji-1,jj )*e1v(ji-1,jj ) ) &
& -vmask(ji,jj-1,1) &
& * ( vb2_b(ji+1,jj-1)*e1v(ji+1,jj-1) &
& -vb2_b(ji-1,jj-1)*e1v(ji-1,jj-1) ) )
+# endif
END DO
END DO
ELSE
@@ -1549,11 +1559,11 @@ CONTAINS
!!----------------------------------------------------------------------
!
IF( before ) THEN
-! IF ( ln_bt_fw ) THEN
+# if defined key_RK3
+ ptab(i1:i2,j1:j2) = e1v(i1:i2,j1:j2) * vn_adv(i1:i2,j1:j2)
+# else
ptab(i1:i2,j1:j2) = e1v(i1:i2,j1:j2) * vb2_b(i1:i2,j1:j2)
-! ELSE
-! ptab(i1:i2,j1:j2) = e1v(i1:i2,j1:j2) * vn_adv(i1:i2,j1:j2)
-! ENDIF
+# endif
ELSE
zrhot = Agrif_rhot()
! Time indexes bounds for integration
@@ -1583,12 +1593,13 @@ CONTAINS
REAL(wp) :: zrhox
!!----------------------------------------------------------------------
IF( before ) THEN
-! IF ( ln_bt_fw ) THEN
+# if defined key_RK3
+ ptab(i1:i2,j1:j2) = e1v(i1:i2,j1:j2) * vn_adv(i1:i2,j1:j2) &
+ * vmask(i1:i2,j1:j2,1)
+# else
ptab(i1:i2,j1:j2) = e1v(i1:i2,j1:j2) * vb2_b(i1:i2,j1:j2) &
* vmask(i1:i2,j1:j2,1)
-! ELSE
-! ptab(i1:i2,j1:j2) = e1v(i1:i2,j1:j2) * vn_adv(i1:i2,j1:j2)
-! ENDIF
+# endif
ELSE
zrhox = Agrif_Rhox()
!
@@ -1618,12 +1629,21 @@ CONTAINS
jmin = MAX(j1, 2) ; jmax = MIN(j2, jpj-1)
DO ji=imin,imax
DO jj=jmin,jmax
+# if defined key_RK3
+ ptab(ji,jj) = 0.25_wp *(umask(ji ,jj,1) &
+ & * ( un_adv(ji ,jj+1)*e2u(ji ,jj+1) &
+ & -un_adv(ji ,jj-1)*e2u(ji ,jj-1) ) &
+ & -umask(ji-1,jj,1) &
+ & * ( un_adv(ji-1,jj+1)*e2u(ji-1,jj+1) &
+ & -un_adv(ji-1,jj-1)*e2u(ji-1,jj-1) ) )
+# else
ptab(ji,jj) = 0.25_wp *(umask(ji ,jj,1) &
& * ( ub2_b(ji ,jj+1)*e2u(ji ,jj+1) &
& -ub2_b(ji ,jj-1)*e2u(ji ,jj-1) ) &
& -umask(ji-1,jj,1) &
& * ( ub2_b(ji-1,jj+1)*e2u(ji-1,jj+1) &
& -ub2_b(ji-1,jj-1)*e2u(ji-1,jj-1) ) )
+# endif
END DO
END DO
ELSE
diff --git a/src/NST/agrif_oce_update.F90 b/src/NST/agrif_oce_update.F90
index 9a96064d7..72fb2eca1 100644
--- a/src/NST/agrif_oce_update.F90
+++ b/src/NST/agrif_oce_update.F90
@@ -96,22 +96,22 @@ CONTAINS
CALL Agrif_Update_Variable(vnb_update_id,locupdate1=(/1+nn_shift_bar,-2/),locupdate2=(/ nn_shift_bar,-2/),procname = updateV2d)
# endif
!
-#if ! defined key_RK3
IF ( ln_dynspg_ts .AND. ln_bt_fw ) THEN
! Update time integrated transports
+# if ! defined key_RK3
# if ! defined DECAL_FEEDBACK_2D
CALL Agrif_Update_Variable(ub2b_update_id,locupdate1=(/ nn_shift_bar,-2/),locupdate2=(/ nn_shift_bar,-2/),procname = updateub2b)
CALL Agrif_Update_Variable(vb2b_update_id,locupdate1=(/ nn_shift_bar,-2/),locupdate2=(/ nn_shift_bar,-2/),procname = updatevb2b)
# else
CALL Agrif_Update_Variable(ub2b_update_id,locupdate1=(/ nn_shift_bar,-2/),locupdate2=(/1+nn_shift_bar,-2/),procname = updateub2b)
CALL Agrif_Update_Variable(vb2b_update_id,locupdate1=(/1+nn_shift_bar,-2/),locupdate2=(/ nn_shift_bar,-2/),procname = updatevb2b)
+# endif
# endif
IF (lk_agrif_fstep) THEN
CALL Agrif_Update_Variable(ub2b_update_id,locupdate1=(/ nn_shift_bar+nn_dist_par_bc-1,-2/),locupdate2=(/ nn_shift_bar+nn_dist_par_bc ,-2/),procname = updateumsk)
CALL Agrif_Update_Variable(vb2b_update_id,locupdate1=(/ nn_shift_bar+nn_dist_par_bc ,-2/),locupdate2=(/ nn_shift_bar+nn_dist_par_bc-1,-2/),procname = updatevmsk)
ENDIF
END IF
-#endif
# if ! defined DECAL_FEEDBACK
CALL Agrif_Update_Variable(un_update_id,procname = updateU)
@@ -989,7 +989,7 @@ CONTAINS
!
END SUBROUTINE updatevmsk
-
+# if ! defined key_RK3
SUBROUTINE updateub2b( tabres, i1, i2, j1, j2, before )
!!----------------------------------------------------------------------
!! *** ROUTINE updateub2b ***
@@ -1025,6 +1025,7 @@ CONTAINS
ENDIF
!
END SUBROUTINE updateub2b
+# endif
SUBROUTINE reflux_sshu( tabres, i1, i2, j1, j2, before, nb, ndir )
!!---------------------------------------------
@@ -1053,7 +1054,11 @@ CONTAINS
!
IF (western_side) THEN
DO jj=j1,j2
+# if defined key_RK3
+ zcor = rn_Dt * r1_e1e2t(i1 ,jj) * e2u(i1,jj) * (un_adv(i1,jj)-tabres(i1,jj))
+# else
zcor = rn_Dt * r1_e1e2t(i1 ,jj) * e2u(i1,jj) * (ub2_b(i1,jj)-tabres(i1,jj))
+# endif
ssh(i1 ,jj,Kmm_a) = ssh(i1 ,jj,Kmm_a) + zcor
#if ! defined key_RK3
IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler))) ssh(i1 ,jj,Kbb_a) = ssh(i1 ,jj,Kbb_a) + rn_atfp * zcor
@@ -1062,7 +1067,11 @@ CONTAINS
ENDIF
IF (eastern_side) THEN
DO jj=j1,j2
+# if defined key_RK3
+ zcor = - rn_Dt * r1_e1e2t(i2+1,jj) * e2u(i2,jj) * (un_adv(i2,jj)-tabres(i2,jj))
+# else
zcor = - rn_Dt * r1_e1e2t(i2+1,jj) * e2u(i2,jj) * (ub2_b(i2,jj)-tabres(i2,jj))
+# endif
ssh(i2+1,jj,Kmm_a) = ssh(i2+1,jj,Kmm_a) + zcor
#if ! defined key_RK3
IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler))) ssh(i2+1,jj,Kbb_a) = ssh(i2+1,jj,Kbb_a) + rn_atfp * zcor
@@ -1074,6 +1083,7 @@ CONTAINS
!
END SUBROUTINE reflux_sshu
+# if ! defined key_RK3
SUBROUTINE updatevb2b( tabres, i1, i2, j1, j2, before )
!!----------------------------------------------------------------------
!! *** ROUTINE updatevb2b ***
@@ -1109,6 +1119,7 @@ CONTAINS
ENDIF
!
END SUBROUTINE updatevb2b
+# endif
SUBROUTINE reflux_sshv( tabres, i1, i2, j1, j2, before, nb, ndir )
!!---------------------------------------------
@@ -1137,7 +1148,11 @@ CONTAINS
!
IF (southern_side) THEN
DO ji=i1,i2
+# if defined key_RK3
+ zcor = rn_Dt * r1_e1e2t(ji,j1 ) * e1v(ji,j1 ) * (vn_adv(ji,j1)-tabres(ji,j1))
+# else
zcor = rn_Dt * r1_e1e2t(ji,j1 ) * e1v(ji,j1 ) * (vb2_b(ji,j1)-tabres(ji,j1))
+# endif
ssh(ji,j1 ,Kmm_a) = ssh(ji,j1 ,Kmm_a) + zcor
#if ! defined key_RK3
IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler))) ssh(ji,j1 ,Kbb_a) = ssh(ji,j1,Kbb_a) + rn_atfp * zcor
@@ -1146,7 +1161,11 @@ CONTAINS
ENDIF
IF (northern_side) THEN
DO ji=i1,i2
+# if defined key_RK3
+ zcor = - rn_Dt * r1_e1e2t(ji,j2+1) * e1v(ji,j2 ) * (vn_adv(ji,j2)-tabres(ji,j2))
+# else
zcor = - rn_Dt * r1_e1e2t(ji,j2+1) * e1v(ji,j2 ) * (vb2_b(ji,j2)-tabres(ji,j2))
+# endif
ssh(ji,j2+1,Kmm_a) = ssh(ji,j2+1,Kmm_a) + zcor
#if ! defined key_RK3
IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler))) ssh(ji,j2+1,Kbb_a) = ssh(ji,j2+1,Kbb_a) + rn_atfp * zcor
diff --git a/src/NST/agrif_user.F90 b/src/NST/agrif_user.F90
index 34be8e3b0..2bdc85c5d 100644
--- a/src/NST/agrif_user.F90
+++ b/src/NST/agrif_user.F90
@@ -414,29 +414,17 @@
hbdy(:,:) = 0._wp
ssh(:,:,Krhs_a) = 0._wp
- IF ( ln_dynspg_ts ) THEN
- Agrif_UseSpecialValue = ln_spc_dyn
- use_sign_north = .TRUE.
- sign_north = -1.
- CALL Agrif_Bc_variable(ub2b_interp_id,calledweight=1.,procname=interpub2b) ! must be called before unb_id to define ubdy
- CALL Agrif_Bc_variable(vb2b_interp_id,calledweight=1.,procname=interpvb2b) ! must be called before vnb_id to define vbdy
- CALL Agrif_Bc_variable( unb_interp_id,calledweight=1.,procname=interpunb )
- CALL Agrif_Bc_variable( vnb_interp_id,calledweight=1.,procname=interpvnb )
- use_sign_north = .FALSE.
- ubdy(:,:) = 0._wp
- vbdy(:,:) = 0._wp
- ELSEIF ( ln_dynspg_EXP ) THEN
- Agrif_UseSpecialValue = ln_spc_dyn
- use_sign_north = .TRUE.
- sign_north = -1.
- ubdy(:,:) = 0._wp
- vbdy(:,:) = 0._wp
- CALL Agrif_Bc_variable( unb_interp_id,calledweight=1.,procname=interpunb )
- CALL Agrif_Bc_variable( vnb_interp_id,calledweight=1.,procname=interpvnb )
- use_sign_north = .FALSE.
- ubdy(:,:) = 0._wp
- vbdy(:,:) = 0._wp
- ENDIF
+ Agrif_UseSpecialValue = ln_spc_dyn
+ use_sign_north = .TRUE.
+ sign_north = -1.
+ ubdy(:,:) = 0._wp
+ vbdy(:,:) = 0._wp
+ CALL Agrif_Bc_variable( unb_interp_id,calledweight=1.,procname=interpunb )
+ CALL Agrif_Bc_variable( vnb_interp_id,calledweight=1.,procname=interpvnb )
+ use_sign_north = .FALSE.
+ ubdy(:,:) = 0._wp
+ vbdy(:,:) = 0._wp
+
Agrif_UseSpecialValue = .FALSE.
l_vremap = .FALSE.
diff --git a/src/OCE/DYN/dynspg_ts.F90 b/src/OCE/DYN/dynspg_ts.F90
index 670277d51..83b2d99ca 100644
--- a/src/OCE/DYN/dynspg_ts.F90
+++ b/src/OCE/DYN/dynspg_ts.F90
@@ -69,8 +69,6 @@ MODULE dynspg_ts
PUBLIC dyn_spg_ts_init ! - - dyn_spg_init
PUBLIC dyn_drg_init ! called by stp2d
- !! Time filtered arrays at baroclinic time step:
- REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: un_adv , vn_adv !: Advection vel. at "now" barocl. step
!
INTEGER, SAVE :: icycle ! Number of barotropic sub-steps for each internal step nn_e <= 2.5 nn_e
REAL(wp),SAVE :: rDt_e ! Barotropic time step
@@ -103,7 +101,7 @@ CONTAINS
!!----------------------------------------------------------------------
!! *** routine dyn_spg_ts_alloc ***
!!----------------------------------------------------------------------
- INTEGER :: ierr(3)
+ INTEGER :: ierr(2)
!!----------------------------------------------------------------------
ierr(:) = 0
!
@@ -111,7 +109,6 @@ CONTAINS
IF( ln_dynvor_een .OR. ln_dynvor_eeT ) &
& ALLOCATE( ftnw(jpi,jpj) , ftne(jpi,jpj) , ftsw(jpi,jpj) , ftse(jpi,jpj), STAT=ierr(2) )
!
- ALLOCATE( un_adv(jpi,jpj), vn_adv(jpi,jpj) , STAT=ierr(3) )
!
dyn_spg_ts_alloc = MAXVAL( ierr(:) )
!
@@ -929,30 +926,29 @@ CONTAINS
!
#if defined key_agrif
-# if defined key_RK3
- ! advective transport from N to N+1 (i.e. Kbb to Kaa)
- ub2_b(:,:) = un_adv(:,:) ! Save integrated transport for AGRIF
- vb2_b(:,:) = vn_adv(:,:)
-# endif
!
! Save time integrated fluxes during child grid integration
! (used to update coarse grid transports at next time step)
!
- IF( .NOT.Agrif_Root() .AND. ln_bt_fw .AND. ln_agrif_2way ) THEN
+ IF( .NOT.Agrif_Root() .AND. ln_agrif_2way ) THEN
IF( Agrif_NbStepint() == 0 ) THEN
ub2_i_b(:,:) = 0._wp
vb2_i_b(:,:) = 0._wp
END IF
!
za1 = 1._wp / REAL(Agrif_rhot(), wp)
+# if defined key_RK3
+ ub2_i_b(:,:) = ub2_i_b(:,:) + za1 * un_adv(:,:)
+ vb2_i_b(:,:) = vb2_i_b(:,:) + za1 * vn_adv(:,:)
+# else
ub2_i_b(:,:) = ub2_i_b(:,:) + za1 * ub2_b(:,:)
vb2_i_b(:,:) = vb2_i_b(:,:) + za1 * vb2_b(:,:)
+
+# endif
ENDIF
#endif
-#if ! defined key_RK3
- ! !* MLF: write time-spliting arrays in the restart
- IF( lrst_oce .AND.ln_bt_fw ) CALL ts_rst_MLF( kt, 'WRITE' )
-#endif
+ ! !: write time-spliting arrays in the restart
+ IF( lrst_oce ) CALL ts_rst( kt, 'WRITE' )
!
IF( ln_wd_il ) DEALLOCATE( zcpx, zcpy )
IF( ln_wd_dl ) DEALLOCATE( ztwdmask, zuwdmask, zvwdmask, zuwdav2, zvwdav2 )
@@ -1041,9 +1037,9 @@ CONTAINS
END SUBROUTINE ts_wgt
- SUBROUTINE ts_rst_MLF( kt, cdrw )
+ SUBROUTINE ts_rst( kt, cdrw )
!!---------------------------------------------------------------------
- !! *** ROUTINE ts_rst_MLF ***
+ !! *** ROUTINE ts_rst ***
!!
!! ** Purpose : Read or write time-splitting arrays in restart file
!!----------------------------------------------------------------------
@@ -1053,11 +1049,15 @@ CONTAINS
!
IF( TRIM(cdrw) == 'READ' ) THEN ! Read/initialise
! ! ---------------
- IF( ln_rstart .AND. ln_bt_fw .AND. .NOT.l_1st_euler ) THEN !* Read the restart file
- CALL iom_get( numror, jpdom_auto, 'ub2_b' , ub2_b (:,:), cd_type = 'U', psgn = -1._wp )
- CALL iom_get( numror, jpdom_auto, 'vb2_b' , vb2_b (:,:), cd_type = 'V', psgn = -1._wp )
- CALL iom_get( numror, jpdom_auto, 'un_bf' , un_bf (:,:), cd_type = 'U', psgn = -1._wp )
- CALL iom_get( numror, jpdom_auto, 'vn_bf' , vn_bf (:,:), cd_type = 'V', psgn = -1._wp )
+ IF( ln_rstart .AND. .NOT.l_1st_euler ) THEN !* Read the restart file
+# if ! defined key_RK3
+ IF ( ln_bt_fw ) THEN
+ CALL iom_get( numror, jpdom_auto, 'ub2_b' , ub2_b (:,:), cd_type = 'U', psgn = -1._wp )
+ CALL iom_get( numror, jpdom_auto, 'vb2_b' , vb2_b (:,:), cd_type = 'V', psgn = -1._wp )
+ CALL iom_get( numror, jpdom_auto, 'un_bf' , un_bf (:,:), cd_type = 'U', psgn = -1._wp )
+ CALL iom_get( numror, jpdom_auto, 'vn_bf' , vn_bf (:,:), cd_type = 'V', psgn = -1._wp )
+ ENDIF
+# endif
IF( .NOT.ln_bt_av ) THEN
CALL iom_get( numror, jpdom_auto, 'sshbb_e' , sshbb_e(:,:), cd_type = 'T', psgn = 1._wp )
CALL iom_get( numror, jpdom_auto, 'ubb_e' , ubb_e(:,:), cd_type = 'U', psgn = -1._wp )
@@ -1074,15 +1074,21 @@ CONTAINS
ELSE
ub2_i_b(:,:) = 0._wp ; vb2_i_b(:,:) = 0._wp ! used in the 1st update of agrif
ENDIF
+# if defined key_RK3
+ CALL iom_get( numror, jpdom_auto, 'un_adv' , un_adv(:,:), cd_type = 'U', psgn = -1._wp )
+ CALL iom_get( numror, jpdom_auto, 'vn_adv' , vn_adv(:,:), cd_type = 'V', psgn = -1._wp )
+# endif
#endif
ELSE
! !* Start from rest or use RK3 time-step
IF(lwp) WRITE(numout,*)
IF(lwp) WRITE(numout,*) ' ==>>> start from rest: set barotropic values to 0'
+# if ! defined key_RK3
ub2_b (:,:) = 0._wp ; vb2_b (:,:) = 0._wp ! used in the 1st interpol of agrif
- un_adv (:,:) = 0._wp ; vn_adv (:,:) = 0._wp ! used in the 1st interpol of agrif
un_bf (:,:) = 0._wp ; vn_bf (:,:) = 0._wp ! used in the 1st update of agrif
-
+#else
+ un_adv (:,:) = 0._wp ; vn_adv (:,:) = 0._wp ! used in the 1st interpol of agrif
+#endif
#if defined key_agrif
ub2_i_b(:,:) = 0._wp ; vb2_i_b(:,:) = 0._wp ! used in the 1st update of agrif
#endif
@@ -1090,11 +1096,15 @@ CONTAINS
!
ELSEIF( TRIM(cdrw) == 'WRITE' ) THEN ! Create restart file
! ! -------------------
- IF(lwp) WRITE(numout,*) '---- ts_rst_MLF ----'
- CALL iom_rstput( kt, nitrst, numrow, 'ub2_b' , ub2_b (:,:) )
- CALL iom_rstput( kt, nitrst, numrow, 'vb2_b' , vb2_b (:,:) )
- CALL iom_rstput( kt, nitrst, numrow, 'un_bf' , un_bf (:,:) )
- CALL iom_rstput( kt, nitrst, numrow, 'vn_bf' , vn_bf (:,:) )
+ IF(lwp) WRITE(numout,*) '---- ts_rst ----'
+# if ! defined key_RK3
+ IF ( ln_bt_fw ) THEN
+ CALL iom_rstput( kt, nitrst, numrow, 'ub2_b' , ub2_b (:,:) )
+ CALL iom_rstput( kt, nitrst, numrow, 'vb2_b' , vb2_b (:,:) )
+ CALL iom_rstput( kt, nitrst, numrow, 'un_bf' , un_bf (:,:) )
+ CALL iom_rstput( kt, nitrst, numrow, 'vn_bf' , vn_bf (:,:) )
+ ENDIF
+# endif
!
IF (.NOT.ln_bt_av) THEN
CALL iom_rstput( kt, nitrst, numrow, 'sshbb_e' , sshbb_e(:,:) )
@@ -1110,10 +1120,14 @@ CONTAINS
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
+# if defined key_RK3
+ CALL iom_rstput( kt, nitrst, numrow, 'un_adv' , un_adv(:,:) )
+ CALL iom_rstput( kt, nitrst, numrow, 'vn_adv' , vn_adv(:,:) )
+# endif
#endif
ENDIF
!
- END SUBROUTINE ts_rst_MLF
+ END SUBROUTINE ts_rst
SUBROUTINE dyn_spg_ts_init
@@ -1160,11 +1174,16 @@ CONTAINS
IF(lwp) WRITE(numout,*) ' ln_bt_av =.false. => No time averaging of barotropic variables '
ENDIF
!
+# if ! defined key_RK3
IF(ln_bt_fw) THEN
IF(lwp) WRITE(numout,*) ' ln_bt_fw=.true. => Forward integration of barotropic variables '
ELSE
IF(lwp) WRITE(numout,*) ' ln_bt_fw =.false.=> Centred integration of barotropic variables '
ENDIF
+# else
+ ! Enforce ln_bt_fw = T with RK3
+ ln_bt_fw = .true.
+# endif
!
#if defined key_agrif
! Restrict the use of Agrif to the forward case only
@@ -1199,11 +1218,8 @@ CONTAINS
! ! Allocate time-splitting arrays
IF( dyn_spg_ts_alloc() /= 0 ) CALL ctl_stop('STOP', 'dyn_spg_init: failed to allocate dynspg_ts arrays' )
!
-#if ! defined key_RK3
- ! !* MLF: restart/initialise
- !
- CALL ts_rst_MLF( nit000, 'READ' )
-#endif
+ ! !: restart/initialise
+ CALL ts_rst( nit000, 'READ' )
!
END SUBROUTINE dyn_spg_ts_init
diff --git a/src/OCE/oce.F90 b/src/OCE/oce.F90
index 8f5e2d451..c180833bb 100644
--- a/src/OCE/oce.F90
+++ b/src/OCE/oce.F90
@@ -45,10 +45,13 @@ MODULE oce
REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hv_e !: external v-depth
REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hur_e !: inverse of u-depth
REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hvr_e !: inverse of v-depth
+ REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: un_adv, vn_adv !: Advective barotropic fluxes
+#if ! defined key_RK3
REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: ub2_b , vb2_b !: Half step fluxes (ln_bt_fw=T)
REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: un_bf , vn_bf !: Asselin filtered half step fluxes (ln_bt_fw=T)
+#endif
#if defined key_agrif
- REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: ub2_i_b, vb2_i_b !: Half step time integrated fluxes
+ REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: ub2_i_b, vb2_i_b !: agrif time integrated fluxes
#endif
!! interpolated gradient (only used in zps case)
@@ -115,23 +118,18 @@ CONTAINS
ALLOCATE( ssha_e(jpi,jpj), sshn_e(jpi,jpj), sshb_e(jpi,jpj), sshbb_e(jpi,jpj), &
& ua_e(jpi,jpj), un_e(jpi,jpj), ub_e(jpi,jpj), ubb_e(jpi,jpj), &
& va_e(jpi,jpj), vn_e(jpi,jpj), vb_e(jpi,jpj), vbb_e(jpi,jpj), &
- & hu_e(jpi,jpj), hur_e(jpi,jpj), hv_e(jpi,jpj), hvr_e(jpi,jpj), STAT=ierr(ii) )
+ & hu_e(jpi,jpj), hur_e(jpi,jpj), hv_e(jpi,jpj), hvr_e(jpi,jpj), &
+ & un_adv(jpi,jpj), vn_adv(jpi,jpj) , STAT=ierr(ii) )
!
-#if defined key_RK3
-# if defined key_agrif
- ii = ii+1 ! RK3+AGRIF: ???
- ALLOCATE( ub2_b(jpi,jpj), vb2_b(jpi,jpj) , STAT=ierr(ii) )
-# endif
-#else
+#if ! defined key_RK3
ii = ii+1 ! MLF: arrays related to Asselin filter + ???
ALLOCATE( un_bf(jpi,jpj), vn_bf(jpi,jpj), ub2_b(jpi,jpj), vb2_b(jpi,jpj) , STAT=ierr(ii) )
#endif
-
#if defined key_agrif
ii = ii+1 ! AGRIF: ???
ALLOCATE( ub2_i_b(jpi,jpj), vb2_i_b(jpi,jpj) , STAT=ierr(ii) )
#endif
- !
+ !
oce_alloc = MAXVAL( ierr )
IF( oce_alloc /= 0 ) CALL ctl_stop( 'STOP', 'oce_alloc: failed to allocate arrays' )
!
diff --git a/src/OCE/stpmlf.F90 b/src/OCE/stpmlf.F90
index 9c59c9e7e..407e38686 100644
--- a/src/OCE/stpmlf.F90
+++ b/src/OCE/stpmlf.F90
@@ -479,7 +479,6 @@ CONTAINS
!! ** Action : puu(Kmm),pvv(Kmm) updated now horizontal velocity (ln_bt_fw=F)
!! puu(Kaa),pvv(Kaa) after horizontal velocity
!!----------------------------------------------------------------------
- USE dynspg_ts, ONLY : un_adv, vn_adv ! updated Kmm barotropic transport
!!
INTEGER , INTENT(in ) :: Kmm, Kaa ! before and after time level indices
REAL(wp), DIMENSION(jpi,jpj,jpk,jpt), INTENT(inout) :: puu, pvv ! velocities
diff --git a/src/OCE/stprk3.F90 b/src/OCE/stprk3.F90
index 087b8e55d..5b1ca4f34 100644
--- a/src/OCE/stprk3.F90
+++ b/src/OCE/stprk3.F90
@@ -21,7 +21,6 @@ MODULE stprk3
USE domqco ! quasi-eulerian coordinate (dom_qco_r3c routine)
USE stprk3_stg ! RK3 stages
USE stp2d ! external mode solver
- USE dynspg_ts, ONLY: un_adv, vn_adv ! updated Kmm barotropic transport
USE trd_oce ! trends: ocean variables
USE diaptr
USE ldftra
@@ -311,7 +310,6 @@ CONTAINS
!! ** Action : puu(Kmm),pvv(Kmm) updated now horizontal velocity (ln_bt_fw=F)
!! puu(Kaa),pvv(Kaa) after horizontal velocity
!!----------------------------------------------------------------------
- USE dynspg_ts, ONLY : un_adv, vn_adv ! updated Kmm barotropic transport
!!
INTEGER , INTENT(in ) :: Kmm, Kaa ! before and after time level indices
REAL(wp), DIMENSION(jpi,jpj,jpk,jpt), INTENT(inout) :: puu, pvv ! velocities
diff --git a/src/OCE/stprk3_stg.F90 b/src/OCE/stprk3_stg.F90
index 38e90efd0..b12a1bebb 100644
--- a/src/OCE/stprk3_stg.F90
+++ b/src/OCE/stprk3_stg.F90
@@ -19,7 +19,6 @@ MODULE stprk3_stg
!!----------------------------------------------------------------------
USE step_oce ! time stepping used modules
USE domqco ! quasi-eulerian coordinate (dom_qco_r3c routine)
- USE dynspg_ts, ONLY: un_adv , vn_adv ! advective transport from N to N+1
USE bdydyn ! ocean open boundary conditions (define bdy_dyn)
USE lbclnk ! ocean lateral boundary conditions (or mpp link)
# if defined key_top
--
GitLab