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 (23)
Showing
with 289 additions and 123 deletions
......@@ -27,7 +27,6 @@
<!-- Files definition -->
<file_definition src="./file_def_nemo-oce.xml"/> <!-- NEMO ocean dynamics -->
<file_definition src="./file_def_nemo-ice.xml"/> <!-- NEMO ocean sea ice -->
<file_definition src="./file_def_nemo-innerttrc.xml"/> <!-- NEMO ocean inert passive tracer -->
<!-- Axis definition -->
......
......@@ -85,7 +85,7 @@
! Type of air-sea fluxes
ln_blk = .true. ! Bulk formulation (T => fill namsbc_blk )
! Sea-ice :
nn_ice = 2 ! =0 no ice boundary condition
nn_ice = 0 ! =0 no ice boundary condition
! ! =1 use observed ice-cover ( => fill namsbc_iif )
! ! =2 or 3 for SI3 and CICE, respectively
! Misc. options of sbc :
......@@ -361,6 +361,10 @@
!-----------------------------------------------------------------------
&namzdf_tke ! turbulent eddy kinetic dependent vertical diffusion (ln_zdftke =T)
!-----------------------------------------------------------------------
nn_mxlice = 0 ! type of scaling under sea-ice
! ! = 0 no scaling under sea-ice
! ! = 1 scaling with constant sea-ice thickness
! ! = 2 scaling with mean sea-ice thickness ( only with SI3 sea-ice model )
/
!!======================================================================
!! *** Diagnostics namelists *** !!
......
......@@ -298,6 +298,10 @@
!-----------------------------------------------------------------------
&namzdf_tke ! turbulent eddy kinetic dependent vertical diffusion (ln_zdftke =T)
!-----------------------------------------------------------------------
nn_mxlice = 2 ! type of scaling under sea-ice
! ! = 0 no scaling under sea-ice
! ! = 1 scaling with constant sea-ice thickness
! ! = 2 scaling with mean sea-ice thickness ( only with SI3 sea-ice model )
/
!!======================================================================
!! *** Diagnostics namelists *** !!
......
......@@ -299,6 +299,10 @@
!-----------------------------------------------------------------------
&namzdf_tke ! turbulent eddy kinetic dependent vertical diffusion (ln_zdftke =T)
!-----------------------------------------------------------------------
nn_mxlice = 2 ! type of scaling under sea-ice
! ! = 0 no scaling under sea-ice
! ! = 1 scaling with constant sea-ice thickness
! ! = 2 scaling with mean sea-ice thickness ( only with SI3 sea-ice model )
/
!!======================================================================
!! *** Diagnostics namelists *** !!
......
......@@ -357,6 +357,10 @@
!-----------------------------------------------------------------------
&namzdf_tke ! turbulent eddy kinetic dependent vertical diffusion (ln_zdftke =T)
!-----------------------------------------------------------------------
nn_mxlice = 2 ! type of scaling under sea-ice
! ! = 0 no scaling under sea-ice
! ! = 1 scaling with constant sea-ice thickness
! ! = 2 scaling with mean sea-ice thickness ( only with SI3 sea-ice model )
/
!!======================================================================
!! *** Diagnostics namelists *** !!
......
......@@ -395,6 +395,10 @@
! ! = 2 add a tke source just at the base of the ML
! ! = 3 as = 1 applied on HF part of the stress (ln_cpl=T)
ln_mxhsw = .false. ! surface mixing length scale = F(wave height)
nn_mxlice = 2 ! type of scaling under sea-ice
! ! = 0 no scaling under sea-ice
! ! = 1 scaling with constant sea-ice thickness
! ! = 2 scaling with mean sea-ice thickness ( only with SI3 sea-ice model )
/
!-----------------------------------------------------------------------
&namzdf_iwm ! internal wave-driven mixing parameterization (ln_zdfiwm =T)
......
......@@ -1199,7 +1199,7 @@
! ! = 2 first vertical derivative of mixing length bounded by 1
! ! = 3 as =2 with distinct dissipative an mixing length scale
ln_mxl0 = .true. ! surface mixing length scale = F(wind stress) (T) or not (F)
nn_mxlice = 2 ! type of scaling under sea-ice
nn_mxlice = 0 ! type of scaling under sea-ice
! ! = 0 no scaling under sea-ice
! ! = 1 scaling with constant sea-ice thickness
! ! = 2 scaling with mean sea-ice thickness ( only with SI3 sea-ice model )
......@@ -1245,7 +1245,7 @@
! ! = 1 roughness uses rn_hsri and is weigthed by 1-TANH(10*fr_i)
! ! = 2 roughness uses rn_hsri and is weighted by 1-fr_i
! ! = 3 roughness uses rn_hsri and is weighted by 1-MIN(1,4*fr_i)
nn_mxlice = 1 ! mixing under sea ice
nn_mxlice = 0 ! mixing under sea ice
! = 0 No scaling under sea-ice
! = 1 scaling with constant Ice-ocean roughness (rn_hsri)
! = 2 scaling with mean sea-ice thickness
......
......@@ -552,6 +552,10 @@
!-----------------------------------------------------------------------
&namzdf_tke ! turbulent eddy kinetic dependent vertical diffusion (ln_zdftke =T)
!-----------------------------------------------------------------------
nn_mxlice = 2 ! type of scaling under sea-ice
! ! = 0 no scaling under sea-ice
! ! = 1 scaling with constant sea-ice thickness
! ! = 2 scaling with mean sea-ice thickness ( only with SI3 sea-ice model )
/
!-----------------------------------------------------------------------
&namzdf_gls ! GLS vertical diffusion (ln_zdfgls =T)
......
......@@ -13,7 +13,7 @@
! ! = 2 first vertical derivative of mixing length bounded by 1
! ! = 3 as =2 with distinct dissipative an mixing length scale
ln_mxl0 = .true. ! surface mixing length scale = F(wind stress) (T) or not (F)
nn_mxlice = 2 ! type of scaling under sea-ice
nn_mxlice = 0 ! type of scaling under sea-ice
! ! = 0 no scaling under sea-ice
! ! = 1 scaling with constant sea-ice thickness
! ! = 2 scaling with mean sea-ice thickness ( only with SI3 sea-ice model )
......
......@@ -691,49 +691,52 @@ CONTAINS
CALL prt_ctl_info(' ========== ')
CALL prt_ctl_info(' - Cell values : ')
CALL prt_ctl_info(' ~~~~~~~~~~~~~ ')
CALL prt_ctl(tab2d_1=e1e2t , clinfo1=' cell area :')
CALL prt_ctl(tab2d_1=at_i , clinfo1=' at_i :')
CALL prt_ctl(tab2d_1=ato_i , clinfo1=' ato_i :')
CALL prt_ctl(tab2d_1=vt_i , clinfo1=' vt_i :')
CALL prt_ctl(tab2d_1=vt_s , clinfo1=' vt_s :')
CALL prt_ctl(tab2d_1=divu_i , clinfo1=' divu_i :')
CALL prt_ctl(tab2d_1=delta_i , clinfo1=' delta_i :')
CALL prt_ctl(tab2d_1=stress1_i , clinfo1=' stress1_i :')
CALL prt_ctl(tab2d_1=stress2_i , clinfo1=' stress2_i :')
CALL prt_ctl(tab2d_1=stress12_i , clinfo1=' stress12_i :')
CALL prt_ctl(tab2d_1=strength , clinfo1=' strength :')
CALL prt_ctl(tab2d_1=delta_i , clinfo1=' delta_i :')
CALL prt_ctl(tab2d_1=u_ice , clinfo1=' u_ice :', tab2d_2=v_ice , clinfo2=' v_ice :')
CALL prt_ctl(tab2d_1=e1e2t , clinfo1=' cell area :', mask1=tmask)
CALL prt_ctl(tab2d_1=at_i , clinfo1=' at_i :', mask1=tmask)
CALL prt_ctl(tab2d_1=ato_i , clinfo1=' ato_i :', mask1=tmask)
CALL prt_ctl(tab2d_1=vt_i , clinfo1=' vt_i :', mask1=tmask)
CALL prt_ctl(tab2d_1=vt_s , clinfo1=' vt_s :', mask1=tmask)
CALL prt_ctl(tab2d_1=divu_i , clinfo1=' divu_i :', mask1=tmask)
CALL prt_ctl(tab2d_1=delta_i , clinfo1=' delta_i :', mask1=tmask)
CALL prt_ctl(tab2d_1=stress1_i , clinfo1=' stress1_i :', mask1=tmask)
CALL prt_ctl(tab2d_1=stress2_i , clinfo1=' stress2_i :', mask1=tmask)
CALL prt_ctl(tab2d_1=stress12_i , clinfo1=' stress12_i :') ! should be fmask
CALL prt_ctl(tab2d_1=strength , clinfo1=' strength :', mask1=tmask)
CALL prt_ctl(tab2d_1=delta_i , clinfo1=' delta_i :', mask1=tmask)
CALL prt_ctl(tab2d_1=u_ice , clinfo1=' u_ice :', mask1=umask, &
& tab2d_2=v_ice , clinfo2=' v_ice :', mask2=vmask)
DO jl = 1, jpl
CALL prt_ctl_info(' ')
CALL prt_ctl_info(' - Category : ', ivar=jl)
CALL prt_ctl_info(' ~~~~~~~~~~')
CALL prt_ctl(tab2d_1=h_i (:,:,jl) , clinfo1= ' h_i : ')
CALL prt_ctl(tab2d_1=h_s (:,:,jl) , clinfo1= ' h_s : ')
CALL prt_ctl(tab2d_1=t_su (:,:,jl) , clinfo1= ' t_su : ')
CALL prt_ctl(tab2d_1=t_s (:,:,1,jl) , clinfo1= ' t_snow : ')
CALL prt_ctl(tab2d_1=s_i (:,:,jl) , clinfo1= ' s_i : ')
CALL prt_ctl(tab2d_1=o_i (:,:,jl) , clinfo1= ' o_i : ')
CALL prt_ctl(tab2d_1=a_i (:,:,jl) , clinfo1= ' a_i : ')
CALL prt_ctl(tab2d_1=v_i (:,:,jl) , clinfo1= ' v_i : ')
CALL prt_ctl(tab2d_1=v_s (:,:,jl) , clinfo1= ' v_s : ')
CALL prt_ctl(tab2d_1=e_s (:,:,1,jl) , clinfo1= ' e_snow : ')
CALL prt_ctl(tab2d_1=sv_i (:,:,jl) , clinfo1= ' sv_i : ')
CALL prt_ctl(tab2d_1=oa_i (:,:,jl) , clinfo1= ' oa_i : ')
CALL prt_ctl(tab2d_1=h_i (:,:,jl) , clinfo1= ' h_i : ', mask1=tmask)
CALL prt_ctl(tab2d_1=h_s (:,:,jl) , clinfo1= ' h_s : ', mask1=tmask)
CALL prt_ctl(tab2d_1=t_su(:,:,jl) , clinfo1= ' t_su : ', mask1=tmask)
CALL prt_ctl(tab2d_1=t_s (:,:,1,jl), clinfo1= ' t_snow : ', mask1=tmask)
CALL prt_ctl(tab2d_1=s_i (:,:,jl) , clinfo1= ' s_i : ', mask1=tmask)
CALL prt_ctl(tab2d_1=o_i (:,:,jl) , clinfo1= ' o_i : ', mask1=tmask)
CALL prt_ctl(tab2d_1=a_i (:,:,jl) , clinfo1= ' a_i : ', mask1=tmask)
CALL prt_ctl(tab2d_1=v_i (:,:,jl) , clinfo1= ' v_i : ', mask1=tmask)
CALL prt_ctl(tab2d_1=v_s (:,:,jl) , clinfo1= ' v_s : ', mask1=tmask)
CALL prt_ctl(tab2d_1=e_s (:,:,1,jl), clinfo1= ' e_snow : ', mask1=tmask)
CALL prt_ctl(tab2d_1=sv_i(:,:,jl) , clinfo1= ' sv_i : ', mask1=tmask)
CALL prt_ctl(tab2d_1=oa_i(:,:,jl) , clinfo1= ' oa_i : ', mask1=tmask)
DO jk = 1, nlay_i
CALL prt_ctl_info(' - Layer : ', ivar=jk)
CALL prt_ctl(tab2d_1=t_i(:,:,jk,jl) , clinfo1= ' t_i : ')
CALL prt_ctl(tab2d_1=e_i(:,:,jk,jl) , clinfo1= ' e_i : ')
CALL prt_ctl(tab2d_1=t_i(:,:,jk,jl) , clinfo1= ' t_i : ', mask1=tmask)
CALL prt_ctl(tab2d_1=e_i(:,:,jk,jl) , clinfo1= ' e_i : ', mask1=tmask)
END DO
END DO
CALL prt_ctl_info(' ')
CALL prt_ctl_info(' - Stresses : ')
CALL prt_ctl_info(' ~~~~~~~~~~ ')
CALL prt_ctl(tab2d_1=utau , clinfo1= ' utau : ', tab2d_2=vtau , clinfo2= ' vtau : ')
CALL prt_ctl(tab2d_1=utau_ice , clinfo1= ' utau_ice : ', tab2d_2=vtau_ice , clinfo2= ' vtau_ice : ')
CALL prt_ctl(tab2d_1=utau , clinfo1= ' utau : ', mask1 = umask, &
& tab2d_2=vtau , clinfo2= ' vtau : ', mask2 = vmask)
CALL prt_ctl(tab2d_1=utau_ice , clinfo1= ' utau_ice : ', mask1 = umask, &
& tab2d_2=vtau_ice , clinfo2= ' vtau_ice : ', mask2 = vmask)
END SUBROUTINE ice_prt3D
......
......@@ -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
......
......@@ -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
......@@ -167,18 +166,36 @@ CONTAINS
END SUBROUTINE Agrif_istate_ssh
SUBROUTINE Agrif_tra
SUBROUTINE Agrif_tra( kt, kstg )
!!----------------------------------------------------------------------
!! *** ROUTINE Agrif_tra ***
!!----------------------------------------------------------------------
INTEGER, INTENT(in) :: kt
INTEGER, OPTIONAL, INTENT(in) :: kstg
REAL(wp) :: ztindex
!
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_UseSpecialValue = l_spc_tra
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.
l_vremap = .FALSE.
......@@ -186,35 +203,52 @@ CONTAINS
END SUBROUTINE Agrif_tra
SUBROUTINE Agrif_dyn( kt )
SUBROUTINE Agrif_dyn( kt, kstg )
!!----------------------------------------------------------------------
!! *** ROUTINE Agrif_DYN ***
!!----------------------------------------------------------------------
INTEGER, INTENT(in) :: kt
INTEGER, OPTIONAL, INTENT(in) :: kstg
!
INTEGER :: ji, jj, jk ! dummy loop indices
INTEGER :: ibdy1, jbdy1, ibdy2, jbdy2
REAL(wp) :: zflag
REAL(wp), DIMENSION(jpi,jpj) :: zub, zvb
REAL(wp), DIMENSION(jpi,jpj) :: zhub, zhvb
REAL(wp) :: ztindex
!!----------------------------------------------------------------------
!
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_UseSpecialValue = ln_spc_dyn
l_vremap = ln_vert_remap
!
use_sign_north = .TRUE.
sign_north = -1.0_wp
CALL Agrif_Bc_variable( un_interp_id, procname=interpun )
CALL Agrif_Bc_variable( vn_interp_id, procname=interpvn )
CALL Agrif_Bc_variable( un_interp_id, calledweight=ztindex, procname=interpun )
CALL Agrif_Bc_variable( vn_interp_id, calledweight=ztindex, procname=interpvn )
IF( .NOT.ln_dynspg_ts ) THEN ! Get transports
ubdy(:,:) = 0._wp ; vbdy(:,:) = 0._wp
utint_stage(:,:) = 0 ; vtint_stage(:,:) = 0
CALL Agrif_Bc_variable( unb_interp_id, procname=interpunb )
CALL Agrif_Bc_variable( vnb_interp_id, procname=interpvnb )
CALL Agrif_Bc_variable( unb_interp_id, calledweight=ztindex, procname=interpunb )
CALL Agrif_Bc_variable( vnb_interp_id, calledweight=ztindex, procname=interpvnb )
ENDIF
use_sign_north = .FALSE.
......@@ -675,6 +709,13 @@ CONTAINS
!
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
!
! Enforce volume conservation if no time refinement:
......@@ -1399,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
......@@ -1431,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()
!
......@@ -1466,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
......@@ -1507,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
......@@ -1541,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()
!
......@@ -1576,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
......
......@@ -36,7 +36,7 @@ MODULE agrif_oce_sponge
# include "do_loop_substitute.h90"
!!----------------------------------------------------------------------
!! 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)
!!----------------------------------------------------------------------
CONTAINS
......@@ -50,8 +50,12 @@ CONTAINS
!!----------------------------------------------------------------------
!
#if defined SPONGE
#if defined key_RK3
zcoef = REAL(Agrif_Nbstepint(), wp)/REAL(Agrif_rhot())
#else
!! Assume persistence:
zcoef = REAL(Agrif_rhot()-1,wp)/REAL(Agrif_rhot())
#endif
Agrif_SpecialValue = 0._wp
Agrif_UseSpecialValue = l_spc_tra
......@@ -78,7 +82,12 @@ CONTAINS
!!----------------------------------------------------------------------
!
#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())
#endif
Agrif_SpecialValue = 0._wp
Agrif_UseSpecialValue = ln_spc_dyn
......
......@@ -41,7 +41,7 @@ MODULE agrif_oce_update
# include "domzgr_substitute.h90"
!!----------------------------------------------------------------------
!! 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)
!!----------------------------------------------------------------------
CONTAINS
......@@ -98,12 +98,14 @@ CONTAINS
!
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)
......@@ -544,6 +546,7 @@ CONTAINS
DO jk=1,jpkm1
DO jj=j1,j2
DO ji=i1,i2
#if ! defined key_RK3
IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler))) THEN ! Add asselin part
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) )
......@@ -553,6 +556,7 @@ CONTAINS
uu(ji,jj,jk,Kbb_a) = ( zub + rn_atfp * ( zunu - zuno) ) &
& * umask(ji,jj,jk) / e3u(ji,jj,jk,Kbb_a)
ENDIF
#endif
!
uu(ji,jj,jk,Kmm_a) = tabres_child(ji,jj,jk) * umask(ji,jj,jk)
END DO
......@@ -693,6 +697,7 @@ CONTAINS
DO jk=1,jpkm1
DO jj=j1,j2
DO ji=i1,i2
#if ! defined key_RK3
IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler))) THEN ! Add asselin part
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) )
......@@ -702,6 +707,7 @@ CONTAINS
vv(ji,jj,jk,Kbb_a) = ( zvb + rn_atfp * ( zvnu - zvno) ) &
& * vmask(ji,jj,jk) / e3v(ji,jj,jk,Kbb_a)
ENDIF
#endif
!
vv(ji,jj,jk,Kmm_a) = tabres_child(ji,jj,jk) * vmask(ji,jj,jk)
END DO
......@@ -768,12 +774,14 @@ CONTAINS
DO ji=i1,i2
!
! Update barotropic velocities:
#if ! defined key_RK3
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
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)
END IF
ENDIF
#endif
uu_b(ji,jj,Kmm_a) = tabres(ji,jj) * r1_hu(ji,jj,Kmm_a) * umask(ji,jj,1)
!
END DO
......@@ -838,12 +846,14 @@ CONTAINS
DO jj=j1,j2
DO ji=i1,i2
! Update barotropic velocities:
#if ! defined key_RK3
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
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)
END IF
ENDIF
#endif
vv_b(ji,jj,Kmm_a) = tabres(ji,jj) * r1_hv(ji,jj,Kmm_a) * vmask(ji,jj,1)
!
END DO
......@@ -903,6 +913,7 @@ CONTAINS
END DO
ELSE
!
#if ! defined key_RK3
IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler))) THEN
DO jj=j1,j2
DO ji=i1,i2
......@@ -911,6 +922,7 @@ CONTAINS
END DO
END DO
ENDIF
#endif
!
DO jj=j1,j2
DO ji=i1,i2
......@@ -977,7 +989,7 @@ CONTAINS
!
END SUBROUTINE updatevmsk
# if ! defined key_RK3
SUBROUTINE updateub2b( tabres, i1, i2, j1, j2, before )
!!----------------------------------------------------------------------
!! *** ROUTINE updateub2b ***
......@@ -1013,6 +1025,7 @@ CONTAINS
ENDIF
!
END SUBROUTINE updateub2b
# endif
SUBROUTINE reflux_sshu( tabres, i1, i2, j1, j2, before, nb, ndir )
!!---------------------------------------------
......@@ -1041,16 +1054,28 @@ 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
#endif
END DO
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
#endif
END DO
ENDIF
!
......@@ -1058,6 +1083,7 @@ CONTAINS
!
END SUBROUTINE reflux_sshu
# if ! defined key_RK3
SUBROUTINE updatevb2b( tabres, i1, i2, j1, j2, before )
!!----------------------------------------------------------------------
!! *** ROUTINE updatevb2b ***
......@@ -1093,6 +1119,7 @@ CONTAINS
ENDIF
!
END SUBROUTINE updatevb2b
# endif
SUBROUTINE reflux_sshv( tabres, i1, i2, j1, j2, before, nb, ndir )
!!---------------------------------------------
......@@ -1121,16 +1148,28 @@ 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
#endif
END DO
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
#endif
END DO
ENDIF
!
......@@ -1232,6 +1271,7 @@ CONTAINS
! of prognostic variables
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
DO jk = 1, jpkm1
DO jj=j1,j2
......@@ -1262,6 +1302,7 @@ CONTAINS
END DO
!
ENDIF
#endif
!
! 2) Updates at NOW time step:
! ----------------------------
......
......@@ -30,23 +30,42 @@ MODULE agrif_top_interp
# include "domzgr_substitute.h90"
!!----------------------------------------------------------------------
!! 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)
!!----------------------------------------------------------------------
CONTAINS
SUBROUTINE Agrif_trc
SUBROUTINE Agrif_trc( kt, kstg )
!!----------------------------------------------------------------------
!! *** ROUTINE Agrif_trc ***
!!----------------------------------------------------------------------
INTEGER, INTENT(in) :: kt
INTEGER, OPTIONAL, INTENT(in) :: kstg
!
REAL(wp) :: ztindex
!
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_UseSpecialValue = l_spc_top
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.
l_vremap = .FALSE.
......
......@@ -33,7 +33,7 @@ MODULE agrif_top_sponge
# include "domzgr_substitute.h90"
!!----------------------------------------------------------------------
!! 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)
!!----------------------------------------------------------------------
CONTAINS
......@@ -46,9 +46,12 @@ CONTAINS
!!----------------------------------------------------------------------
!
#if defined SPONGE_TOP
#if defined key_RK3
zcoef = REAL(Agrif_Nbstepint(), wp)/REAL(Agrif_rhot())
#else
!! Assume persistence:
zcoef = REAL(Agrif_rhot()-1,wp)/REAL(Agrif_rhot())
#endif
Agrif_SpecialValue = 0._wp
Agrif_UseSpecialValue = l_spc_top
l_vremap = ln_vert_remap
......
......@@ -29,7 +29,7 @@ MODULE agrif_top_update
# include "domzgr_substitute.h90"
!!----------------------------------------------------------------------
!! 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)
!!----------------------------------------------------------------------
CONTAINS
......@@ -121,7 +121,7 @@ CONTAINS
ENDIF
ENDDO
ENDDO
#if ! defined key_RK3
IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler))) THEN
! Add asselin part
DO jn = 1,jptra
......@@ -142,6 +142,7 @@ CONTAINS
END DO
END DO
ENDIF
#endif
DO jn = 1,jptra
DO jk = 1, jpkm1
DO jj = j1, j2
......@@ -160,6 +161,7 @@ CONTAINS
& * tmask(i1:i2,j1:j2,jk)
END DO
ENDDO
#if ! defined key_RK3
IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler))) THEN
! Add asselin part
DO jn = 1,jptra
......@@ -180,6 +182,7 @@ CONTAINS
END DO
END DO
ENDIF
#endif
DO jn = 1,jptra
DO jk=k1,k2
DO jj=j1,j2
......
......@@ -27,7 +27,11 @@
!!----------------------------------------------------------------------
!
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
# endif
!
! !* Agrif initialization
CALL Agrif_InitValues_cont
......@@ -410,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.
......
......@@ -617,10 +617,10 @@ CONTAINS
!!gm
IF( ln_zps .AND. .NOT. ln_c1d .AND. .NOT. ln_isfcav) &
& CALL zps_hde ( kt, Kmm, jpts, pts(:,:,:,:,Kbb), gtsu, gtsv, & ! Partial steps: before horizontal gradient
& CALL zps_hde ( kt, jpts, pts(:,:,:,:,Kbb), gtsu, gtsv, & ! Partial steps: before horizontal gradient
& rhd, gru , grv ) ! of t, s, rd at the last ocean level
IF( ln_zps .AND. .NOT. ln_c1d .AND. ln_isfcav) &
& CALL zps_hde_isf( nit000, Kmm, jpts, pts(:,:,:,:,Kbb), gtsu, gtsv, gtui, gtvi, & ! Partial steps for top cell (ISF)
& CALL zps_hde_isf( nit000, jpts, pts(:,:,:,:,Kbb), gtsu, gtsv, gtui, gtvi, & ! Partial steps for top cell (ISF)
& rhd, gru , grv , grui, grvi ) ! of t, s, rd at the last ocean level
DEALLOCATE( t_bkginc )
......
......@@ -106,12 +106,13 @@ CONTAINS
IF( ln_timing ) CALL timing_start('dia_hth')
IF( kt == nit000 ) THEN
l_hth = .FALSE.
IF( iom_use( 'mlddzt' ) .OR. iom_use( 'mldr0_3' ) .OR. iom_use( 'mldr0_1' ) .OR. &
& iom_use( 'mld_dt02' ) .OR. iom_use( 'topthdep' ) .OR. iom_use( 'mldr10_3' ) .OR. &
& iom_use( '20d' ) .OR. iom_use( '26d' ) .OR. iom_use( '28d' ) .OR. &
& iom_use( 'hc300' ) .OR. iom_use( 'hc700' ) .OR. iom_use( 'hc2000' ) .OR. &
& iom_use( 'pycndep' ) .OR. iom_use( 'tinv' ) .OR. iom_use( 'depti' ) ) l_hth = .TRUE.
!
l_hth = iom_use( 'mlddzt' ) .OR. iom_use( 'mldr0_3' ) .OR. iom_use( 'mldr0_1' ) .OR. &
& iom_use( 'mld_dt02' ) .OR. iom_use( 'topthdep' ) .OR. iom_use( 'mldr10_3' ) .OR. &
& iom_use( '20d' ) .OR. iom_use( '26d' ) .OR. iom_use( '28d' ) .OR. &
& iom_use( 'hc300' ) .OR. iom_use( 'hc700' ) .OR. iom_use( 'hc2000' ) .OR. &
& iom_use( 'pycndep' ) .OR. iom_use( 'tinv' ) .OR. iom_use( 'depti' )
!
! ! allocate dia_hth array
IF( l_hth ) THEN
IF( dia_hth_alloc() /= 0 ) CALL ctl_stop( 'STOP', 'dia_hth : unable to allocate standard arrays' )
......@@ -124,11 +125,12 @@ CONTAINS
IF( l_hth ) THEN
!
IF( iom_use( 'mlddzt' ) .OR. iom_use( 'mldr0_3' ) .OR. iom_use( 'mldr0_1' ) ) THEN
! initialization
ztinv (:,:) = 0._wp
zdepinv(:,:) = 0._wp
zmaxdzT(:,:) = 0._wp
! initialization
IF( iom_use( 'tinv' ) ) ztinv (:,:) = 0._wp
IF( iom_use( 'depti' ) ) zdepinv(:,:) = 0._wp
IF( iom_use( 'mlddzt' ) ) zmaxdzT(:,:) = 0._wp
IF( iom_use( 'mlddzt' ) .OR. iom_use( 'mld_dt02' ) .OR. iom_use( 'topthdep' ) &
& .OR. iom_use( 'mldr10_3' ) .OR. iom_use( 'pycndep' ) ) THEN
DO_2D( 1, 1, 1, 1 )
zztmp = gdepw(ji,jj,mbkt(ji,jj)+1,Kmm)
hth (ji,jj) = zztmp
......@@ -137,6 +139,8 @@ CONTAINS
zrho10_3(ji,jj) = zztmp
zpycn (ji,jj) = zztmp
END_2D
ENDIF
IF( iom_use( 'mldr0_3' ) .OR. iom_use( 'mldr0_1' ) ) THEN
IF( nla10 > 1 ) THEN
DO_2D( 1, 1, 1, 1 )
zztmp = gdepw(ji,jj,mbkt(ji,jj)+1,Kmm)
......@@ -144,25 +148,9 @@ CONTAINS
zrho0_1(ji,jj) = zztmp
END_2D
ENDIF
! Preliminary computation
! computation of zdelr = (dr/dT)(T,S,10m)*(-0.2 degC)
DO_2D( 1, 1, 1, 1 )
IF( tmask(ji,jj,nla10) == 1. ) THEN
zu = 1779.50 + 11.250 * ts(ji,jj,nla10,jp_tem,Kmm) - 3.80 * ts(ji,jj,nla10,jp_sal,Kmm) &
& - 0.0745 * ts(ji,jj,nla10,jp_tem,Kmm) * ts(ji,jj,nla10,jp_tem,Kmm) &
& - 0.0100 * ts(ji,jj,nla10,jp_tem,Kmm) * ts(ji,jj,nla10,jp_sal,Kmm)
zv = 5891.00 + 38.000 * ts(ji,jj,nla10,jp_tem,Kmm) + 3.00 * ts(ji,jj,nla10,jp_sal,Kmm) &
& - 0.3750 * ts(ji,jj,nla10,jp_tem,Kmm) * ts(ji,jj,nla10,jp_tem,Kmm)
zut = 11.25 - 0.149 * ts(ji,jj,nla10,jp_tem,Kmm) - 0.01 * ts(ji,jj,nla10,jp_sal,Kmm)
zvt = 38.00 - 0.750 * ts(ji,jj,nla10,jp_tem,Kmm)
zw = (zu + 0.698*zv) * (zu + 0.698*zv)
zdelr(ji,jj) = ztem2 * (1000.*(zut*zv - zvt*zu)/zw)
ELSE
zdelr(ji,jj) = 0._wp
ENDIF
END_2D
ENDIF
IF( iom_use( 'mlddzt' ) .OR. iom_use( 'mldr0_3' ) .OR. iom_use( 'mldr0_1' ) ) THEN
! ------------------------------------------------------------- !
! thermocline depth: strongest vertical gradient of temperature !
! turbocline depth (mixing layer depth): avt = zavt5 !
......@@ -198,6 +186,25 @@ CONTAINS
!
IF( iom_use( 'mld_dt02' ) .OR. iom_use( 'topthdep' ) .OR. iom_use( 'mldr10_3' ) .OR. &
& iom_use( 'pycndep' ) .OR. iom_use( 'tinv' ) .OR. iom_use( 'depti' ) ) THEN
!
! Preliminary computation
! computation of zdelr = (dr/dT)(T,S,10m)*(-0.2 degC)
DO_2D( 1, 1, 1, 1 )
IF( tmask(ji,jj,nla10) == 1. ) THEN
zu = 1779.50 + 11.250 * ts(ji,jj,nla10,jp_tem,Kmm) - 3.80 * ts(ji,jj,nla10,jp_sal,Kmm) &
& - 0.0745 * ts(ji,jj,nla10,jp_tem,Kmm) * ts(ji,jj,nla10,jp_tem,Kmm) &
& - 0.0100 * ts(ji,jj,nla10,jp_tem,Kmm) * ts(ji,jj,nla10,jp_sal,Kmm)
zv = 5891.00 + 38.000 * ts(ji,jj,nla10,jp_tem,Kmm) + 3.00 * ts(ji,jj,nla10,jp_sal,Kmm) &
& - 0.3750 * ts(ji,jj,nla10,jp_tem,Kmm) * ts(ji,jj,nla10,jp_tem,Kmm)
zut = 11.25 - 0.149 * ts(ji,jj,nla10,jp_tem,Kmm) - 0.01 * ts(ji,jj,nla10,jp_sal,Kmm)
zvt = 38.00 - 0.750 * ts(ji,jj,nla10,jp_tem,Kmm)
zw = (zu + 0.698*zv) * (zu + 0.698*zv)
zdelr(ji,jj) = ztem2 * (1000.*(zut*zv - zvt*zu)/zw)
ELSE
zdelr(ji,jj) = 0._wp
ENDIF
END_2D
!
! ------------------------------------------------------------- !
! MLD: abs( tn - tn(10m) ) = ztem2 !
! Top of thermocline: tn = tn(10m) - ztem2 !
......