Skip to content
Snippets Groups Projects

Compare revisions

Changes are shown as if the source revision was being merged into the target revision. Learn more about comparing revisions.

Source

Select target project
No results found

Target

Select target project
No results found
Show changes
Commits on Source (41)
Hide whitespace changes
Inline Side-by-side
Showing
with 439 additions and 130 deletions
arch/arch-osx_gfortran.fcm
View file @ da31f855
...@@ -35,7 +35,7 @@ ...@@ -35,7 +35,7 @@
%CPP cpp -Dkey_nosignedzero %CPP cpp -Dkey_nosignedzero
%FC mpif90 %FC mpif90
%FCFLAGS -fdefault-real-8 -O3 -funroll-all-loops -fcray-pointer -ffree-line-length-none %FCFLAGS -fdefault-real-8 -O3 -funroll-all-loops -fcray-pointer -ffree-line-length-none -fallow-argument-mismatch
%FFLAGS %FCFLAGS %FFLAGS %FCFLAGS
%LD %FC %LD %FC
%LDFLAGS %LDFLAGS
......
cfgs/SHARED/grid_def_nemo.xml
View file @ da31f855
...@@ -12,7 +12,7 @@ ...@@ -12,7 +12,7 @@
<domain domain_ref="grid_T" /> <domain domain_ref="grid_T" />
</grid> </grid>
<grid id="grid_T_2D_inner" > <grid id="grid_T_2D_inner" >
<domain domain_ref="grid_T_inner" /> <domain domain_ref="grid_T_inner" name="grid_T" />
</grid> </grid>
<!-- --> <!-- -->
<grid id="grid_T_ncatice" > <grid id="grid_T_ncatice" >
...@@ -20,7 +20,7 @@ ...@@ -20,7 +20,7 @@
<axis axis_ref="ncatice" /> <axis axis_ref="ncatice" />
</grid> </grid>
<grid id="grid_T_ncatice_inner" > <grid id="grid_T_ncatice_inner" >
<domain domain_ref="grid_T_inner" /> <domain domain_ref="grid_T_inner" name="grid_T" />
<axis axis_ref="ncatice" /> <axis axis_ref="ncatice" />
</grid> </grid>
<!-- --> <!-- -->
...@@ -29,7 +29,7 @@ ...@@ -29,7 +29,7 @@
<axis axis_ref="deptht" /> <axis axis_ref="deptht" />
</grid> </grid>
<grid id="grid_T_3D_inner" > <grid id="grid_T_3D_inner" >
<domain domain_ref="grid_T_inner" /> <domain domain_ref="grid_T_inner" name="grid_T" />
<axis axis_ref="deptht" /> <axis axis_ref="deptht" />
</grid> </grid>
<!-- --> <!-- -->
...@@ -38,7 +38,7 @@ ...@@ -38,7 +38,7 @@
<axis axis_ref="profsed" /> <axis axis_ref="profsed" />
</grid> </grid>
<grid id="grid_T_3DS_inner" > <grid id="grid_T_3DS_inner" >
<domain domain_ref="grid_T_inner" /> <domain domain_ref="grid_T_inner" name="grid_T" />
<axis axis_ref="profsed" /> <axis axis_ref="profsed" />
</grid> </grid>
<!-- --> <!-- -->
...@@ -46,7 +46,7 @@ ...@@ -46,7 +46,7 @@
<domain domain_ref="grid_U" /> <domain domain_ref="grid_U" />
</grid> </grid>
<grid id="grid_U_2D_inner" > <grid id="grid_U_2D_inner" >
<domain domain_ref="grid_U_inner" /> <domain domain_ref="grid_U_inner" name="grid_U" />
</grid> </grid>
<!-- --> <!-- -->
<grid id="grid_U_3D" > <grid id="grid_U_3D" >
...@@ -54,7 +54,7 @@ ...@@ -54,7 +54,7 @@
<axis axis_ref="depthu" /> <axis axis_ref="depthu" />
</grid> </grid>
<grid id="grid_U_3D_inner" > <grid id="grid_U_3D_inner" >
<domain domain_ref="grid_U_inner" /> <domain domain_ref="grid_U_inner" name="grid_U" />
<axis axis_ref="depthu" /> <axis axis_ref="depthu" />
</grid> </grid>
<!-- --> <!-- -->
...@@ -62,7 +62,7 @@ ...@@ -62,7 +62,7 @@
<domain domain_ref="grid_V" /> <domain domain_ref="grid_V" />
</grid> </grid>
<grid id="grid_V_2D_inner" > <grid id="grid_V_2D_inner" >
<domain domain_ref="grid_V_inner" /> <domain domain_ref="grid_V_inner" name="grid_V" />
</grid> </grid>
<!-- --> <!-- -->
<grid id="grid_V_3D" > <grid id="grid_V_3D" >
...@@ -70,7 +70,7 @@ ...@@ -70,7 +70,7 @@
<axis axis_ref="depthv" /> <axis axis_ref="depthv" />
</grid> </grid>
<grid id="grid_V_3D_inner" > <grid id="grid_V_3D_inner" >
<domain domain_ref="grid_V_inner" /> <domain domain_ref="grid_V_inner" name="grid_V" />
<axis axis_ref="depthv" /> <axis axis_ref="depthv" />
</grid> </grid>
<!-- --> <!-- -->
...@@ -78,7 +78,7 @@ ...@@ -78,7 +78,7 @@
<domain domain_ref="grid_W" /> <domain domain_ref="grid_W" />
</grid> </grid>
<grid id="grid_W_2D_inner" > <grid id="grid_W_2D_inner" >
<domain domain_ref="grid_W_inner" /> <domain domain_ref="grid_W_inner" name="grid_W" />
</grid> </grid>
<!-- --> <!-- -->
<grid id="grid_W_3D" > <grid id="grid_W_3D" >
...@@ -86,7 +86,7 @@ ...@@ -86,7 +86,7 @@
<axis axis_ref="depthw" /> <axis axis_ref="depthw" />
</grid> </grid>
<grid id="grid_W_3D_inner" > <grid id="grid_W_3D_inner" >
<domain domain_ref="grid_W_inner" /> <domain domain_ref="grid_W_inner" name="grid_W" />
<axis axis_ref="depthw" /> <axis axis_ref="depthw" />
</grid> </grid>
<!-- --> <!-- -->
...@@ -94,7 +94,7 @@ ...@@ -94,7 +94,7 @@
<domain domain_ref="grid_F" /> <domain domain_ref="grid_F" />
</grid> </grid>
<grid id="grid_F_2D_inner" > <grid id="grid_F_2D_inner" >
<domain domain_ref="grid_F_inner" /> <domain domain_ref="grid_F_inner" name="grid_F" />
</grid> </grid>
<!-- --> <!-- -->
<grid id="grid_F_3D" > <grid id="grid_F_3D" >
...@@ -102,7 +102,7 @@ ...@@ -102,7 +102,7 @@
<axis axis_ref="depthf" /> <axis axis_ref="depthf" />
</grid> </grid>
<grid id="grid_F_3D_inner" > <grid id="grid_F_3D_inner" >
<domain domain_ref="grid_F_inner" /> <domain domain_ref="grid_F_inner" name="grid_F" />
<axis axis_ref="depthf" /> <axis axis_ref="depthf" />
</grid> </grid>
<!-- --> <!-- -->
......
sette/sette.sh
View file @ da31f855
...@@ -221,7 +221,7 @@ if [ ${USING_LOOP_FUSION} == "no" ] ; then export DEL_KEYS="${DEL_KEYS}key_loop ...@@ -221,7 +221,7 @@ if [ ${USING_LOOP_FUSION} == "no" ] ; then export DEL_KEYS="${DEL_KEYS}key_loop
if [ ${USING_QCO} == "yes" ] ; then export ADD_KEYS="${ADD_KEYS}key_qco " ; fi if [ ${USING_QCO} == "yes" ] ; then export ADD_KEYS="${ADD_KEYS}key_qco " ; fi
if [ ${USING_QCO} == "no" ] ; then export DEL_KEYS="${DEL_KEYS}key_qco key_linssh " ; fi if [ ${USING_QCO} == "no" ] ; then export DEL_KEYS="${DEL_KEYS}key_qco key_linssh " ; fi
# #
if [ ${USING_RK3} == "yes" ] ; then export ADD_KEYS="${ADD_KEYS}key_qco key_RK3 " ; fi if [ ${USING_RK3} == "yes" ] ; then export ADD_KEYS="${ADD_KEYS}key_RK3 " ; fi
if [ ${USING_RK3} == "no" ] ; then export DEL_KEYS="${DEL_KEYS}key_RK3 " ; fi if [ ${USING_RK3} == "no" ] ; then export DEL_KEYS="${DEL_KEYS}key_RK3 " ; fi
# #
......
src/NST/agrif_all_update.F90
View file @ da31f855
...@@ -99,9 +99,13 @@ CONTAINS ...@@ -99,9 +99,13 @@ CONTAINS
& uu_b(:,:, Kbb_a), 'U',-1._wp, & & uu_b(:,:, Kbb_a), 'U',-1._wp, &
& vv_b(:,:, Kmm_a), 'V',-1._wp, & & vv_b(:,:, Kmm_a), 'V',-1._wp, &
& vv_b(:,:, Kbb_a), 'V',-1._wp, & & vv_b(:,:, Kbb_a), 'V',-1._wp, &
# if ! defined key_RK3
& ub2_b(:,:), 'U',-1._wp, & & ub2_b(:,:), 'U',-1._wp, &
& ub2_i_b(:,:), 'U',-1._wp, & & un_bf(:,:), 'U',-1._wp, &
& vb2_b(:,:), 'V',-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 ) & vb2_i_b(:,:), 'V',-1._wp )
#if defined key_qco #if defined key_qco
......
src/NST/agrif_oce_interp.F90
View file @ da31f855
...@@ -28,7 +28,6 @@ MODULE agrif_oce_interp ...@@ -28,7 +28,6 @@ MODULE agrif_oce_interp
USE zdf_oce USE zdf_oce
USE agrif_oce USE agrif_oce
USE phycst USE phycst
!!! USE dynspg_ts, ONLY: un_adv, vn_adv
! !
USE in_out_manager USE in_out_manager
USE agrif_oce_sponge USE agrif_oce_sponge
...@@ -167,18 +166,36 @@ CONTAINS ...@@ -167,18 +166,36 @@ CONTAINS
END SUBROUTINE Agrif_istate_ssh END SUBROUTINE Agrif_istate_ssh
SUBROUTINE Agrif_tra SUBROUTINE Agrif_tra( kt, kstg )
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
!! *** ROUTINE Agrif_tra *** !! *** ROUTINE Agrif_tra ***
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
INTEGER, INTENT(in) :: kt
INTEGER, OPTIONAL, INTENT(in) :: kstg
REAL(wp) :: ztindex
! !
IF( Agrif_Root() ) RETURN IF( Agrif_Root() ) RETURN
! !
! Set time index depending on stage in case of RK3 time stepping:
IF ( PRESENT( kstg ) ) THEN
ztindex = REAL(Agrif_Nbstepint(), wp)
IF ( kstg == 1 ) THEN
ztindex = ztindex + 1._wp / 3._wp
ELSEIF ( kstg == 2 ) THEN
ztindex = ztindex + 1._wp / 2._wp
ELSEIF ( kstg == 3 ) THEN
ztindex = ztindex + 1._wp
ENDIF
ztindex = ztindex / Agrif_Rhot()
ELSE
ztindex = REAL(Agrif_Nbstepint()+1, wp) / Agrif_Rhot()
ENDIF
!
Agrif_SpecialValue = 0._wp Agrif_SpecialValue = 0._wp
Agrif_UseSpecialValue = l_spc_tra Agrif_UseSpecialValue = l_spc_tra
l_vremap = ln_vert_remap l_vremap = ln_vert_remap
! !
CALL Agrif_Bc_variable( ts_interp_id, procname=interptsn ) CALL Agrif_Bc_variable( ts_interp_id, calledweight=ztindex, procname=interptsn )
! !
Agrif_UseSpecialValue = .FALSE. Agrif_UseSpecialValue = .FALSE.
l_vremap = .FALSE. l_vremap = .FALSE.
...@@ -186,35 +203,52 @@ CONTAINS ...@@ -186,35 +203,52 @@ CONTAINS
END SUBROUTINE Agrif_tra END SUBROUTINE Agrif_tra
SUBROUTINE Agrif_dyn( kt ) SUBROUTINE Agrif_dyn( kt, kstg )
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
!! *** ROUTINE Agrif_DYN *** !! *** ROUTINE Agrif_DYN ***
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
INTEGER, INTENT(in) :: kt INTEGER, INTENT(in) :: kt
INTEGER, OPTIONAL, INTENT(in) :: kstg
! !
INTEGER :: ji, jj, jk ! dummy loop indices INTEGER :: ji, jj, jk ! dummy loop indices
INTEGER :: ibdy1, jbdy1, ibdy2, jbdy2 INTEGER :: ibdy1, jbdy1, ibdy2, jbdy2
REAL(wp) :: zflag REAL(wp) :: zflag
REAL(wp), DIMENSION(jpi,jpj) :: zub, zvb REAL(wp), DIMENSION(jpi,jpj) :: zub, zvb
REAL(wp), DIMENSION(jpi,jpj) :: zhub, zhvb REAL(wp), DIMENSION(jpi,jpj) :: zhub, zhvb
REAL(wp) :: ztindex
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
! !
IF( Agrif_Root() ) RETURN IF( Agrif_Root() ) RETURN
! !
! Set time index depending on stage in case of RK3 time stepping:
IF ( PRESENT( kstg ) ) THEN
ztindex = REAL(Agrif_Nbstepint(), wp)
IF ( kstg == 1 ) THEN
ztindex = ztindex + 1._wp / 3._wp
ELSEIF ( kstg == 2 ) THEN
ztindex = ztindex + 1._wp / 2._wp
ELSEIF ( kstg == 3 ) THEN
ztindex = ztindex + 1._wp
ENDIF
ztindex = ztindex / Agrif_Rhot()
ELSE
ztindex = REAL(Agrif_Nbstepint()+1, wp) / Agrif_Rhot()
ENDIF
!
Agrif_SpecialValue = 0.0_wp Agrif_SpecialValue = 0.0_wp
Agrif_UseSpecialValue = ln_spc_dyn Agrif_UseSpecialValue = ln_spc_dyn
l_vremap = ln_vert_remap l_vremap = ln_vert_remap
! !
use_sign_north = .TRUE. use_sign_north = .TRUE.
sign_north = -1.0_wp sign_north = -1.0_wp
CALL Agrif_Bc_variable( un_interp_id, procname=interpun ) CALL Agrif_Bc_variable( un_interp_id, calledweight=ztindex, procname=interpun )
CALL Agrif_Bc_variable( vn_interp_id, procname=interpvn ) CALL Agrif_Bc_variable( vn_interp_id, calledweight=ztindex, procname=interpvn )
IF( .NOT.ln_dynspg_ts ) THEN ! Get transports IF( .NOT.ln_dynspg_ts ) THEN ! Get transports
ubdy(:,:) = 0._wp ; vbdy(:,:) = 0._wp ubdy(:,:) = 0._wp ; vbdy(:,:) = 0._wp
utint_stage(:,:) = 0 ; vtint_stage(:,:) = 0 utint_stage(:,:) = 0 ; vtint_stage(:,:) = 0
CALL Agrif_Bc_variable( unb_interp_id, procname=interpunb ) CALL Agrif_Bc_variable( unb_interp_id, calledweight=ztindex, procname=interpunb )
CALL Agrif_Bc_variable( vnb_interp_id, procname=interpvnb ) CALL Agrif_Bc_variable( vnb_interp_id, calledweight=ztindex, procname=interpvnb )
ENDIF ENDIF
use_sign_north = .FALSE. use_sign_north = .FALSE.
...@@ -675,6 +709,13 @@ CONTAINS ...@@ -675,6 +709,13 @@ CONTAINS
! !
IF( Agrif_Root() ) RETURN IF( Agrif_Root() ) RETURN
! !
#if defined key_RK3
Agrif_SpecialValue = 0._wp
Agrif_UseSpecialValue = .TRUE.
CALL Agrif_Bc_variable(sshn_id, procname=interpsshn )
Agrif_UseSpecialValue = .FALSE.
#endif
!
ll_int_cons = ln_bt_fw ! Assume conservative temporal integration in the forward case only ll_int_cons = ln_bt_fw ! Assume conservative temporal integration in the forward case only
! !
! Enforce volume conservation if no time refinement: ! Enforce volume conservation if no time refinement:
...@@ -1399,10 +1440,11 @@ CONTAINS ...@@ -1399,10 +1440,11 @@ CONTAINS
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
IF( before ) THEN IF( before ) THEN
! IF ( ln_bt_fw ) 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) ptab(i1:i2,j1:j2) = e2u(i1:i2,j1:j2) * ub2_b(i1:i2,j1:j2)
! ELSE # endif
! ptab(i1:i2,j1:j2) = e2u(i1:i2,j1:j2) * un_adv(i1:i2,j1:j2)
! ENDIF
ELSE ELSE
zrhot = Agrif_rhot() zrhot = Agrif_rhot()
! Time indexes bounds for integration ! Time indexes bounds for integration
...@@ -1431,12 +1473,13 @@ CONTAINS ...@@ -1431,12 +1473,13 @@ CONTAINS
REAL(wp) :: zrhoy REAL(wp) :: zrhoy
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
IF( before ) THEN 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) & ptab(i1:i2,j1:j2) = e2u(i1:i2,j1:j2) * ub2_b(i1:i2,j1:j2) &
* umask(i1:i2,j1:j2,1) * umask(i1:i2,j1:j2,1)
! ELSE # endif
! ptab(i1:i2,j1:j2) = e2u(i1:i2,j1:j2) * un_adv(i1:i2,j1:j2)
! ENDIF
ELSE ELSE
zrhoy = Agrif_Rhoy() zrhoy = Agrif_Rhoy()
! !
...@@ -1466,12 +1509,21 @@ CONTAINS ...@@ -1466,12 +1509,21 @@ CONTAINS
jmin = MAX(j1, 2) ; jmax = MIN(j2, jpj-1) jmin = MAX(j1, 2) ; jmax = MIN(j2, jpj-1)
DO ji=imin,imax DO ji=imin,imax
DO jj=jmin,jmax 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) & 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 ) &
& -vb2_b(ji-1,jj )*e1v(ji-1,jj ) ) & & -vb2_b(ji-1,jj )*e1v(ji-1,jj ) ) &
& -vmask(ji,jj-1,1) & & -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) &
& -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
END DO END DO
ELSE ELSE
...@@ -1507,11 +1559,11 @@ CONTAINS ...@@ -1507,11 +1559,11 @@ CONTAINS
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
! !
IF( before ) THEN 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) ptab(i1:i2,j1:j2) = e1v(i1:i2,j1:j2) * vb2_b(i1:i2,j1:j2)
! ELSE # endif
! ptab(i1:i2,j1:j2) = e1v(i1:i2,j1:j2) * vn_adv(i1:i2,j1:j2)
! ENDIF
ELSE ELSE
zrhot = Agrif_rhot() zrhot = Agrif_rhot()
! Time indexes bounds for integration ! Time indexes bounds for integration
...@@ -1541,12 +1593,13 @@ CONTAINS ...@@ -1541,12 +1593,13 @@ CONTAINS
REAL(wp) :: zrhox REAL(wp) :: zrhox
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
IF( before ) THEN 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) & ptab(i1:i2,j1:j2) = e1v(i1:i2,j1:j2) * vb2_b(i1:i2,j1:j2) &
* vmask(i1:i2,j1:j2,1) * vmask(i1:i2,j1:j2,1)
! ELSE # endif
! ptab(i1:i2,j1:j2) = e1v(i1:i2,j1:j2) * vn_adv(i1:i2,j1:j2)
! ENDIF
ELSE ELSE
zrhox = Agrif_Rhox() zrhox = Agrif_Rhox()
! !
...@@ -1576,12 +1629,21 @@ CONTAINS ...@@ -1576,12 +1629,21 @@ CONTAINS
jmin = MAX(j1, 2) ; jmax = MIN(j2, jpj-1) jmin = MAX(j1, 2) ; jmax = MIN(j2, jpj-1)
DO ji=imin,imax DO ji=imin,imax
DO jj=jmin,jmax 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) & 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) &
& -ub2_b(ji ,jj-1)*e2u(ji ,jj-1) ) & & -ub2_b(ji ,jj-1)*e2u(ji ,jj-1) ) &
& -umask(ji-1,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) &
& -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
END DO END DO
ELSE ELSE
......
src/NST/agrif_oce_sponge.F90
View file @ da31f855
...@@ -36,7 +36,7 @@ MODULE agrif_oce_sponge ...@@ -36,7 +36,7 @@ MODULE agrif_oce_sponge
# include "do_loop_substitute.h90" # include "do_loop_substitute.h90"
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
!! NEMO/NST 4.0 , NEMO Consortium (2018) !! NEMO/NST 4.0 , NEMO Consortium (2018)
!! $Id: agrif_oce_sponge.F90 15437 2021-10-22 12:21:20Z jchanut $ !! $Id: agrif_oce_sponge.F90 14800 2021-05-06 15:42:46Z jchanut $
!! Software governed by the CeCILL license (see ./LICENSE) !! Software governed by the CeCILL license (see ./LICENSE)
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
CONTAINS CONTAINS
...@@ -50,8 +50,12 @@ CONTAINS ...@@ -50,8 +50,12 @@ CONTAINS
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
! !
#if defined SPONGE #if defined SPONGE
#if defined key_RK3
zcoef = REAL(Agrif_Nbstepint(), wp)/REAL(Agrif_rhot())
#else
!! Assume persistence: !! Assume persistence:
zcoef = REAL(Agrif_rhot()-1,wp)/REAL(Agrif_rhot()) zcoef = REAL(Agrif_rhot()-1,wp)/REAL(Agrif_rhot())
#endif
Agrif_SpecialValue = 0._wp Agrif_SpecialValue = 0._wp
Agrif_UseSpecialValue = l_spc_tra Agrif_UseSpecialValue = l_spc_tra
...@@ -78,7 +82,12 @@ CONTAINS ...@@ -78,7 +82,12 @@ CONTAINS
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
! !
#if defined SPONGE #if defined SPONGE
#if defined key_RK3
zcoef = REAL(Agrif_Nbstepint(), wp)/REAL(Agrif_rhot())
#else
zcoef = REAL(Agrif_rhot()-1,wp)/REAL(Agrif_rhot()) zcoef = REAL(Agrif_rhot()-1,wp)/REAL(Agrif_rhot())
#endif
Agrif_SpecialValue = 0._wp Agrif_SpecialValue = 0._wp
Agrif_UseSpecialValue = ln_spc_dyn Agrif_UseSpecialValue = ln_spc_dyn
......
src/NST/agrif_oce_update.F90
View file @ da31f855
...@@ -41,7 +41,7 @@ MODULE agrif_oce_update ...@@ -41,7 +41,7 @@ MODULE agrif_oce_update
# include "domzgr_substitute.h90" # include "domzgr_substitute.h90"
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
!! NEMO/NST 4.0 , NEMO Consortium (2018) !! NEMO/NST 4.0 , NEMO Consortium (2018)
!! $Id: agrif_oce_update.F90 15317 2021-10-01 16:09:36Z jchanut $ !! $Id: agrif_oce_update.F90 14800 2021-05-06 15:42:46Z jchanut $
!! Software governed by the CeCILL license (see ./LICENSE) !! Software governed by the CeCILL license (see ./LICENSE)
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
CONTAINS CONTAINS
...@@ -98,12 +98,14 @@ CONTAINS ...@@ -98,12 +98,14 @@ CONTAINS
! !
IF ( ln_dynspg_ts .AND. ln_bt_fw ) THEN IF ( ln_dynspg_ts .AND. ln_bt_fw ) THEN
! Update time integrated transports ! Update time integrated transports
# if ! defined key_RK3
# if ! defined DECAL_FEEDBACK_2D # 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(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) CALL Agrif_Update_Variable(vb2b_update_id,locupdate1=(/ nn_shift_bar,-2/),locupdate2=(/ nn_shift_bar,-2/),procname = updatevb2b)
# else # 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(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) CALL Agrif_Update_Variable(vb2b_update_id,locupdate1=(/1+nn_shift_bar,-2/),locupdate2=(/ nn_shift_bar,-2/),procname = updatevb2b)
# endif
# endif # endif
IF (lk_agrif_fstep) THEN 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(ub2b_update_id,locupdate1=(/ nn_shift_bar+nn_dist_par_bc-1,-2/),locupdate2=(/ nn_shift_bar+nn_dist_par_bc ,-2/),procname = updateumsk)
...@@ -544,6 +546,7 @@ CONTAINS ...@@ -544,6 +546,7 @@ CONTAINS
DO jk=1,jpkm1 DO jk=1,jpkm1
DO jj=j1,j2 DO jj=j1,j2
DO ji=i1,i2 DO ji=i1,i2
#if ! defined key_RK3
IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler))) THEN ! Add asselin part IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler))) THEN ! Add asselin part
ze3b = e3u(ji,jj,jk,Kbb_a) & ! Recover e3ub before update ze3b = e3u(ji,jj,jk,Kbb_a) & ! Recover e3ub before update
& - rn_atfp * ( e3u(ji,jj,jk,Kmm_a) - e3u(ji,jj,jk,Krhs_a) ) & - rn_atfp * ( e3u(ji,jj,jk,Kmm_a) - e3u(ji,jj,jk,Krhs_a) )
...@@ -553,6 +556,7 @@ CONTAINS ...@@ -553,6 +556,7 @@ CONTAINS
uu(ji,jj,jk,Kbb_a) = ( zub + rn_atfp * ( zunu - zuno) ) & uu(ji,jj,jk,Kbb_a) = ( zub + rn_atfp * ( zunu - zuno) ) &
& * umask(ji,jj,jk) / e3u(ji,jj,jk,Kbb_a) & * umask(ji,jj,jk) / e3u(ji,jj,jk,Kbb_a)
ENDIF ENDIF
#endif
! !
uu(ji,jj,jk,Kmm_a) = tabres_child(ji,jj,jk) * umask(ji,jj,jk) uu(ji,jj,jk,Kmm_a) = tabres_child(ji,jj,jk) * umask(ji,jj,jk)
END DO END DO
...@@ -693,6 +697,7 @@ CONTAINS ...@@ -693,6 +697,7 @@ CONTAINS
DO jk=1,jpkm1 DO jk=1,jpkm1
DO jj=j1,j2 DO jj=j1,j2
DO ji=i1,i2 DO ji=i1,i2
#if ! defined key_RK3
IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler))) THEN ! Add asselin part IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler))) THEN ! Add asselin part
ze3b = e3v(ji,jj,jk,Kbb_a) & ! Recover e3vb before update ze3b = e3v(ji,jj,jk,Kbb_a) & ! Recover e3vb before update
& - rn_atfp * ( e3v(ji,jj,jk,Kmm_a) - e3v(ji,jj,jk,Krhs_a) ) & - rn_atfp * ( e3v(ji,jj,jk,Kmm_a) - e3v(ji,jj,jk,Krhs_a) )
...@@ -702,6 +707,7 @@ CONTAINS ...@@ -702,6 +707,7 @@ CONTAINS
vv(ji,jj,jk,Kbb_a) = ( zvb + rn_atfp * ( zvnu - zvno) ) & vv(ji,jj,jk,Kbb_a) = ( zvb + rn_atfp * ( zvnu - zvno) ) &
& * vmask(ji,jj,jk) / e3v(ji,jj,jk,Kbb_a) & * vmask(ji,jj,jk) / e3v(ji,jj,jk,Kbb_a)
ENDIF ENDIF
#endif
! !
vv(ji,jj,jk,Kmm_a) = tabres_child(ji,jj,jk) * vmask(ji,jj,jk) vv(ji,jj,jk,Kmm_a) = tabres_child(ji,jj,jk) * vmask(ji,jj,jk)
END DO END DO
...@@ -768,12 +774,14 @@ CONTAINS ...@@ -768,12 +774,14 @@ CONTAINS
DO ji=i1,i2 DO ji=i1,i2
! !
! Update barotropic velocities: ! Update barotropic velocities:
#if ! defined key_RK3
IF ( .NOT.ln_dynspg_ts .OR. (ln_dynspg_ts.AND.(.NOT.ln_bt_fw)) ) THEN IF ( .NOT.ln_dynspg_ts .OR. (ln_dynspg_ts.AND.(.NOT.ln_bt_fw)) ) THEN
IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler))) THEN ! Add asselin part IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler))) THEN ! Add asselin part
zcorr = (tabres(ji,jj) - uu_b(ji,jj,Kmm_a) * hu(ji,jj,Krhs_a)) * r1_hu(ji,jj,Kbb_a) zcorr = (tabres(ji,jj) - uu_b(ji,jj,Kmm_a) * hu(ji,jj,Krhs_a)) * r1_hu(ji,jj,Kbb_a)
uu_b(ji,jj,Kbb_a) = uu_b(ji,jj,Kbb_a) + rn_atfp * zcorr * umask(ji,jj,1) uu_b(ji,jj,Kbb_a) = uu_b(ji,jj,Kbb_a) + rn_atfp * zcorr * umask(ji,jj,1)
END IF END IF
ENDIF ENDIF
#endif
uu_b(ji,jj,Kmm_a) = tabres(ji,jj) * r1_hu(ji,jj,Kmm_a) * umask(ji,jj,1) uu_b(ji,jj,Kmm_a) = tabres(ji,jj) * r1_hu(ji,jj,Kmm_a) * umask(ji,jj,1)
! !
END DO END DO
...@@ -838,12 +846,14 @@ CONTAINS ...@@ -838,12 +846,14 @@ CONTAINS
DO jj=j1,j2 DO jj=j1,j2
DO ji=i1,i2 DO ji=i1,i2
! Update barotropic velocities: ! Update barotropic velocities:
#if ! defined key_RK3
IF ( .NOT.ln_dynspg_ts .OR. (ln_dynspg_ts.AND.(.NOT.ln_bt_fw)) ) THEN IF ( .NOT.ln_dynspg_ts .OR. (ln_dynspg_ts.AND.(.NOT.ln_bt_fw)) ) THEN
IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler))) THEN ! Add asselin part IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler))) THEN ! Add asselin part
zcorr = (tabres(ji,jj) - vv_b(ji,jj,Kmm_a) * hv(ji,jj,Krhs_a)) * r1_hv(ji,jj,Kbb_a) zcorr = (tabres(ji,jj) - vv_b(ji,jj,Kmm_a) * hv(ji,jj,Krhs_a)) * r1_hv(ji,jj,Kbb_a)
vv_b(ji,jj,Kbb_a) = vv_b(ji,jj,Kbb_a) + rn_atfp * zcorr * vmask(ji,jj,1) vv_b(ji,jj,Kbb_a) = vv_b(ji,jj,Kbb_a) + rn_atfp * zcorr * vmask(ji,jj,1)
END IF END IF
ENDIF ENDIF
#endif
vv_b(ji,jj,Kmm_a) = tabres(ji,jj) * r1_hv(ji,jj,Kmm_a) * vmask(ji,jj,1) vv_b(ji,jj,Kmm_a) = tabres(ji,jj) * r1_hv(ji,jj,Kmm_a) * vmask(ji,jj,1)
! !
END DO END DO
...@@ -903,6 +913,7 @@ CONTAINS ...@@ -903,6 +913,7 @@ CONTAINS
END DO END DO
ELSE ELSE
! !
#if ! defined key_RK3
IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler))) THEN IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler))) THEN
DO jj=j1,j2 DO jj=j1,j2
DO ji=i1,i2 DO ji=i1,i2
...@@ -911,6 +922,7 @@ CONTAINS ...@@ -911,6 +922,7 @@ CONTAINS
END DO END DO
END DO END DO
ENDIF ENDIF
#endif
! !
DO jj=j1,j2 DO jj=j1,j2
DO ji=i1,i2 DO ji=i1,i2
...@@ -977,7 +989,7 @@ CONTAINS ...@@ -977,7 +989,7 @@ CONTAINS
! !
END SUBROUTINE updatevmsk END SUBROUTINE updatevmsk
# if ! defined key_RK3
SUBROUTINE updateub2b( tabres, i1, i2, j1, j2, before ) SUBROUTINE updateub2b( tabres, i1, i2, j1, j2, before )
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
!! *** ROUTINE updateub2b *** !! *** ROUTINE updateub2b ***
...@@ -1013,6 +1025,7 @@ CONTAINS ...@@ -1013,6 +1025,7 @@ CONTAINS
ENDIF ENDIF
! !
END SUBROUTINE updateub2b END SUBROUTINE updateub2b
# endif
SUBROUTINE reflux_sshu( tabres, i1, i2, j1, j2, before, nb, ndir ) SUBROUTINE reflux_sshu( tabres, i1, i2, j1, j2, before, nb, ndir )
!!--------------------------------------------- !!---------------------------------------------
...@@ -1041,16 +1054,28 @@ CONTAINS ...@@ -1041,16 +1054,28 @@ CONTAINS
! !
IF (western_side) THEN IF (western_side) THEN
DO jj=j1,j2 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)) 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 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 IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler))) ssh(i1 ,jj,Kbb_a) = ssh(i1 ,jj,Kbb_a) + rn_atfp * zcor
#endif
END DO END DO
ENDIF ENDIF
IF (eastern_side) THEN IF (eastern_side) THEN
DO jj=j1,j2 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)) 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 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 IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler))) ssh(i2+1,jj,Kbb_a) = ssh(i2+1,jj,Kbb_a) + rn_atfp * zcor
#endif
END DO END DO
ENDIF ENDIF
! !
...@@ -1058,6 +1083,7 @@ CONTAINS ...@@ -1058,6 +1083,7 @@ CONTAINS
! !
END SUBROUTINE reflux_sshu END SUBROUTINE reflux_sshu
# if ! defined key_RK3
SUBROUTINE updatevb2b( tabres, i1, i2, j1, j2, before ) SUBROUTINE updatevb2b( tabres, i1, i2, j1, j2, before )
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
!! *** ROUTINE updatevb2b *** !! *** ROUTINE updatevb2b ***
...@@ -1093,6 +1119,7 @@ CONTAINS ...@@ -1093,6 +1119,7 @@ CONTAINS
ENDIF ENDIF
! !
END SUBROUTINE updatevb2b END SUBROUTINE updatevb2b
# endif
SUBROUTINE reflux_sshv( tabres, i1, i2, j1, j2, before, nb, ndir ) SUBROUTINE reflux_sshv( tabres, i1, i2, j1, j2, before, nb, ndir )
!!--------------------------------------------- !!---------------------------------------------
...@@ -1121,16 +1148,28 @@ CONTAINS ...@@ -1121,16 +1148,28 @@ CONTAINS
! !
IF (southern_side) THEN IF (southern_side) THEN
DO ji=i1,i2 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)) 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 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 IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler))) ssh(ji,j1 ,Kbb_a) = ssh(ji,j1,Kbb_a) + rn_atfp * zcor
#endif
END DO END DO
ENDIF ENDIF
IF (northern_side) THEN IF (northern_side) THEN
DO ji=i1,i2 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)) 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 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 IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler))) ssh(ji,j2+1,Kbb_a) = ssh(ji,j2+1,Kbb_a) + rn_atfp * zcor
#endif
END DO END DO
ENDIF ENDIF
! !
...@@ -1232,6 +1271,7 @@ CONTAINS ...@@ -1232,6 +1271,7 @@ CONTAINS
! of prognostic variables ! of prognostic variables
e3t(i1:i2,j1:j2,1:jpkm1,Krhs_a) = e3t(i1:i2,j1:j2,1:jpkm1,Kmm_a) e3t(i1:i2,j1:j2,1:jpkm1,Krhs_a) = e3t(i1:i2,j1:j2,1:jpkm1,Kmm_a)
#if ! defined key_RK3
IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler) )) THEN IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler) )) THEN
DO jk = 1, jpkm1 DO jk = 1, jpkm1
DO jj=j1,j2 DO jj=j1,j2
...@@ -1262,6 +1302,7 @@ CONTAINS ...@@ -1262,6 +1302,7 @@ CONTAINS
END DO END DO
! !
ENDIF ENDIF
#endif
! !
! 2) Updates at NOW time step: ! 2) Updates at NOW time step:
! ---------------------------- ! ----------------------------
......
src/NST/agrif_top_interp.F90
View file @ da31f855
...@@ -30,23 +30,42 @@ MODULE agrif_top_interp ...@@ -30,23 +30,42 @@ MODULE agrif_top_interp
# include "domzgr_substitute.h90" # include "domzgr_substitute.h90"
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
!! NEMO/NST 4.0 , NEMO Consortium (2018) !! NEMO/NST 4.0 , NEMO Consortium (2018)
!! $Id: agrif_top_interp.F90 14218 2020-12-18 16:44:52Z jchanut $ !! $Id: agrif_top_interp.F90 14800 2021-05-06 15:42:46Z jchanut $
!! Software governed by the CeCILL license (see ./LICENSE) !! Software governed by the CeCILL license (see ./LICENSE)
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
CONTAINS CONTAINS
SUBROUTINE Agrif_trc SUBROUTINE Agrif_trc( kt, kstg )
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
!! *** ROUTINE Agrif_trc *** !! *** ROUTINE Agrif_trc ***
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
INTEGER, INTENT(in) :: kt
INTEGER, OPTIONAL, INTENT(in) :: kstg
!
REAL(wp) :: ztindex
! !
IF( Agrif_Root() ) RETURN IF( Agrif_Root() ) RETURN
! !
! Set time index depending on stage in case of RK3 time stepping:
IF ( PRESENT( kstg ) ) THEN
ztindex = REAL(Agrif_Nbstepint(), wp)
IF ( kstg == 1 ) THEN
ztindex = ztindex + 1._wp / 3._wp
ELSEIF ( kstg == 2 ) THEN
ztindex = ztindex + 1._wp / 2._wp
ELSEIF ( kstg == 3 ) THEN
ztindex = ztindex + 1._wp
ENDIF
ztindex = ztindex / Agrif_Rhot()
ELSE
ztindex = REAL(Agrif_Nbstepint()+1, wp) / Agrif_Rhot()
ENDIF
!
Agrif_SpecialValue = 0._wp Agrif_SpecialValue = 0._wp
Agrif_UseSpecialValue = l_spc_top Agrif_UseSpecialValue = l_spc_top
l_vremap = ln_vert_remap l_vremap = ln_vert_remap
! !
CALL Agrif_Bc_variable( trn_id, procname=interptrn ) CALL Agrif_Bc_variable( trn_id,calledweight=ztindex, procname=interptrn )
! !
Agrif_UseSpecialValue = .FALSE. Agrif_UseSpecialValue = .FALSE.
l_vremap = .FALSE. l_vremap = .FALSE.
......
src/NST/agrif_top_sponge.F90
View file @ da31f855
...@@ -33,7 +33,7 @@ MODULE agrif_top_sponge ...@@ -33,7 +33,7 @@ MODULE agrif_top_sponge
# include "domzgr_substitute.h90" # include "domzgr_substitute.h90"
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
!! NEMO/NST 4.0 , NEMO Consortium (2018) !! NEMO/NST 4.0 , NEMO Consortium (2018)
!! $Id: agrif_top_sponge.F90 15437 2021-10-22 12:21:20Z jchanut $ !! $Id: agrif_top_sponge.F90 14800 2021-05-06 15:42:46Z jchanut $
!! Software governed by the CeCILL license (see ./LICENSE) !! Software governed by the CeCILL license (see ./LICENSE)
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
CONTAINS CONTAINS
...@@ -46,9 +46,12 @@ CONTAINS ...@@ -46,9 +46,12 @@ CONTAINS
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
! !
#if defined SPONGE_TOP #if defined SPONGE_TOP
#if defined key_RK3
zcoef = REAL(Agrif_Nbstepint(), wp)/REAL(Agrif_rhot())
#else
!! Assume persistence: !! Assume persistence:
zcoef = REAL(Agrif_rhot()-1,wp)/REAL(Agrif_rhot()) zcoef = REAL(Agrif_rhot()-1,wp)/REAL(Agrif_rhot())
#endif
Agrif_SpecialValue = 0._wp Agrif_SpecialValue = 0._wp
Agrif_UseSpecialValue = l_spc_top Agrif_UseSpecialValue = l_spc_top
l_vremap = ln_vert_remap l_vremap = ln_vert_remap
......
src/NST/agrif_top_update.F90
View file @ da31f855
...@@ -29,7 +29,7 @@ MODULE agrif_top_update ...@@ -29,7 +29,7 @@ MODULE agrif_top_update
# include "domzgr_substitute.h90" # include "domzgr_substitute.h90"
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
!! NEMO/NST 4.0 , NEMO Consortium (2018) !! NEMO/NST 4.0 , NEMO Consortium (2018)
!! $Id: agrif_top_update.F90 15265 2021-09-16 11:13:13Z jchanut $ !! $Id: agrif_top_update.F90 14800 2021-05-06 15:42:46Z jchanut $
!! Software governed by the CeCILL license (see ./LICENSE) !! Software governed by the CeCILL license (see ./LICENSE)
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
CONTAINS CONTAINS
...@@ -121,7 +121,7 @@ CONTAINS ...@@ -121,7 +121,7 @@ CONTAINS
ENDIF ENDIF
ENDDO ENDDO
ENDDO ENDDO
#if ! defined key_RK3
IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler))) THEN IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler))) THEN
! Add asselin part ! Add asselin part
DO jn = 1,jptra DO jn = 1,jptra
...@@ -142,6 +142,7 @@ CONTAINS ...@@ -142,6 +142,7 @@ CONTAINS
END DO END DO
END DO END DO
ENDIF ENDIF
#endif
DO jn = 1,jptra DO jn = 1,jptra
DO jk = 1, jpkm1 DO jk = 1, jpkm1
DO jj = j1, j2 DO jj = j1, j2
...@@ -160,6 +161,7 @@ CONTAINS ...@@ -160,6 +161,7 @@ CONTAINS
& * tmask(i1:i2,j1:j2,jk) & * tmask(i1:i2,j1:j2,jk)
END DO END DO
ENDDO ENDDO
#if ! defined key_RK3
IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler))) THEN IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler))) THEN
! Add asselin part ! Add asselin part
DO jn = 1,jptra DO jn = 1,jptra
...@@ -180,6 +182,7 @@ CONTAINS ...@@ -180,6 +182,7 @@ CONTAINS
END DO END DO
END DO END DO
ENDIF ENDIF
#endif
DO jn = 1,jptra DO jn = 1,jptra
DO jk=k1,k2 DO jk=k1,k2
DO jj=j1,j2 DO jj=j1,j2
......
src/NST/agrif_user.F90
View file @ da31f855
...@@ -27,7 +27,11 @@ ...@@ -27,7 +27,11 @@
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
! !
CALL nemo_init !* Initializations of each fine grid CALL nemo_init !* Initializations of each fine grid
# if defined key_RK3
Kbb_a = Nbb; Kmm_a = Nbb; Krhs_a = Nrhs
# else
Kbb_a = Nbb; Kmm_a = Nnn; Krhs_a = Nrhs ! agrif_oce module copies of time level indices Kbb_a = Nbb; Kmm_a = Nnn; Krhs_a = Nrhs ! agrif_oce module copies of time level indices
# endif
! !
! !* Agrif initialization ! !* Agrif initialization
CALL Agrif_InitValues_cont CALL Agrif_InitValues_cont
...@@ -410,29 +414,17 @@ ...@@ -410,29 +414,17 @@
hbdy(:,:) = 0._wp hbdy(:,:) = 0._wp
ssh(:,:,Krhs_a) = 0._wp ssh(:,:,Krhs_a) = 0._wp
IF ( ln_dynspg_ts ) THEN Agrif_UseSpecialValue = ln_spc_dyn
Agrif_UseSpecialValue = ln_spc_dyn use_sign_north = .TRUE.
use_sign_north = .TRUE. sign_north = -1.
sign_north = -1. ubdy(:,:) = 0._wp
CALL Agrif_Bc_variable(ub2b_interp_id,calledweight=1.,procname=interpub2b) ! must be called before unb_id to define ubdy vbdy(:,:) = 0._wp
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( unb_interp_id,calledweight=1.,procname=interpunb ) CALL Agrif_Bc_variable( vnb_interp_id,calledweight=1.,procname=interpvnb )
CALL Agrif_Bc_variable( vnb_interp_id,calledweight=1.,procname=interpvnb ) use_sign_north = .FALSE.
use_sign_north = .FALSE. ubdy(:,:) = 0._wp
ubdy(:,:) = 0._wp vbdy(:,:) = 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 = .FALSE. Agrif_UseSpecialValue = .FALSE.
l_vremap = .FALSE. l_vremap = .FALSE.
......
src/OCE/DIA/diaar5.F90
View file @ da31f855
...@@ -79,7 +79,6 @@ CONTAINS ...@@ -79,7 +79,6 @@ CONTAINS
REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: z2d, zpe ! 2D workspace REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: z2d, zpe ! 2D workspace
REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: z3d, zrhd, ztpot, zgdept ! 3D workspace (zgdept: needed to use the substitute) REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: z3d, zrhd, ztpot, zgdept ! 3D workspace (zgdept: needed to use the substitute)
REAL(wp), ALLOCATABLE, DIMENSION(:,:,:,:) :: ztsn ! 4D workspace REAL(wp), ALLOCATABLE, DIMENSION(:,:,:,:) :: ztsn ! 4D workspace
!!-------------------------------------------------------------------- !!--------------------------------------------------------------------
IF( ln_timing ) CALL timing_start('dia_ar5') IF( ln_timing ) CALL timing_start('dia_ar5')
...@@ -318,6 +317,7 @@ CONTAINS ...@@ -318,6 +317,7 @@ CONTAINS
! !
INTEGER :: ji, jj, jk INTEGER :: ji, jj, jk
REAL(wp), DIMENSION(A2D(nn_hls)) :: z2d REAL(wp), DIMENSION(A2D(nn_hls)) :: z2d
!!----------------------------------------------------------------------
z2d(:,:) = puflx(:,:,1) z2d(:,:) = puflx(:,:,1)
DO_3D( 0, 0, 0, 0, 1, jpkm1 ) DO_3D( 0, 0, 0, 0, 1, jpkm1 )
...@@ -355,7 +355,7 @@ CONTAINS ...@@ -355,7 +355,7 @@ CONTAINS
!! ** Purpose : initialization for AR5 diagnostic computation !! ** Purpose : initialization for AR5 diagnostic computation
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
INTEGER :: inum INTEGER :: inum
INTEGER :: ik, idep INTEGER :: ik
INTEGER :: ji, jj, jk ! dummy loop indices INTEGER :: ji, jj, jk ! dummy loop indices
REAL(wp) :: zztmp REAL(wp) :: zztmp
REAL(wp), ALLOCATABLE, DIMENSION(:,:,:,:) :: zsaldta ! Jan/Dec levitus salinity REAL(wp), ALLOCATABLE, DIMENSION(:,:,:,:) :: zsaldta ! Jan/Dec levitus salinity
...@@ -384,9 +384,9 @@ CONTAINS ...@@ -384,9 +384,9 @@ CONTAINS
zvol0 (:,:) = 0._wp zvol0 (:,:) = 0._wp
thick0(:,:) = 0._wp thick0(:,:) = 0._wp
DO_3D( 1, 1, 1, 1, 1, jpkm1 ) ! interpolation of salinity at the last ocean level (i.e. the partial step) DO_3D( 1, 1, 1, 1, 1, jpkm1 ) ! interpolation of salinity at the last ocean level (i.e. the partial step)
idep = tmask(ji,jj,jk) * e3t_0(ji,jj,jk) zztmp = tmask(ji,jj,jk) * e3t_0(ji,jj,jk)
zvol0 (ji,jj) = zvol0 (ji,jj) + idep * e1e2t(ji,jj) zvol0 (ji,jj) = zvol0 (ji,jj) + zztmp * e1e2t(ji,jj)
thick0(ji,jj) = thick0(ji,jj) + idep thick0(ji,jj) = thick0(ji,jj) + zztmp
END_3D END_3D
vol0 = glob_sum( 'diaar5', zvol0 ) vol0 = glob_sum( 'diaar5', zvol0 )
DEALLOCATE( zvol0 ) DEALLOCATE( zvol0 )
......
src/OCE/DIA/diaptr.F90
View file @ da31f855
...@@ -41,6 +41,7 @@ MODULE diaptr ...@@ -41,6 +41,7 @@ MODULE diaptr
END INTERFACE END INTERFACE
PUBLIC dia_ptr ! call in step module PUBLIC dia_ptr ! call in step module
PUBLIC dia_ptr_init ! call in stprk3 module
PUBLIC dia_ptr_hst ! called from tra_ldf/tra_adv routines PUBLIC dia_ptr_hst ! called from tra_ldf/tra_adv routines
REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: hstr_adv, hstr_ldf, hstr_eiv !: Heat/Salt TRansports(adv, diff, Bolus.) REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: hstr_adv, hstr_ldf, hstr_eiv !: Heat/Salt TRansports(adv, diff, Bolus.)
...@@ -65,7 +66,7 @@ MODULE diaptr ...@@ -65,7 +66,7 @@ MODULE diaptr
# include "domzgr_substitute.h90" # include "domzgr_substitute.h90"
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
!! NEMO/OCE 4.0 , NEMO Consortium (2018) !! NEMO/OCE 4.0 , NEMO Consortium (2018)
!! $Id: diaptr.F90 14834 2021-05-11 09:24:44Z hadcv $ !! $Id: diaptr.F90 15513 2021-11-15 17:31:29Z techene $
!! Software governed by the CeCILL license (see ./LICENSE) !! Software governed by the CeCILL license (see ./LICENSE)
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
CONTAINS CONTAINS
...@@ -82,7 +83,9 @@ CONTAINS ...@@ -82,7 +83,9 @@ CONTAINS
! !
IF( ln_timing ) CALL timing_start('dia_ptr') IF( ln_timing ) CALL timing_start('dia_ptr')
#if ! defined key_RK3
IF( kt == nit000 .AND. ll_init ) CALL dia_ptr_init ! -> will define l_diaptr and nbasin IF( kt == nit000 .AND. ll_init ) CALL dia_ptr_init ! -> will define l_diaptr and nbasin
#endif
! !
IF( l_diaptr ) THEN IF( l_diaptr ) THEN
! Calculate zonal integrals ! Calculate zonal integrals
......
src/OCE/DOM/dom_oce.F90
View file @ da31f855
...@@ -152,13 +152,12 @@ MODULE dom_oce ...@@ -152,13 +152,12 @@ MODULE dom_oce
! ! reference depths of cells ! ! reference depths of cells
REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: gdept_0 !: t- depth [m] REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: gdept_0 !: t- depth [m]
REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: gdepw_0 !: w- depth [m] REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: gdepw_0 !: w- depth [m]
REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: gde3w_0 !: w- depth (sum of e3w) [m]
! ! time-dependent depths of cells (domvvl) ! ! time-dependent depths of cells (domvvl)
#if defined key_qco || defined key_linssh #if defined key_qco || defined key_linssh
#else #else
REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) :: gdept, gdepw REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: gde3w_0, gde3w !: w- depth (sum of e3w) [m]
REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: gde3w REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) :: gdept, gdepw
#endif #endif
! ! reference heights of ocean water column and its inverse ! ! reference heights of ocean water column and its inverse
REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: ht_0, r1_ht_0 !: t-depth [m] and [1/m] REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: ht_0, r1_ht_0 !: t-depth [m] and [1/m]
...@@ -286,7 +285,7 @@ CONTAINS ...@@ -286,7 +285,7 @@ CONTAINS
& ff_f (jpi,jpj) , ff_t (jpi,jpj) , STAT=ierr(ii) ) & ff_f (jpi,jpj) , ff_t (jpi,jpj) , STAT=ierr(ii) )
! !
ii = ii+1 ii = ii+1
ALLOCATE( gdept_0 (jpi,jpj,jpk) , gdepw_0 (jpi,jpj,jpk) , gde3w_0(jpi,jpj,jpk) , & ALLOCATE( gdept_0 (jpi,jpj,jpk) , gdepw_0 (jpi,jpj,jpk) , &
& gdept_1d( jpk) , gdepw_1d( jpk) , STAT=ierr(ii) ) & gdept_1d( jpk) , gdepw_1d( jpk) , STAT=ierr(ii) )
! !
ii = ii+1 ii = ii+1
...@@ -311,7 +310,8 @@ CONTAINS ...@@ -311,7 +310,8 @@ CONTAINS
#else #else
! vvl : time varation for all vertical coordinate variables ! vvl : time varation for all vertical coordinate variables
ii = ii+1 ii = ii+1
ALLOCATE( gdept (jpi,jpj,jpk,jpt) , gdepw (jpi,jpj,jpk,jpt) , gde3w (jpi,jpj,jpk) , STAT=ierr(ii) ) ALLOCATE( gdept (jpi,jpj,jpk,jpt) , gdepw (jpi,jpj,jpk,jpt) , &
& gde3w_0(jpi,jpj,jpk) , gde3w (jpi,jpj,jpk) , STAT=ierr(ii) )
! !
ii = ii+1 ii = ii+1
ALLOCATE( e3t(jpi,jpj,jpk,jpt) , e3u (jpi,jpj,jpk,jpt) , e3v (jpi,jpj,jpk,jpt) , e3f(jpi,jpj,jpk) , & ALLOCATE( e3t(jpi,jpj,jpk,jpt) , e3u (jpi,jpj,jpk,jpt) , e3v (jpi,jpj,jpk,jpt) , e3f(jpi,jpj,jpk) , &
......
src/OCE/DOM/domain.F90
View file @ da31f855
...@@ -64,7 +64,7 @@ MODULE domain ...@@ -64,7 +64,7 @@ MODULE domain
# include "do_loop_substitute.h90" # include "do_loop_substitute.h90"
!!------------------------------------------------------------------------- !!-------------------------------------------------------------------------
!! NEMO/OCE 4.0 , NEMO Consortium (2018) !! NEMO/OCE 4.0 , NEMO Consortium (2018)
!! $Id: domain.F90 15270 2021-09-17 14:27:55Z smasson $ !! $Id: domain.F90 14547 2021-02-25 17:07:15Z techene $
!! Software governed by the CeCILL license (see ./LICENSE) !! Software governed by the CeCILL license (see ./LICENSE)
!!------------------------------------------------------------------------- !!-------------------------------------------------------------------------
CONTAINS CONTAINS
...@@ -88,7 +88,6 @@ CONTAINS ...@@ -88,7 +88,6 @@ CONTAINS
! !
INTEGER :: ji, jj, jk, jt ! dummy loop indices INTEGER :: ji, jj, jk, jt ! dummy loop indices
INTEGER :: iconf = 0 ! local integers INTEGER :: iconf = 0 ! local integers
REAL(wp):: zrdt
CHARACTER (len=64) :: cform = "(A12, 3(A13, I7))" CHARACTER (len=64) :: cform = "(A12, 3(A13, I7))"
INTEGER , DIMENSION(jpi,jpj) :: ik_top , ik_bot ! top and bottom ocean level INTEGER , DIMENSION(jpi,jpj) :: ik_top , ik_bot ! top and bottom ocean level
REAL(wp), DIMENSION(jpi,jpj) :: z1_hu_0, z1_hv_0 REAL(wp), DIMENSION(jpi,jpj) :: z1_hu_0, z1_hv_0
...@@ -310,8 +309,27 @@ CONTAINS ...@@ -310,8 +309,27 @@ CONTAINS
ENDIF ENDIF
! !
! set current model timestep rDt = 2*rn_Dt if MLF or rDt = rn_Dt if RK3 ! set current model timestep rDt = 2*rn_Dt if MLF or rDt = rn_Dt if RK3
#if defined key_RK3
rDt = rn_Dt
r1_Dt = 1._wp / rDt
!
IF(lwp) THEN
WRITE(numout,*)
WRITE(numout,*) ' ===>>> Runge Kutta 3rd order (RK3) : rDt = ', rDt
WRITE(numout,*)
ENDIF
!
#else
rDt = 2._wp * rn_Dt rDt = 2._wp * rn_Dt
r1_Dt = 1._wp / rDt r1_Dt = 1._wp / rDt
!
IF(lwp) THEN
WRITE(numout,*)
WRITE(numout,*) ' ===>>> Modified Leap-Frog (MLF) : rDt = ', rDt
WRITE(numout,*)
ENDIF
!
#endif
! !
IF( l_SAS .AND. .NOT.ln_linssh ) THEN IF( l_SAS .AND. .NOT.ln_linssh ) THEN
CALL ctl_warn( 'SAS requires linear ssh : force ln_linssh = T' ) CALL ctl_warn( 'SAS requires linear ssh : force ln_linssh = T' )
...@@ -392,7 +410,16 @@ CONTAINS ...@@ -392,7 +410,16 @@ CONTAINS
IF( nn_wxios > 0 ) lwxios = .TRUE. !* set output file type for XIOS based on NEMO namelist IF( nn_wxios > 0 ) lwxios = .TRUE. !* set output file type for XIOS based on NEMO namelist
nxioso = nn_wxios nxioso = nn_wxios
ENDIF ENDIF
! !== Check consistency between ln_rstart and ln_1st_euler ==! (i.e. set l_1st_euler) !
#if defined key_RK3
! !== RK3: Open the restart file ==!
IF( ln_rstart ) THEN
IF(lwp) WRITE(numout,*)
IF(lwp) WRITE(numout,*) ' open the restart file'
CALL rst_read_open
ENDIF
#else
! !== MLF: Check consistency between ln_rstart and ln_1st_euler ==! (i.e. set l_1st_euler)
l_1st_euler = ln_1st_euler l_1st_euler = ln_1st_euler
! !
IF( ln_rstart ) THEN !* Restart case IF( ln_rstart ) THEN !* Restart case
...@@ -427,6 +454,7 @@ CONTAINS ...@@ -427,6 +454,7 @@ CONTAINS
IF(lwp) WRITE(numout,*)' an Euler initial time step is used : l_1st_euler is forced to .true. ' IF(lwp) WRITE(numout,*)' an Euler initial time step is used : l_1st_euler is forced to .true. '
l_1st_euler = .TRUE. l_1st_euler = .TRUE.
ENDIF ENDIF
#endif
! !
! !== control of output frequency ==! ! !== control of output frequency ==!
! !
......
src/OCE/DOM/domqco.F90
View file @ da31f855
...@@ -39,6 +39,7 @@ MODULE domqco ...@@ -39,6 +39,7 @@ MODULE domqco
PUBLIC dom_qco_init ! called by domain.F90 PUBLIC dom_qco_init ! called by domain.F90
PUBLIC dom_qco_zgr ! called by isfcpl.F90 PUBLIC dom_qco_zgr ! called by isfcpl.F90
PUBLIC dom_qco_r3c ! called by steplf.F90 PUBLIC dom_qco_r3c ! called by steplf.F90
PUBLIC dom_qco_r3c_RK3 ! called by stprk3_stg.F90
! !!* Namelist nam_vvl ! !!* Namelist nam_vvl
LOGICAL , PUBLIC :: ln_vvl_zstar = .FALSE. ! zstar vertical coordinate LOGICAL , PUBLIC :: ln_vvl_zstar = .FALSE. ! zstar vertical coordinate
...@@ -123,7 +124,9 @@ CONTAINS ...@@ -123,7 +124,9 @@ CONTAINS
! ! Horizontal interpolation of e3t ! ! Horizontal interpolation of e3t
#if defined key_RK3 #if defined key_RK3
CALL dom_qco_r3c( ssh(:,:,Kbb), r3t(:,:,Kbb), r3u(:,:,Kbb), r3v(:,:,Kbb), r3f(:,:) ) CALL dom_qco_r3c( ssh(:,:,Kbb), r3t(:,:,Kbb), r3u(:,:,Kbb), r3v(:,:,Kbb), r3f(:,:) )
CALL dom_qco_r3c( ssh(:,:,Kmm), r3t(:,:,Kmm), r3u(:,:,Kmm), r3v(:,:,Kmm) ) r3t(:,:,Kmm) = r3t(:,:,Kbb) !!st r3 at Kmm needed to be initialised for Agrid_Grid call in nemo_gcm
r3u(:,:,Kmm) = r3u(:,:,Kbb) !! maybe we only need zeros ???
r3v(:,:,Kmm) = r3v(:,:,Kbb)
#else #else
CALL dom_qco_r3c( ssh(:,:,Kbb), r3t(:,:,Kbb), r3u(:,:,Kbb), r3v(:,:,Kbb) ) CALL dom_qco_r3c( ssh(:,:,Kbb), r3t(:,:,Kbb), r3u(:,:,Kbb), r3v(:,:,Kbb) )
CALL dom_qco_r3c( ssh(:,:,Kmm), r3t(:,:,Kmm), r3u(:,:,Kmm), r3v(:,:,Kmm), r3f(:,:) ) CALL dom_qco_r3c( ssh(:,:,Kmm), r3t(:,:,Kmm), r3u(:,:,Kmm), r3v(:,:,Kmm), r3f(:,:) )
...@@ -136,6 +139,61 @@ CONTAINS ...@@ -136,6 +139,61 @@ CONTAINS
SUBROUTINE dom_qco_r3c( pssh, pr3t, pr3u, pr3v, pr3f ) SUBROUTINE dom_qco_r3c( pssh, pr3t, pr3u, pr3v, pr3f )
!!---------------------------------------------------------------------
!! *** ROUTINE r3c ***
!!
!! ** Purpose : compute the filtered ratio ssh/h_0 at t-,u-,v-,f-points
!!
!! ** Method : - compute the ssh at u- and v-points (f-point optional)
!! Vector Form : surface weighted averaging
!! Flux Form : simple averaging
!! - compute the ratio ssh/h_0 at t-,u-,v-pts, (f-pt optional)
!!----------------------------------------------------------------------
REAL(wp), DIMENSION(:,:) , INTENT(in ) :: pssh ! sea surface height [m]
REAL(wp), DIMENSION(:,:) , INTENT( out) :: pr3t, pr3u, pr3v ! ssh/h0 ratio at t-, u-, v-,points [-]
REAL(wp), DIMENSION(:,:), OPTIONAL, INTENT( out) :: pr3f ! ssh/h0 ratio at f-point [-]
!
INTEGER :: ji, jj ! dummy loop indices
!!----------------------------------------------------------------------
!
!
DO_2D_OVR( nn_hls, nn_hls, nn_hls, nn_hls )
pr3t(ji,jj) = pssh(ji,jj) * r1_ht_0(ji,jj) !== ratio at t-point ==!
END_2D
!
! !== ratio at u-,v-point ==!
!
DO_2D_OVR( nn_hls, nn_hls-1, nn_hls, nn_hls-1 )
pr3u(ji,jj) = 0.5_wp * ( e1e2t(ji ,jj) * pssh(ji ,jj) &
& + e1e2t(ji+1,jj) * pssh(ji+1,jj) ) * r1_hu_0(ji,jj) * r1_e1e2u(ji,jj)
pr3v(ji,jj) = 0.5_wp * ( e1e2t(ji,jj ) * pssh(ji,jj ) &
& + e1e2t(ji,jj+1) * pssh(ji,jj+1) ) * r1_hv_0(ji,jj) * r1_e1e2v(ji,jj)
END_2D
!
IF( .NOT.PRESENT( pr3f ) ) THEN !- lbc on ratio at u-, v-points only
IF (nn_hls==1) CALL lbc_lnk( 'dom_qco_r3c', pr3u, 'U', 1._wp, pr3v, 'V', 1._wp )
!
!
ELSE !== ratio at f-point ==!
!
DO_2D_OVR( nn_hls, nn_hls-1, nn_hls, nn_hls-1 )
! round brackets added to fix the order of floating point operations
! needed to ensure halo 1 - halo 2 compatibility
pr3f(ji,jj) = 0.25_wp * ( ( e1e2t(ji ,jj ) * pssh(ji ,jj ) &
& + e1e2t(ji+1,jj ) * pssh(ji+1,jj ) ) & ! bracket for halo 1 - halo 2 compatibility
& + ( e1e2t(ji ,jj+1) * pssh(ji ,jj+1) &
& + e1e2t(ji+1,jj+1) * pssh(ji+1,jj+1) ) & ! bracket for halo 1 - halo 2 compatibility
& ) * r1_hf_0(ji,jj) * r1_e1e2f(ji,jj)
END_2D
! ! lbc on ratio at u-,v-,f-points
IF (nn_hls==1) CALL lbc_lnk( 'dom_qco_r3c', pr3u, 'U', 1._wp, pr3v, 'V', 1._wp, pr3f, 'F', 1._wp )
!
ENDIF
!
END SUBROUTINE dom_qco_r3c
SUBROUTINE dom_qco_r3c_RK3( pssh, pr3t, pr3u, pr3v, pr3f )
!!--------------------------------------------------------------------- !!---------------------------------------------------------------------
!! *** ROUTINE r3c *** !! *** ROUTINE r3c ***
!! !!
...@@ -164,7 +222,7 @@ CONTAINS ...@@ -164,7 +222,7 @@ CONTAINS
!!st IF( ln_dynadv_vec ) THEN !- Vector Form (thickness weighted averaging) !!st IF( ln_dynadv_vec ) THEN !- Vector Form (thickness weighted averaging)
#if ! defined key_qcoTest_FluxForm #if ! defined key_qcoTest_FluxForm
! ! no 'key_qcoTest_FluxForm' : surface weighted ssh average ! ! no 'key_qcoTest_FluxForm' : surface weighted ssh average
DO_2D_OVR( nn_hls, nn_hls-1, nn_hls, nn_hls-1 ) DO_2D( 0, 0, 0, 0 )
pr3u(ji,jj) = 0.5_wp * ( e1e2t(ji ,jj) * pssh(ji ,jj) & pr3u(ji,jj) = 0.5_wp * ( e1e2t(ji ,jj) * pssh(ji ,jj) &
& + e1e2t(ji+1,jj) * pssh(ji+1,jj) ) * r1_hu_0(ji,jj) * r1_e1e2u(ji,jj) & + e1e2t(ji+1,jj) * pssh(ji+1,jj) ) * r1_hu_0(ji,jj) * r1_e1e2u(ji,jj)
pr3v(ji,jj) = 0.5_wp * ( e1e2t(ji,jj ) * pssh(ji,jj ) & pr3v(ji,jj) = 0.5_wp * ( e1e2t(ji,jj ) * pssh(ji,jj ) &
...@@ -172,52 +230,43 @@ CONTAINS ...@@ -172,52 +230,43 @@ CONTAINS
END_2D END_2D
!!st ELSE !- Flux Form (simple averaging) !!st ELSE !- Flux Form (simple averaging)
#else #else
DO_2D_OVR( nn_hls, nn_hls-1, nn_hls, nn_hls-1 ) DO_2D( 0, 0, 0, 0 )
pr3u(ji,jj) = 0.5_wp * ( pssh(ji,jj) + pssh(ji+1,jj ) ) * r1_hu_0(ji,jj) pr3u(ji,jj) = 0.5_wp * ( pssh(ji,jj) + pssh(ji+1,jj ) ) * r1_hu_0(ji,jj)
pr3v(ji,jj) = 0.5_wp * ( pssh(ji,jj) + pssh(ji ,jj+1) ) * r1_hv_0(ji,jj) pr3v(ji,jj) = 0.5_wp * ( pssh(ji,jj) + pssh(ji ,jj+1) ) * r1_hv_0(ji,jj)
END_2D END_2D
!!st ENDIF !!st ENDIF
#endif #endif
! !
IF( .NOT.PRESENT( pr3f ) ) THEN !- lbc on ratio at u-, v-points only IF( PRESENT( pr3f ) ) THEN !== ratio at f-point ==!
IF (nn_hls==1) CALL lbc_lnk( 'dom_qco_r3c', pr3u, 'U', 1._wp, pr3v, 'V', 1._wp )
!
!
ELSE !== ratio at f-point ==!
! !
!!st IF( ln_dynadv_vec ) THEN !- Vector Form (thickness weighted averaging) !!st IF( ln_dynadv_vec ) THEN !- Vector Form (thickness weighted averaging)
#if ! defined key_qcoTest_FluxForm #if ! defined key_qcoTest_FluxForm
! ! no 'key_qcoTest_FluxForm' : surface weighted ssh average ! ! no 'key_qcoTest_FluxForm' : surface weighted ssh average
DO_2D_OVR( nn_hls, nn_hls-1, nn_hls, nn_hls-1 ) DO_2D( 0, 0, 0, 0 )
! round brackets added to fix the order of floating point operations ! round brackets added to fix the order of floating point operations
! needed to ensure halo 1 - halo 2 compatibility ! needed to ensure halo 1 - halo 2 compatibility
pr3f(ji,jj) = 0.25_wp * ( ( e1e2t(ji ,jj ) * pssh(ji ,jj ) & pr3f(ji,jj) = 0.25_wp * ( ( e1e2t(ji ,jj ) * pssh(ji ,jj ) &
& + e1e2t(ji+1,jj ) * pssh(ji+1,jj ) & & + e1e2t(ji+1,jj ) * pssh(ji+1,jj ) ) & ! bracket for halo 1 - halo 2 compatibility
& ) & ! bracket for halo 1 - halo 2 compatibility & + ( e1e2t(ji ,jj+1) * pssh(ji ,jj+1) &
& + ( e1e2t(ji ,jj+1) * pssh(ji ,jj+1) & & + e1e2t(ji+1,jj+1) * pssh(ji+1,jj+1) ) & ! bracket for halo 1 - halo 2 compatibility
& + e1e2t(ji+1,jj+1) * pssh(ji+1,jj+1) &
& ) & ! bracket for halo 1 - halo 2 compatibility
& ) * r1_hf_0(ji,jj) * r1_e1e2f(ji,jj) & ) * r1_hf_0(ji,jj) * r1_e1e2f(ji,jj)
END_2D END_2D
!!st ELSE !- Flux Form (simple averaging) !!st ELSE !- Flux Form (simple averaging)
#else #else
DO_2D_OVR( nn_hls, nn_hls-1, nn_hls, nn_hls-1 ) DO_2D( 0, 0, 0, 0 )
! round brackets added to fix the order of floating point operations ! round brackets added to fix the order of floating point operations
! needed to ensure halo 1 - halo 2 compatibility ! needed to ensure halo 1 - halo 2 compatibility
pr3f(ji,jj) = 0.25_wp * ( ( pssh(ji,jj ) + pssh(ji+1,jj ) ) & pr3f(ji,jj) = 0.25_wp * ( ( pssh(ji,jj ) + pssh(ji+1,jj ) ) &
& + ( pssh(ji,jj+1) + pssh(ji+1,jj+1) & & + ( pssh(ji,jj+1) + pssh(ji+1,jj+1) ) & ! bracket for halo 1 - halo 2 compatibility
& ) & ! bracket for halo 1 - halo 2 compatibility
& ) * r1_hf_0(ji,jj) & ) * r1_hf_0(ji,jj)
END_2D END_2D
!!st ENDIF !!st ENDIF
#endif #endif
! ! lbc on ratio at u-,v-,f-points
IF (nn_hls==1) CALL lbc_lnk( 'dom_qco_r3c', pr3u, 'U', 1._wp, pr3v, 'V', 1._wp, pr3f, 'F', 1._wp )
! !
ENDIF ENDIF
! !
END SUBROUTINE dom_qco_r3c END SUBROUTINE dom_qco_r3c_RK3
SUBROUTINE qco_ctl SUBROUTINE qco_ctl
......
src/OCE/DOM/domzgr.F90
View file @ da31f855
...@@ -46,7 +46,7 @@ MODULE domzgr ...@@ -46,7 +46,7 @@ MODULE domzgr
# include "do_loop_substitute.h90" # include "do_loop_substitute.h90"
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
!! NEMO/OCE 4.0 , NEMO Consortium (2018) !! NEMO/OCE 4.0 , NEMO Consortium (2018)
!! $Id: domzgr.F90 15556 2021-11-29 15:23:06Z jchanut $ !! $Id: domzgr.F90 15157 2021-07-29 08:28:32Z techene $
!! Software governed by the CeCILL license (see ./LICENSE) !! Software governed by the CeCILL license (see ./LICENSE)
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
CONTAINS CONTAINS
...@@ -142,12 +142,17 @@ CONTAINS ...@@ -142,12 +142,17 @@ CONTAINS
CALL lbc_lnk( 'usrdef_zgr', zmsk, 'T', 1. ) ! set halos CALL lbc_lnk( 'usrdef_zgr', zmsk, 'T', 1. ) ! set halos
k_top(:,:) = k_top(:,:) * NINT( zmsk(:,:) ) k_top(:,:) = k_top(:,:) * NINT( zmsk(:,:) )
! !
!!gm to be remove when removing the OLD definition of e3 scale factors so that gde3w disappears #if ! defined key_qco && ! defined key_linssh
! OLD implementation of coordinate (not with 'key_qco' or 'key_linssh')
! gde3w_0 has to be defined
!!gm to be remove when removing the OLD definition of e3 scale factors so that gde3w_0=gdept_0
!!gm therefore gde3w_0 disappears
! Compute gde3w_0 (vertical sum of e3w) ! Compute gde3w_0 (vertical sum of e3w)
gde3w_0(:,:,1) = 0.5_wp * e3w_0(:,:,1) gde3w_0(:,:,1) = 0.5_wp * e3w_0(:,:,1)
DO jk = 2, jpk DO jk = 2, jpk
gde3w_0(:,:,jk) = gde3w_0(:,:,jk-1) + e3w_0(:,:,jk) gde3w_0(:,:,jk) = gde3w_0(:,:,jk-1) + e3w_0(:,:,jk)
END DO END DO
#endif
! !
! Any closed seas (defined by closea_mask > 0 in domain_cfg file) to be filled ! Any closed seas (defined by closea_mask > 0 in domain_cfg file) to be filled
! in at runtime if ln_closea=.false. ! in at runtime if ln_closea=.false.
...@@ -200,14 +205,20 @@ CONTAINS ...@@ -200,14 +205,20 @@ CONTAINS
WRITE(numout,*) ' MIN val k_top ', MINVAL( k_top(:,:) ), ' MAX ', MAXVAL( k_top(:,:) ) WRITE(numout,*) ' MIN val k_top ', MINVAL( k_top(:,:) ), ' MAX ', MAXVAL( k_top(:,:) )
WRITE(numout,*) ' MIN val k_bot ', MINVAL( k_bot(:,:) ), ' MAX ', MAXVAL( k_bot(:,:) ) WRITE(numout,*) ' MIN val k_bot ', MINVAL( k_bot(:,:) ), ' MAX ', MAXVAL( k_bot(:,:) )
WRITE(numout,*) ' MIN val depth t ', MINVAL( gdept_0(:,:,:) ), & WRITE(numout,*) ' MIN val depth t ', MINVAL( gdept_0(:,:,:) ), &
& ' w ', MINVAL( gdepw_0(:,:,:) ), '3w ', MINVAL( gde3w_0(:,:,:) ) #if ! defined key_qco && ! defined key_linssh
& '3w ', MINVAL( gde3w_0(:,:,:) ), &
#endif
& ' w ', MINVAL( gdepw_0(:,:,:) )
WRITE(numout,*) ' MIN val e3 t ', MINVAL( e3t_0(:,:,:) ), ' f ', MINVAL( e3f_0(:,:,:) ), & WRITE(numout,*) ' MIN val e3 t ', MINVAL( e3t_0(:,:,:) ), ' f ', MINVAL( e3f_0(:,:,:) ), &
& ' u ', MINVAL( e3u_0(:,:,:) ), ' u ', MINVAL( e3v_0(:,:,:) ), & & ' u ', MINVAL( e3u_0(:,:,:) ), ' u ', MINVAL( e3v_0(:,:,:) ), &
& ' uw', MINVAL( e3uw_0(:,:,:) ), ' vw', MINVAL( e3vw_0(:,:,:)), & & ' uw', MINVAL( e3uw_0(:,:,:) ), ' vw', MINVAL( e3vw_0(:,:,:)), &
& ' w ', MINVAL( e3w_0(:,:,:) ) & ' w ', MINVAL( e3w_0(:,:,:) )
WRITE(numout,*) ' MAX val depth t ', MAXVAL( gdept_0(:,:,:) ), & WRITE(numout,*) ' MAX val depth t ', MAXVAL( gdept_0(:,:,:) ), &
& ' w ', MAXVAL( gdepw_0(:,:,:) ), '3w ', MAXVAL( gde3w_0(:,:,:) ) #if ! defined key_qco && ! defined key_linssh
& '3w ', MINVAL( gde3w_0(:,:,:) ), &
#endif
& ' w ', MINVAL( gdepw_0(:,:,:) )
WRITE(numout,*) ' MAX val e3 t ', MAXVAL( e3t_0(:,:,:) ), ' f ', MAXVAL( e3f_0(:,:,:) ), & WRITE(numout,*) ' MAX val e3 t ', MAXVAL( e3t_0(:,:,:) ), ' f ', MAXVAL( e3f_0(:,:,:) ), &
& ' u ', MAXVAL( e3u_0(:,:,:) ), ' u ', MAXVAL( e3v_0(:,:,:) ), & & ' u ', MAXVAL( e3u_0(:,:,:) ), ' u ', MAXVAL( e3v_0(:,:,:) ), &
& ' uw', MAXVAL( e3uw_0(:,:,:) ), ' vw', MAXVAL( e3vw_0(:,:,:) ), & & ' uw', MAXVAL( e3uw_0(:,:,:) ), ' vw', MAXVAL( e3vw_0(:,:,:) ), &
......
src/OCE/DOM/domzgr_substitute.h90
View file @ da31f855
...@@ -51,4 +51,3 @@ ...@@ -51,4 +51,3 @@
# define gde3w(i,j,k) gdept_0(i,j,k) # define gde3w(i,j,k) gdept_0(i,j,k)
#endif #endif
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
src/OCE/DOM/istate.F90
View file @ da31f855
...@@ -50,7 +50,7 @@ MODULE istate ...@@ -50,7 +50,7 @@ MODULE istate
# include "domzgr_substitute.h90" # include "domzgr_substitute.h90"
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
!! NEMO/OCE 4.0 , NEMO Consortium (2018) !! NEMO/OCE 4.0 , NEMO Consortium (2018)
!! $Id: istate.F90 15052 2021-06-24 14:39:14Z smasson $ !! $Id: istate.F90 14991 2021-06-14 19:52:31Z techene $
!! Software governed by the CeCILL license (see ./LICENSE) !! Software governed by the CeCILL license (see ./LICENSE)
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
CONTAINS CONTAINS
...@@ -138,32 +138,47 @@ CONTAINS ...@@ -138,32 +138,47 @@ CONTAINS
ts (:,:,:,:,Kmm) = ts (:,:,:,:,Kbb) ! set now values from to before ones ts (:,:,:,:,Kmm) = ts (:,:,:,:,Kbb) ! set now values from to before ones
uu (:,:,:,Kmm) = uu (:,:,:,Kbb) uu (:,:,:,Kmm) = uu (:,:,:,Kbb)
vv (:,:,:,Kmm) = vv (:,:,:,Kbb) vv (:,:,:,Kmm) = vv (:,:,:,Kbb)
ENDIF ENDIF
#if defined key_agrif #if defined key_agrif
ENDIF ENDIF
#endif #endif
! !
! Initialize "now" and "before" barotropic velocities: #if defined key_RK3
! Do it whatever the free surface method, these arrays being eventually used IF( .NOT. ln_rstart ) THEN
#endif
! Initialize "before" barotropic velocities. "now" values are always set but
! "before" values may have been read from a restart to ensure restartability.
! In the non-restart or non-RK3 cases they need to be initialised here:
uu_b(:,:,Kbb) = 0._wp ; vv_b(:,:,Kbb) = 0._wp
DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpkm1 )
uu_b(ji,jj,Kbb) = uu_b(ji,jj,Kbb) + e3u(ji,jj,jk,Kbb) * uu(ji,jj,jk,Kbb) * umask(ji,jj,jk)
vv_b(ji,jj,Kbb) = vv_b(ji,jj,Kbb) + e3v(ji,jj,jk,Kbb) * vv(ji,jj,jk,Kbb) * vmask(ji,jj,jk)
END_3D
uu_b(:,:,Kbb) = uu_b(:,:,Kbb) * r1_hu(:,:,Kbb)
vv_b(:,:,Kbb) = vv_b(:,:,Kbb) * r1_hv(:,:,Kbb)
!
#if defined key_RK3
ENDIF
#endif
! !
uu_b(:,:,Kmm) = 0._wp ; vv_b(:,:,Kmm) = 0._wp ! Initialize "now" barotropic velocities:
uu_b(:,:,Kbb) = 0._wp ; vv_b(:,:,Kbb) = 0._wp ! Do it whatever the free surface method, these arrays being used eventually
! !
!!gm the use of umsak & vmask is not necessary below as uu(:,:,:,Kmm), vv(:,:,:,Kmm), uu(:,:,:,Kbb), vv(:,:,:,Kbb) are always masked #if defined key_RK3
IF( .NOT. ln_rstart ) THEN
uu_b(:,:,Kmm) = uu_b(:,:,Kbb) ! Kmm value set to Kbb for initialisation in Agrif_Regrid in namo_gcm
vv_b(:,:,Kmm) = vv_b(:,:,Kbb)
ENDIF
#else
!!gm the use of umask & vmask is not necessary below as uu(:,:,:,Kmm), vv(:,:,:,Kmm), uu(:,:,:,Kbb), vv(:,:,:,Kbb) are always masked
uu_b(:,:,Kmm) = 0._wp ; vv_b(:,:,Kmm) = 0._wp
DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpkm1 ) DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpkm1 )
uu_b(ji,jj,Kmm) = uu_b(ji,jj,Kmm) + e3u(ji,jj,jk,Kmm) * uu(ji,jj,jk,Kmm) * umask(ji,jj,jk) uu_b(ji,jj,Kmm) = uu_b(ji,jj,Kmm) + e3u(ji,jj,jk,Kmm) * uu(ji,jj,jk,Kmm) * umask(ji,jj,jk)
vv_b(ji,jj,Kmm) = vv_b(ji,jj,Kmm) + e3v(ji,jj,jk,Kmm) * vv(ji,jj,jk,Kmm) * vmask(ji,jj,jk) vv_b(ji,jj,Kmm) = vv_b(ji,jj,Kmm) + e3v(ji,jj,jk,Kmm) * vv(ji,jj,jk,Kmm) * vmask(ji,jj,jk)
!
uu_b(ji,jj,Kbb) = uu_b(ji,jj,Kbb) + e3u(ji,jj,jk,Kbb) * uu(ji,jj,jk,Kbb) * umask(ji,jj,jk)
vv_b(ji,jj,Kbb) = vv_b(ji,jj,Kbb) + e3v(ji,jj,jk,Kbb) * vv(ji,jj,jk,Kbb) * vmask(ji,jj,jk)
END_3D END_3D
!
uu_b(:,:,Kmm) = uu_b(:,:,Kmm) * r1_hu(:,:,Kmm) uu_b(:,:,Kmm) = uu_b(:,:,Kmm) * r1_hu(:,:,Kmm)
vv_b(:,:,Kmm) = vv_b(:,:,Kmm) * r1_hv(:,:,Kmm) vv_b(:,:,Kmm) = vv_b(:,:,Kmm) * r1_hv(:,:,Kmm)
! #endif
uu_b(:,:,Kbb) = uu_b(:,:,Kbb) * r1_hu(:,:,Kbb)
vv_b(:,:,Kbb) = vv_b(:,:,Kbb) * r1_hv(:,:,Kbb)
! !
END SUBROUTINE istate_init END SUBROUTINE istate_init
......
src/OCE/DYN/divhor.F90
View file @ da31f855
...@@ -12,6 +12,7 @@ MODULE divhor ...@@ -12,6 +12,7 @@ MODULE divhor
!! 3.7 ! 2014-01 (G. Madec) suppression of velocity curl from in-core memory !! 3.7 ! 2014-01 (G. Madec) suppression of velocity curl from in-core memory
!! - ! 2014-12 (G. Madec) suppression of cross land advection option !! - ! 2014-12 (G. Madec) suppression of cross land advection option
!! - ! 2015-10 (G. Madec) add velocity and rnf flag in argument of div_hor !! - ! 2015-10 (G. Madec) add velocity and rnf flag in argument of div_hor
!! 4.5 ! 2015-10 (S. Techene, G. Madec) hdiv replaced by e3divUh
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
...@@ -35,19 +36,89 @@ MODULE divhor ...@@ -35,19 +36,89 @@ MODULE divhor
IMPLICIT NONE IMPLICIT NONE
PRIVATE PRIVATE
PUBLIC div_hor ! routine called by step.F90 and istate.F90 ! !! * Interface
INTERFACE div_hor
MODULE PROCEDURE div_hor_RK3, div_hor_old
END INTERFACE
PUBLIC div_hor ! routine called by ssh_nxt.F90 and istate.F90
!! * Substitutions !! * Substitutions
# include "do_loop_substitute.h90" # include "do_loop_substitute.h90"
# include "domzgr_substitute.h90" # include "domzgr_substitute.h90"
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
!! NEMO/OCE 4.0 , NEMO Consortium (2018) !! NEMO/OCE 4.0 , NEMO Consortium (2018)
!! $Id: divhor.F90 15150 2021-07-27 10:38:24Z smasson $ !! $Id: divhor.F90 14808 2021-05-07 12:00:45Z jchanut $
!! Software governed by the CeCILL license (see ./LICENSE) !! Software governed by the CeCILL license (see ./LICENSE)
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
CONTAINS CONTAINS
SUBROUTINE div_hor( kt, Kbb, Kmm ) SUBROUTINE div_hor_RK3( kt, Kbb, Kmm, puu, pvv, pe3divUh )
!!----------------------------------------------------------------------
!! *** ROUTINE div_hor_RK3 ***
!!
!! ** Purpose : compute the horizontal divergence at now time-step
!!
!! ** Method : the now divergence is computed as :
!! hdiv = 1/(e1e2t*e3t) ( di[e2u*e3u un] + dj[e1v*e3v vn] )
!! and correct with runoff inflow (div_rnf) and cross land flow (div_cla)
!!
!! ** Action : - thickness weighted horizontal divergence of in input velocity (puu,pvv)
!!----------------------------------------------------------------------
INTEGER , INTENT(in ) :: kt, Kbb, Kmm ! ocean time-step & time-level indices
REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in ) :: puu, pvv ! horizontal velocity
REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT( out) :: pe3divUh ! e3t*div[Uh]
!
INTEGER :: ji, jj, jk ! dummy loop indices
!!----------------------------------------------------------------------
!
IF( ln_timing ) CALL timing_start('div_hor_RK3')
!
IF( kt == nit000 ) THEN
IF(lwp) WRITE(numout,*)
IF(lwp) WRITE(numout,*) 'div_hor_RK3 : thickness weighted horizontal divergence '
IF(lwp) WRITE(numout,*) '~~~~~~~~~~~'
hdiv (:,:,:) = 0._wp ! initialize hdiv & pe3divUh for the halos and jpk level at the first time step
ENDIF
!
pe3divUh(:,:,:) = 0._wp !!gm to be applied to the halos only
!
DO_3D_OVR( nn_hls-1, nn_hls, nn_hls-1, nn_hls, 1, jpkm1 )
hdiv(ji,jj,jk) = ( e2u(ji ,jj) * e3u(ji ,jj,jk,Kmm) * puu(ji ,jj,jk) &
& - e2u(ji-1,jj) * e3u(ji-1,jj,jk,Kmm) * puu(ji-1,jj,jk) &
& + e1v(ji,jj ) * e3v(ji,jj ,jk,Kmm) * pvv(ji,jj ,jk) &
& - e1v(ji,jj-1) * e3v(ji,jj-1,jk,Kmm) * pvv(ji,jj-1,jk) ) &
& * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm)
END_3D
!
IF( ln_rnf ) CALL sbc_rnf_div( hdiv, Kmm ) !== + runoffs divergence ==!
!
#if defined key_asminc
IF( ln_sshinc .AND. ln_asmiau ) & !== + SSH assimilation increment ==!
& CALL ssh_asm_div( kt, Kbb, Kmm, hdiv )
#endif
!
IF( ln_isf ) CALL isf_hdiv( kt, Kmm, hdiv ) !== + ice-shelf mass exchange ==!
!
IF( nn_hls==1 ) CALL lbc_lnk( 'divhor', hdiv, 'T', 1._wp ) ! (no sign change)
!
!!gm Patch before suppression of hdiv from all modules that use it
! DO_3D( 0, 0, 0, 0, 1, jpkm1 ) !== e3t * Horizontal divergence ==!
! pe3divUh(ji,jj,jk) = hdiv(ji,jj,jk) * e3t(ji,jj,jk,Kmm)
! END_3D
!JC: over whole domain, and after lbclnk on hdiv to prevent from reproducibility issues
DO_3D_OVR( nn_hls-1, nn_hls, nn_hls-1, nn_hls, 1, jpkm1 )
pe3divUh(ji,jj,jk) = hdiv(ji,jj,jk) * e3t(ji,jj,jk,Kmm)
END_3D
!!gm end
!
!
IF( ln_timing ) CALL timing_stop('div_hor_RK3')
!
END SUBROUTINE div_hor_RK3
SUBROUTINE div_hor_old( kt, Kbb, Kmm )
!!---------------------------------------------------------------------- !!----------------------------------------------------------------------
!! *** ROUTINE div_hor *** !! *** ROUTINE div_hor ***
!! !!
...@@ -102,7 +173,7 @@ CONTAINS ...@@ -102,7 +173,7 @@ CONTAINS
! ! needed for ww in sshwzv ! ! needed for ww in sshwzv
IF( ln_timing ) CALL timing_stop('div_hor') IF( ln_timing ) CALL timing_stop('div_hor')
! !
END SUBROUTINE div_hor END SUBROUTINE div_hor_old
!!====================================================================== !!======================================================================
END MODULE divhor END MODULE divhor