From 4389904baaceaf51a0d1465868ebf66231b404fc Mon Sep 17 00:00:00 2001
From: sibylle <sibylle.techene@locean.ipsl.fr>
Date: Tue, 12 Apr 2022 17:20:29 +0200
Subject: [PATCH] #2 RK3 branch merged with main + update stprk3 wrt zps_hde
change of interface
---
cfgs/AGRIF_DEMO/EXPREF/1_context_nemo.xml | 1 -
cfgs/AGRIF_DEMO/EXPREF/1_namelist_cfg | 6 +-
cfgs/AGRIF_DEMO/EXPREF/2_namelist_cfg | 4 +
cfgs/AGRIF_DEMO/EXPREF/3_namelist_cfg | 4 +
cfgs/AGRIF_DEMO/EXPREF/namelist_cfg | 4 +
cfgs/ORCA2_ICE_PISCES/EXPREF/namelist_cfg | 4 +
cfgs/SHARED/namelist_ref | 4 +-
cfgs/WED025/EXPREF/namelist_cfg | 4 +
doc/namelists/namzdf_tke | 2 +-
src/ICE/icectl.F90 | 61 +++++-----
src/OCE/ASM/asminc.F90 | 4 +-
src/OCE/DIA/diahth.F90 | 67 ++++++-----
src/OCE/DOM/domzgr.F90 | 4 +-
src/OCE/DOM/dtatsd.F90 | 57 ++++-----
src/OCE/DYN/dynhpg.F90 | 2 +-
src/OCE/DYN/dynvor.F90 | 2 +-
src/OCE/IOM/iom.F90 | 2 +-
src/OCE/IOM/iom_nf90.F90 | 6 +-
src/OCE/LBC/lbc_lnk_pt2pt_generic.h90 | 2 +-
src/OCE/LBC/mppini.F90 | 32 +++--
src/OCE/SBC/geo2ocean.F90 | 82 ++++++-------
src/OCE/SBC/sbcblk.F90 | 54 ++++++---
src/OCE/SBC/sbcfwb.F90 | 10 +-
src/OCE/SBC/sbcmod.F90 | 11 +-
src/OCE/SBC/sbcrnf.F90 | 6 +-
src/OCE/SBC/sbcssm.F90 | 2 +-
src/OCE/TRA/eosbn2.F90 | 4 +-
src/OCE/TRA/traldf_lap_blp.F90 | 4 +-
src/OCE/TRA/traldf_triad.F90 | 1 -
src/OCE/TRA/zpshde.F90 | 138 +++++++++++-----------
src/OCE/ZDF/zdfgls.F90 | 7 +-
src/OCE/ZDF/zdftke.F90 | 5 +
src/OCE/step.F90 | 4 +-
src/OCE/stpmlf.F90 | 4 +-
src/OCE/stprk3.F90 | 16 +--
src/OFF/dtadyn.F90 | 4 +-
src/TOP/TRP/trctrp.F90 | 4 +-
src/TOP/trcini.F90 | 1 -
38 files changed, 351 insertions(+), 278 deletions(-)
diff --git a/cfgs/AGRIF_DEMO/EXPREF/1_context_nemo.xml b/cfgs/AGRIF_DEMO/EXPREF/1_context_nemo.xml
index 87b0d3b2c..ee852b604 100644
--- a/cfgs/AGRIF_DEMO/EXPREF/1_context_nemo.xml
+++ b/cfgs/AGRIF_DEMO/EXPREF/1_context_nemo.xml
@@ -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 -->
diff --git a/cfgs/AGRIF_DEMO/EXPREF/1_namelist_cfg b/cfgs/AGRIF_DEMO/EXPREF/1_namelist_cfg
index fb7cc8224..6dacc95db 100644
--- a/cfgs/AGRIF_DEMO/EXPREF/1_namelist_cfg
+++ b/cfgs/AGRIF_DEMO/EXPREF/1_namelist_cfg
@@ -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 *** !!
diff --git a/cfgs/AGRIF_DEMO/EXPREF/2_namelist_cfg b/cfgs/AGRIF_DEMO/EXPREF/2_namelist_cfg
index ed652bd07..cb10389e9 100644
--- a/cfgs/AGRIF_DEMO/EXPREF/2_namelist_cfg
+++ b/cfgs/AGRIF_DEMO/EXPREF/2_namelist_cfg
@@ -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 *** !!
diff --git a/cfgs/AGRIF_DEMO/EXPREF/3_namelist_cfg b/cfgs/AGRIF_DEMO/EXPREF/3_namelist_cfg
index 4b74ffed5..5a082933d 100644
--- a/cfgs/AGRIF_DEMO/EXPREF/3_namelist_cfg
+++ b/cfgs/AGRIF_DEMO/EXPREF/3_namelist_cfg
@@ -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 *** !!
diff --git a/cfgs/AGRIF_DEMO/EXPREF/namelist_cfg b/cfgs/AGRIF_DEMO/EXPREF/namelist_cfg
index 63973baf2..94bb5c091 100644
--- a/cfgs/AGRIF_DEMO/EXPREF/namelist_cfg
+++ b/cfgs/AGRIF_DEMO/EXPREF/namelist_cfg
@@ -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 *** !!
diff --git a/cfgs/ORCA2_ICE_PISCES/EXPREF/namelist_cfg b/cfgs/ORCA2_ICE_PISCES/EXPREF/namelist_cfg
index 7f302e9a0..c8353d3ad 100644
--- a/cfgs/ORCA2_ICE_PISCES/EXPREF/namelist_cfg
+++ b/cfgs/ORCA2_ICE_PISCES/EXPREF/namelist_cfg
@@ -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)
diff --git a/cfgs/SHARED/namelist_ref b/cfgs/SHARED/namelist_ref
index cde6a6b16..ae9d42ee7 100644
--- a/cfgs/SHARED/namelist_ref
+++ b/cfgs/SHARED/namelist_ref
@@ -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
diff --git a/cfgs/WED025/EXPREF/namelist_cfg b/cfgs/WED025/EXPREF/namelist_cfg
index b4293c931..71602bf7b 100644
--- a/cfgs/WED025/EXPREF/namelist_cfg
+++ b/cfgs/WED025/EXPREF/namelist_cfg
@@ -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)
diff --git a/doc/namelists/namzdf_tke b/doc/namelists/namzdf_tke
index 85990fe70..d5d10df2c 100644
--- a/doc/namelists/namzdf_tke
+++ b/doc/namelists/namzdf_tke
@@ -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 )
diff --git a/src/ICE/icectl.F90 b/src/ICE/icectl.F90
index e66fe7798..ad20e6734 100644
--- a/src/ICE/icectl.F90
+++ b/src/ICE/icectl.F90
@@ -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
diff --git a/src/OCE/ASM/asminc.F90 b/src/OCE/ASM/asminc.F90
index aba92a15d..54371a62c 100644
--- a/src/OCE/ASM/asminc.F90
+++ b/src/OCE/ASM/asminc.F90
@@ -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 )
diff --git a/src/OCE/DIA/diahth.F90 b/src/OCE/DIA/diahth.F90
index 2cee430d3..0b8650acf 100644
--- a/src/OCE/DIA/diahth.F90
+++ b/src/OCE/DIA/diahth.F90
@@ -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 !
diff --git a/src/OCE/DOM/domzgr.F90 b/src/OCE/DOM/domzgr.F90
index 76fc34641..85f89ed2c 100644
--- a/src/OCE/DOM/domzgr.F90
+++ b/src/OCE/DOM/domzgr.F90
@@ -324,9 +324,9 @@ CONTAINS
#if defined key_qco && key_isf
DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 2, jpk ) ! vertical sum at partial cell xxxx other level
IF( jk == k_top(ji,jj) ) THEN ! first ocean point : partial cell
- gdept_0(ji,jj,jk) = gdepw_0(ji,jj,jk ) + 0.5_wp * e3w_0(ji,jj,jk) ! = risfdep + 1/2 e3w_0(mikt)
+ pdept(ji,jj,jk) = pdepw(ji,jj,jk ) + 0.5_wp * pe3w(ji,jj,jk) ! = risfdep + 1/2 e3w_0(mikt)
ELSE ! other levels
- gdept_0(ji,jj,jk) = gdept_0(ji,jj,jk-1) + e3w_0(ji,jj,jk)
+ pdept(ji,jj,jk) = pdept(ji,jj,jk-1) + pe3w(ji,jj,jk)
ENDIF
END_3D
#endif
diff --git a/src/OCE/DOM/dtatsd.F90 b/src/OCE/DOM/dtatsd.F90
index 7ffccc4cf..5863789cc 100644
--- a/src/OCE/DOM/dtatsd.F90
+++ b/src/OCE/DOM/dtatsd.F90
@@ -199,8 +199,7 @@ CONTAINS
ptsd(ji,jj,jk,jp_sal) = sf_tsd(jp_sal)%fnow(ji,jj,jk)
END_3D
!
-! JC I think it's more convenient to consider the general sco case as the rule
-! IF( ln_sco ) THEN !== s- or mixed s-zps-coordinate ==!
+ IF( ln_sco ) THEN !== s- or mixed s-zps-coordinate ==!
!
IF( .NOT. l_istiled .OR. ntile == 1 ) THEN ! Do only on the first tile
IF( kt == nit000 .AND. lwp )THEN
@@ -236,31 +235,35 @@ CONTAINS
ptsd(ji,jj,jpk,jp_sal) = 0._wp
END_2D
!
-! ELSE !== z- or zps- coordinate ==!
-! !
-! DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpk )
-! ptsd(ji,jj,jk,jp_tem) = ptsd(ji,jj,jk,jp_tem) * tmask(ji,jj,jk) ! Mask
-! ptsd(ji,jj,jk,jp_sal) = ptsd(ji,jj,jk,jp_sal) * tmask(ji,jj,jk)
-! END_3D
-! !
-! IF( ln_zps ) THEN ! zps-coordinate (partial steps) interpolation at the last ocean level
-! DO_2D( nn_hls, nn_hls, nn_hls, nn_hls )
-! ik = mbkt(ji,jj)
-! IF( ik > 1 ) THEN
-! zl = ( gdept_1d(ik) - gdept_0(ji,jj,ik) ) / ( gdept_1d(ik) - gdept_1d(ik-1) )
-! ptsd(ji,jj,ik,jp_tem) = (1.-zl) * ptsd(ji,jj,ik,jp_tem) + zl * ptsd(ji,jj,ik-1,jp_tem)
-! ptsd(ji,jj,ik,jp_sal) = (1.-zl) * ptsd(ji,jj,ik,jp_sal) + zl * ptsd(ji,jj,ik-1,jp_sal)
-! ENDIF
-! ik = mikt(ji,jj)
-! IF( ik > 1 ) THEN
-! zl = ( gdept_0(ji,jj,ik) - gdept_1d(ik) ) / ( gdept_1d(ik+1) - gdept_1d(ik) )
-! ptsd(ji,jj,ik,jp_tem) = (1.-zl) * ptsd(ji,jj,ik,jp_tem) + zl * ptsd(ji,jj,ik+1,jp_tem)
-! ptsd(ji,jj,ik,jp_sal) = (1.-zl) * ptsd(ji,jj,ik,jp_sal) + zl * ptsd(ji,jj,ik+1,jp_sal)
-! END IF
-! END_2D
-! ENDIF
-! !
-! ENDIF
+ ELSE !== z- or zps- coordinate ==!
+ !
+ ! We must keep this definition in a case different from the general case of s-coordinate as we don't
+ ! want to use "underground" values (levels below ocean bottom) to be able to start the model from
+ ! masked temp and sal (read for example in a restart or in output.init)
+ !
+ DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpk )
+ ptsd(ji,jj,jk,jp_tem) = ptsd(ji,jj,jk,jp_tem) * tmask(ji,jj,jk) ! Mask
+ ptsd(ji,jj,jk,jp_sal) = ptsd(ji,jj,jk,jp_sal) * tmask(ji,jj,jk)
+ END_3D
+ !
+ IF( ln_zps ) THEN ! zps-coordinate (partial steps) interpolation at the last ocean level
+ DO_2D( nn_hls, nn_hls, nn_hls, nn_hls )
+ ik = mbkt(ji,jj)
+ IF( ik > 1 ) THEN
+ zl = ( gdept_1d(ik) - gdept_0(ji,jj,ik) ) / ( gdept_1d(ik) - gdept_1d(ik-1) )
+ ptsd(ji,jj,ik,jp_tem) = (1.-zl) * ptsd(ji,jj,ik,jp_tem) + zl * ptsd(ji,jj,ik-1,jp_tem)
+ ptsd(ji,jj,ik,jp_sal) = (1.-zl) * ptsd(ji,jj,ik,jp_sal) + zl * ptsd(ji,jj,ik-1,jp_sal)
+ ENDIF
+ ik = mikt(ji,jj)
+ IF( ik > 1 ) THEN
+ zl = ( gdept_0(ji,jj,ik) - gdept_1d(ik) ) / ( gdept_1d(ik+1) - gdept_1d(ik) )
+ ptsd(ji,jj,ik,jp_tem) = (1.-zl) * ptsd(ji,jj,ik,jp_tem) + zl * ptsd(ji,jj,ik+1,jp_tem)
+ ptsd(ji,jj,ik,jp_sal) = (1.-zl) * ptsd(ji,jj,ik,jp_sal) + zl * ptsd(ji,jj,ik+1,jp_sal)
+ END IF
+ END_2D
+ ENDIF
+ !
+ ENDIF
!
IF( .NOT.ln_tsd_dmp ) THEN !== deallocate T & S structure ==!
! (data used only for initialisation)
diff --git a/src/OCE/DYN/dynhpg.F90 b/src/OCE/DYN/dynhpg.F90
index 4ca2b7a46..145c7f9e1 100644
--- a/src/OCE/DYN/dynhpg.F90
+++ b/src/OCE/DYN/dynhpg.F90
@@ -328,7 +328,7 @@ CONTAINS
ENDIF
! Partial steps: Compute NOW horizontal gradient of t, s, rd at the last ocean level
- CALL zps_hde( kt, Kmm, jpts, ts(:,:,:,:,Kmm), zgtsu, zgtsv, rhd, zgru , zgrv )
+ CALL zps_hde( kt, jpts, ts(:,:,:,:,Kmm), zgtsu, zgtsv, rhd, zgru , zgrv )
! Local constant initialization
zcoef0 = - grav * 0.5_wp
diff --git a/src/OCE/DYN/dynvor.F90 b/src/OCE/DYN/dynvor.F90
index 0dde705b2..03dda78ce 100644
--- a/src/OCE/DYN/dynvor.F90
+++ b/src/OCE/DYN/dynvor.F90
@@ -1098,7 +1098,7 @@ CONTAINS
!
CASE DEFAULT !* F-point metric term : pre-compute di(e2u)/(2*e1e2f) and dj(e1v)/(2*e1e2f)
ALLOCATE( di_e2v_2e1e2f(jpi,jpj), dj_e1u_2e1e2f(jpi,jpj) )
- DO_2D( 1, 0, 1, 0 )
+ DO_2D( 0, 0, 0, 0 )
di_e2v_2e1e2f(ji,jj) = ( e2v(ji+1,jj ) - e2v(ji,jj) ) * 0.5 * r1_e1e2f(ji,jj)
dj_e1u_2e1e2f(ji,jj) = ( e1u(ji ,jj+1) - e1u(ji,jj) ) * 0.5 * r1_e1e2f(ji,jj)
END_2D
diff --git a/src/OCE/IOM/iom.F90 b/src/OCE/IOM/iom.F90
index c137d5737..d66b13d19 100644
--- a/src/OCE/IOM/iom.F90
+++ b/src/OCE/IOM/iom.F90
@@ -752,7 +752,7 @@ CONTAINS
! =============
clname = trim(cdname)
IF ( .NOT. Agrif_Root() .AND. .NOT. lliof ) THEN
- iln = INDEX(clname,'/')
+ iln = INDEX(clname,'/', BACK=.TRUE.)
cltmpn = clname(1:iln)
clname = clname(iln+1:LEN_TRIM(clname))
clname=TRIM(cltmpn)//TRIM(Agrif_CFixed())//'_'//TRIM(clname)
diff --git a/src/OCE/IOM/iom_nf90.F90 b/src/OCE/IOM/iom_nf90.F90
index 5f664536b..c49d9538e 100644
--- a/src/OCE/IOM/iom_nf90.F90
+++ b/src/OCE/IOM/iom_nf90.F90
@@ -556,6 +556,7 @@ CONTAINS
LOGICAL :: lchunk ! logical switch to activate chunking and compression
! ! when appropriate (currently chunking is applied to 4d fields only)
INTEGER :: idlv ! local variable
+ CHARACTER(LEN=256) :: ccname ! local variable
!---------------------------------------------------------------------
!
clinfo = ' iom_nf90_rp0123d, file: '//TRIM(iom_file(kiomid)%name)//', var: '//TRIM(cdvar)
@@ -571,7 +572,10 @@ CONTAINS
! define the dimension variables if it is not already done
DO jd = 1, 2
CALL iom_nf90_check(NF90_INQUIRE_DIMENSION(if90id,jd,iom_file(kiomid)%cn_var(jd),iom_file(kiomid)%dimsz(jd,jd)),clinfo)
- CALL iom_nf90_check(NF90_DEF_VAR( if90id, TRIM(iom_file(kiomid)%cn_var(jd)), NF90_FLOAT , (/ 1, 2 /), &
+ ccname = TRIM(iom_file(kiomid)%cn_var(jd))
+ IF ( ccname == 'x') ccname = 'nav_lon'
+ IF ( ccname == 'y') ccname = 'nav_lat'
+ CALL iom_nf90_check(NF90_DEF_VAR( if90id, ccname, NF90_FLOAT , (/ 1, 2 /), &
& iom_file(kiomid)%nvid(jd) ), clinfo)
END DO
iom_file(kiomid)%dimsz(2,1) = iom_file(kiomid)%dimsz(2,2) ! second dim of first variable
diff --git a/src/OCE/LBC/lbc_lnk_pt2pt_generic.h90 b/src/OCE/LBC/lbc_lnk_pt2pt_generic.h90
index 8d336f7d2..395a3e93c 100644
--- a/src/OCE/LBC/lbc_lnk_pt2pt_generic.h90
+++ b/src/OCE/LBC/lbc_lnk_pt2pt_generic.h90
@@ -120,7 +120,7 @@
! ! ip0i ip1i im1i im0i
!
iwewe(:) = (/ jpwe,jpea,jpwe,jpea /) ; issnn(:) = (/ jpso,jpso,jpno,jpno /)
- !cd sides: west east south north ; corners: so-we, so-ea, no-we, no-ea
+ ! sides: west east south north ; corners: so-we, so-ea, no-we, no-ea
isizei(1:4) = (/ ihls, ihls, ipi, ipi /) ; isizei(5:8) = ihls ! i- count
isizej(1:4) = (/ Nj_0, Nj_0, ihls, ihls /) ; isizej(5:8) = ihls ! j- count
ishtSi(1:4) = (/ ip1i, im1i, ip0i, ip0i /) ; ishtSi(5:8) = ishtSi( iwewe ) ! i- shift send data
diff --git a/src/OCE/LBC/mppini.F90 b/src/OCE/LBC/mppini.F90
index 31efed848..e5cd21c35 100644
--- a/src/OCE/LBC/mppini.F90
+++ b/src/OCE/LBC/mppini.F90
@@ -175,6 +175,7 @@ CONTAINS
ENDIF
WRITE(numout,*) ' avoid use of mpi_allgather at the north fold ln_nnogather = ', ln_nnogather
WRITE(numout,*) ' halo width (applies to both rows and columns) nn_hls = ', nn_hls
+ WRITE(numout,*) ' choice of communication method nn_comm = ', nn_comm
ENDIF
!
IF(lwm) WRITE( numond, nammpp )
@@ -854,6 +855,7 @@ CONTAINS
IF(lwp) THEN
WRITE(numout,*)
WRITE(numout,*) ' -----------------------------------------------------------'
+ CALL FLUSH(numout)
ENDIF
CALL mppsync
CALL mppstop( ld_abort = .TRUE. )
@@ -1152,11 +1154,11 @@ CONTAINS
!!
!!----------------------------------------------------------------------
INTEGER :: inumsave
- INTEGER :: jh
+ INTEGER :: ji,jj,jh
INTEGER :: ipi, ipj
INTEGER :: iiwe, iiea, iist, iisz
INTEGER :: ijso, ijno, ijst, ijsz
- REAL(wp), DIMENSION(:,:), ALLOCATABLE :: zmsk
+ REAL(wp), DIMENSION(:,:), ALLOCATABLE :: zmsk0, zmsk
LOGICAL , DIMENSION(Ni_0,Nj_0,1) :: lloce
!!----------------------------------------------------------------------
!
@@ -1175,13 +1177,21 @@ CONTAINS
ipi = Ni_0 + 2*jh ! local domain size
ipj = Nj_0 + 2*jh
!
- ALLOCATE( zmsk(ipi,ipj) )
- zmsk(jh+1:jh+Ni_0,jh+1:jh+Nj_0) = REAL(COUNT(lloce, dim = 3), wp) ! define inner domain -> need REAL to use lbclnk
- CALL lbc_lnk('mppini', zmsk, 'T', 1._wp, khls = jh) ! fill halos
- ! Beware, coastal F points can be used in the code -> we may need communications for these points F points even if tmask = 0
- ! -> the mask we must use here is equal to 1 as soon as one of the 4 neighbours is oce (sum of the mask, not multiplication)
- zmsk(jh+1:jh+Ni_0,jh+1:jh+Nj_0) = zmsk(jh+1:jh+Ni_0,jh+1 :jh+Nj_0 ) + zmsk(jh+1+1:jh+Ni_0+1,jh+1 :jh+Nj_0 ) &
- & + zmsk(jh+1:jh+Ni_0,jh+1+1:jh+Nj_0+1) + zmsk(jh+1+1:jh+Ni_0+1,jh+1+1:jh+Nj_0+1)
+ ALLOCATE( zmsk0(ipi,ipj), zmsk(ipi,ipj) )
+ zmsk0(jh+1:jh+Ni_0,jh+1:jh+Nj_0) = REAL(COUNT(lloce, dim = 3), wp) ! define inner domain -> need REAL to use lbclnk
+ CALL lbc_lnk('mppini', zmsk0, 'T', 1._wp, khls = jh) ! fill halos
+ ! Beware about the mask we must use here :
+ DO jj = jh+1, jh+Nj_0
+ DO ji = jh+1, jh+Ni_0
+ zmsk(ji,jj) = zmsk0(ji,jj) &
+ ! 1) dynvor may use scale factors on i+1 (e2v for di_e2v_2e1e2f) and j+1 (e1u for dj_e1u_2e1e2f) even if land
+ ! -> the mask must be > 1 if south/west neighbours is oce as we may need to send these arrays to these neighbours
+ & + zmsk0(ji-1,jj) + zmsk0(ji,jj-1) &
+ ! 2) coastal F points can be used, so we may need communications for these points F points even IF tmask = 0
+ ! -> the mask must be > 1 as soon as one of the 3 neighbours is oce: (i,j+1) (i+1,j) (i+1,j+1)
+ & + zmsk0(ji+1,jj) + zmsk0(ji,jj+1) + zmsk0(ji+1,jj+1)
+ END DO
+ END DO
CALL lbc_lnk('mppini', zmsk, 'T', 1._wp, khls = jh) ! fill halos again!
!
iiwe = jh ; iiea = Ni_0 ! bottom-left corner - 1 of the sent data
@@ -1206,7 +1216,7 @@ IF( nn_comm == 1 ) THEN ! SM: NOT WORKING FOR NEIGHBOURHOOD COLLECTIVE COM
!
iiwe = iiwe-jh ; iiea = iiea+jh ! bottom-left corner - 1 of the received data
ijso = ijso-jh ; ijno = ijno+jh
- ! do not send if we send only land points
+ ! do not recv if we recv only land points
IF( NINT(SUM( zmsk(iiwe+1:iiwe+jh ,ijst+1:ijst+ijsz) )) == 0 ) mpiRnei(jh,jpwe) = -1
IF( NINT(SUM( zmsk(iiea+1:iiea+jh ,ijst+1:ijst+ijsz) )) == 0 ) mpiRnei(jh,jpea) = -1
IF( NINT(SUM( zmsk(iist+1:iist+iisz,ijso+1:ijso+jh ) )) == 0 ) mpiRnei(jh,jpso) = -1
@@ -1225,7 +1235,7 @@ ENDIF
IF( mpiRnei(jh,jpno) > -1 ) mpiRnei(jh, (/jpnw,jpne/) ) = -1 ! NW and NE corners will be received through North nei
ENDIF
!
- DEALLOCATE( zmsk )
+ DEALLOCATE( zmsk0, zmsk )
!
CALL mpp_ini_nc(jh) ! Initialize/Update communicator for neighbourhood collective communications
!
diff --git a/src/OCE/SBC/geo2ocean.F90 b/src/OCE/SBC/geo2ocean.F90
index 21e80952c..182e5e030 100644
--- a/src/OCE/SBC/geo2ocean.F90
+++ b/src/OCE/SBC/geo2ocean.F90
@@ -163,71 +163,71 @@ CONTAINS
!
zlam = plamt(ji,jj) ! north pole direction & modulous (at t-point)
zphi = pphit(ji,jj)
- zxnpt = 0. - 2. * COS( rad*zlam ) * TAN( rpi/4. - rad*zphi/2. )
- zynpt = 0. - 2. * SIN( rad*zlam ) * TAN( rpi/4. - rad*zphi/2. )
+ zxnpt = 0._wp- 2._wp* COS( rad*zlam ) * TAN( rpi/4._wp- rad*zphi/2._wp)
+ zynpt = 0._wp- 2._wp* SIN( rad*zlam ) * TAN( rpi/4._wp- rad*zphi/2._wp)
znnpt = zxnpt*zxnpt + zynpt*zynpt
!
zlam = plamu(ji,jj) ! north pole direction & modulous (at u-point)
zphi = pphiu(ji,jj)
- zxnpu = 0. - 2. * COS( rad*zlam ) * TAN( rpi/4. - rad*zphi/2. )
- zynpu = 0. - 2. * SIN( rad*zlam ) * TAN( rpi/4. - rad*zphi/2. )
+ zxnpu = 0._wp- 2._wp* COS( rad*zlam ) * TAN( rpi/4._wp- rad*zphi/2._wp)
+ zynpu = 0._wp- 2._wp* SIN( rad*zlam ) * TAN( rpi/4._wp- rad*zphi/2._wp)
znnpu = zxnpu*zxnpu + zynpu*zynpu
!
zlam = plamv(ji,jj) ! north pole direction & modulous (at v-point)
zphi = pphiv(ji,jj)
- zxnpv = 0. - 2. * COS( rad*zlam ) * TAN( rpi/4. - rad*zphi/2. )
- zynpv = 0. - 2. * SIN( rad*zlam ) * TAN( rpi/4. - rad*zphi/2. )
+ zxnpv = 0._wp- 2._wp* COS( rad*zlam ) * TAN( rpi/4._wp- rad*zphi/2._wp)
+ zynpv = 0._wp- 2._wp* SIN( rad*zlam ) * TAN( rpi/4._wp- rad*zphi/2._wp)
znnpv = zxnpv*zxnpv + zynpv*zynpv
!
zlam = plamf(ji,jj) ! north pole direction & modulous (at f-point)
zphi = pphif(ji,jj)
- zxnpf = 0. - 2. * COS( rad*zlam ) * TAN( rpi/4. - rad*zphi/2. )
- zynpf = 0. - 2. * SIN( rad*zlam ) * TAN( rpi/4. - rad*zphi/2. )
+ zxnpf = 0._wp- 2._wp* COS( rad*zlam ) * TAN( rpi/4._wp- rad*zphi/2._wp)
+ zynpf = 0._wp- 2._wp* SIN( rad*zlam ) * TAN( rpi/4._wp- rad*zphi/2._wp)
znnpf = zxnpf*zxnpf + zynpf*zynpf
!
zlam = plamv(ji,jj ) ! j-direction: v-point segment direction (around t-point)
zphi = pphiv(ji,jj )
zlan = plamv(ji,jj-1)
zphh = pphiv(ji,jj-1)
- zxvvt = 2. * COS( rad*zlam ) * TAN( rpi/4. - rad*zphi/2. ) &
- & - 2. * COS( rad*zlan ) * TAN( rpi/4. - rad*zphh/2. )
- zyvvt = 2. * SIN( rad*zlam ) * TAN( rpi/4. - rad*zphi/2. ) &
- & - 2. * SIN( rad*zlan ) * TAN( rpi/4. - rad*zphh/2. )
+ zxvvt = 2._wp* COS( rad*zlam ) * TAN( rpi/4._wp- rad*zphi/2._wp) &
+ & - 2._wp* COS( rad*zlan ) * TAN( rpi/4._wp- rad*zphh/2._wp)
+ zyvvt = 2._wp* SIN( rad*zlam ) * TAN( rpi/4._wp- rad*zphi/2._wp) &
+ & - 2._wp* SIN( rad*zlan ) * TAN( rpi/4._wp- rad*zphh/2._wp)
znvvt = SQRT( znnpt * ( zxvvt*zxvvt + zyvvt*zyvvt ) )
- znvvt = MAX( znvvt, 1.e-14 )
+ znvvt = MAX( znvvt, 1.e-14_wp )
!
zlam = plamf(ji,jj ) ! j-direction: f-point segment direction (around u-point)
zphi = pphif(ji,jj )
zlan = plamf(ji,jj-1)
zphh = pphif(ji,jj-1)
- zxffu = 2. * COS( rad*zlam ) * TAN( rpi/4. - rad*zphi/2. ) &
- & - 2. * COS( rad*zlan ) * TAN( rpi/4. - rad*zphh/2. )
- zyffu = 2. * SIN( rad*zlam ) * TAN( rpi/4. - rad*zphi/2. ) &
- & - 2. * SIN( rad*zlan ) * TAN( rpi/4. - rad*zphh/2. )
+ zxffu = 2._wp* COS( rad*zlam ) * TAN( rpi/4._wp- rad*zphi/2._wp) &
+ & - 2._wp* COS( rad*zlan ) * TAN( rpi/4._wp- rad*zphh/2._wp)
+ zyffu = 2._wp* SIN( rad*zlam ) * TAN( rpi/4._wp- rad*zphi/2._wp) &
+ & - 2._wp* SIN( rad*zlan ) * TAN( rpi/4._wp- rad*zphh/2._wp)
znffu = SQRT( znnpu * ( zxffu*zxffu + zyffu*zyffu ) )
- znffu = MAX( znffu, 1.e-14 )
+ znffu = MAX( znffu, 1.e-14_wp )
!
zlam = plamf(ji ,jj) ! i-direction: f-point segment direction (around v-point)
zphi = pphif(ji ,jj)
zlan = plamf(ji-1,jj)
zphh = pphif(ji-1,jj)
- zxffv = 2. * COS( rad*zlam ) * TAN( rpi/4. - rad*zphi/2. ) &
- & - 2. * COS( rad*zlan ) * TAN( rpi/4. - rad*zphh/2. )
- zyffv = 2. * SIN( rad*zlam ) * TAN( rpi/4. - rad*zphi/2. ) &
- & - 2. * SIN( rad*zlan ) * TAN( rpi/4. - rad*zphh/2. )
+ zxffv = 2._wp* COS( rad*zlam ) * TAN( rpi/4._wp- rad*zphi/2._wp) &
+ & - 2._wp* COS( rad*zlan ) * TAN( rpi/4._wp- rad*zphh/2._wp)
+ zyffv = 2._wp* SIN( rad*zlam ) * TAN( rpi/4._wp- rad*zphi/2._wp) &
+ & - 2._wp* SIN( rad*zlan ) * TAN( rpi/4._wp- rad*zphh/2._wp)
znffv = SQRT( znnpv * ( zxffv*zxffv + zyffv*zyffv ) )
- znffv = MAX( znffv, 1.e-14 )
+ znffv = MAX( znffv, 1.e-14_wp )
!
zlam = plamu(ji,jj+1) ! j-direction: u-point segment direction (around f-point)
zphi = pphiu(ji,jj+1)
zlan = plamu(ji,jj )
zphh = pphiu(ji,jj )
- zxuuf = 2. * COS( rad*zlam ) * TAN( rpi/4. - rad*zphi/2. ) &
- & - 2. * COS( rad*zlan ) * TAN( rpi/4. - rad*zphh/2. )
- zyuuf = 2. * SIN( rad*zlam ) * TAN( rpi/4. - rad*zphi/2. ) &
- & - 2. * SIN( rad*zlan ) * TAN( rpi/4. - rad*zphh/2. )
+ zxuuf = 2._wp* COS( rad*zlam ) * TAN( rpi/4._wp- rad*zphi/2._wp) &
+ & - 2._wp* COS( rad*zlan ) * TAN( rpi/4._wp- rad*zphh/2._wp)
+ zyuuf = 2._wp* SIN( rad*zlam ) * TAN( rpi/4._wp- rad*zphi/2._wp) &
+ & - 2._wp* SIN( rad*zlan ) * TAN( rpi/4._wp- rad*zphh/2._wp)
znuuf = SQRT( znnpf * ( zxuuf*zxuuf + zyuuf*zyuuf ) )
- znuuf = MAX( znuuf, 1.e-14 )
+ znuuf = MAX( znuuf, 1.e-14_wp )
!
! ! cosinus and sinus using dot and cross products
gsint(ji,jj) = ( zxnpt*zyvvt - zynpt*zxvvt ) / znvvt
@@ -249,21 +249,21 @@ CONTAINS
! =============== !
DO_2D( 0, 1, 0, 0 )
- IF( MOD( ABS( plamv(ji,jj) - plamv(ji,jj-1) ), 360. ) < 1.e-8 ) THEN
- gsint(ji,jj) = 0.
- gcost(ji,jj) = 1.
+ IF( MOD( ABS( plamv(ji,jj) - plamv(ji,jj-1) ), 360._wp) < 1.e-8_wp ) THEN
+ gsint(ji,jj) = 0._wp
+ gcost(ji,jj) = 1._wp
ENDIF
- IF( MOD( ABS( plamf(ji,jj) - plamf(ji,jj-1) ), 360. ) < 1.e-8 ) THEN
- gsinu(ji,jj) = 0.
- gcosu(ji,jj) = 1.
+ IF( MOD( ABS( plamf(ji,jj) - plamf(ji,jj-1) ), 360._wp) < 1.e-8_wp ) THEN
+ gsinu(ji,jj) = 0._wp
+ gcosu(ji,jj) = 1._wp
ENDIF
- IF( ABS( pphif(ji,jj) - pphif(ji-1,jj) ) < 1.e-8 ) THEN
- gsinv(ji,jj) = 0.
+ IF( ABS( pphif(ji,jj) - pphif(ji-1,jj) ) < 1.e-8_wp ) THEN
+ gsinv(ji,jj) = 0._wp
gcosv(ji,jj) = 1.
ENDIF
- IF( MOD( ABS( plamu(ji,jj) - plamu(ji,jj+1) ), 360. ) < 1.e-8 ) THEN
- gsinf(ji,jj) = 0.
- gcosf(ji,jj) = 1.
+ IF( MOD( ABS( plamu(ji,jj) - plamu(ji,jj+1) ), 360._wp) < 1.e-8_wp ) THEN
+ gsinf(ji,jj) = 0._wp
+ gcosf(ji,jj) = 1._wp
ENDIF
END_2D
@@ -367,8 +367,8 @@ CONTAINS
CHARACTER(len=1) , INTENT( IN ) :: cgrid
REAL(wp), DIMENSION(jpi,jpj), INTENT( OUT ) :: pxx , pyy , pzz
!!
- REAL(wp), PARAMETER :: rpi = 3.141592653E0
- REAL(wp), PARAMETER :: rad = rpi / 180.e0
+ REAL(wp), PARAMETER :: rpi = 3.141592653e0_wp
+ REAL(wp), PARAMETER :: rad = rpi / 180.e0_wp
INTEGER :: ig !
INTEGER :: ierr ! local integer
!!----------------------------------------------------------------------
diff --git a/src/OCE/SBC/sbcblk.F90 b/src/OCE/SBC/sbcblk.F90
index a68659642..aa3668c33 100644
--- a/src/OCE/SBC/sbcblk.F90
+++ b/src/OCE/SBC/sbcblk.F90
@@ -369,6 +369,8 @@ CONTAINS
sf(jfpr)%fnow(:,:,1:ipka) = 101325._wp ! use standard pressure in Pa
ELSEIF( jfpr == jp_prec .OR. jfpr == jp_snow .OR. jfpr == jp_uoatm .OR. jfpr == jp_voatm ) THEN
sf(jfpr)%fnow(:,:,1:ipka) = 0._wp ! no precip or no snow or no surface currents
+ ELSEIF( jfpr == jp_wndi .OR. jfpr == jp_wndj ) THEN
+ sf(jfpr)%fnow(:,:,1:ipka) = 0._wp
ELSEIF( jfpr == jp_hpgi .OR. jfpr == jp_hpgj ) THEN
IF( .NOT. ln_abl ) THEN
DEALLOCATE( sf(jfpr)%fnow ) ! deallocate as not used in this case
@@ -866,11 +868,11 @@ CONTAINS
CALL iom_put( "vtau_oce", ztau_j(:,:)*tmask(:,:,1) ) ! vtau at T-points!
IF(sn_cfctl%l_prtctl) THEN
- CALL prt_ctl( tab2d_1=pssq , clinfo1=' blk_oce_1: pssq : ')
- CALL prt_ctl( tab2d_1=wndm , clinfo1=' blk_oce_1: wndm : ')
+ CALL prt_ctl( tab2d_1=pssq , clinfo1=' blk_oce_1: pssq : ', mask1=tmask )
+ CALL prt_ctl( tab2d_1=wndm , clinfo1=' blk_oce_1: wndm : ', mask1=tmask )
CALL prt_ctl( tab2d_1=utau , clinfo1=' blk_oce_1: utau : ', mask1=umask, &
& tab2d_2=vtau , clinfo2=' vtau : ', mask2=vmask )
- CALL prt_ctl( tab2d_1=zcd_oce, clinfo1=' blk_oce_1: Cd : ')
+ CALL prt_ctl( tab2d_1=zcd_oce, clinfo1=' blk_oce_1: Cd : ', mask1=tmask )
ENDIF
!
ENDIF ! ln_blk / ln_abl
@@ -943,8 +945,10 @@ CONTAINS
qns(:,:) = qns(:,:) * tmask(:,:,1)
!
#if defined key_si3
- qns_oce(:,:) = zqlw(:,:) + psen(:,:) + plat(:,:) ! non solar without emp (only needed by SI3)
- qsr_oce(:,:) = qsr(:,:)
+ IF ( nn_ice == 2 ) THEN
+ qns_oce(:,:) = zqlw(:,:) + psen(:,:) + plat(:,:) ! non solar without emp (only needed by SI3)
+ qsr_oce(:,:) = qsr(:,:)
+ ENDIF
#endif
!
CALL iom_put( "rho_air" , rhoa*tmask(:,:,1) ) ! output air density [kg/m^3]
@@ -965,11 +969,11 @@ CONTAINS
ENDIF
!
IF(sn_cfctl%l_prtctl) THEN
- CALL prt_ctl(tab2d_1=zqlw , clinfo1=' blk_oce_2: zqlw : ')
- CALL prt_ctl(tab2d_1=psen , clinfo1=' blk_oce_2: psen : ' )
- CALL prt_ctl(tab2d_1=plat , clinfo1=' blk_oce_2: plat : ' )
- CALL prt_ctl(tab2d_1=qns , clinfo1=' blk_oce_2: qns : ' )
- CALL prt_ctl(tab2d_1=emp , clinfo1=' blk_oce_2: emp : ')
+ CALL prt_ctl(tab2d_1=zqlw , clinfo1=' blk_oce_2: zqlw : ', mask1=tmask )
+ CALL prt_ctl(tab2d_1=psen , clinfo1=' blk_oce_2: psen : ', mask1=tmask )
+ CALL prt_ctl(tab2d_1=plat , clinfo1=' blk_oce_2: plat : ', mask1=tmask )
+ CALL prt_ctl(tab2d_1=qns , clinfo1=' blk_oce_2: qns : ', mask1=tmask )
+ CALL prt_ctl(tab2d_1=emp , clinfo1=' blk_oce_2: emp : ', mask1=tmask )
ENDIF
!
END SUBROUTINE blk_oce_2
@@ -1092,8 +1096,8 @@ CONTAINS
END_2D
CALL lbc_lnk( 'sbcblk', putaui, 'U', -1._wp, pvtaui, 'V', -1._wp )
!
- IF(sn_cfctl%l_prtctl) CALL prt_ctl( tab2d_1=putaui , clinfo1=' blk_ice: putaui : ' &
- & , tab2d_2=pvtaui , clinfo2=' pvtaui : ' )
+ IF(sn_cfctl%l_prtctl) CALL prt_ctl( tab2d_1=putaui , clinfo1=' blk_ice: putaui : ', mask1=umask &
+ & , tab2d_2=pvtaui , clinfo2=' pvtaui : ', mask2=vmask )
ELSE ! ln_abl
DO_2D( nn_hls, nn_hls, nn_hls, nn_hls )
@@ -1105,7 +1109,7 @@ CONTAINS
ENDIF ! ln_blk / ln_abl
!
- IF(sn_cfctl%l_prtctl) CALL prt_ctl(tab2d_1=wndm_ice , clinfo1=' blk_ice: wndm_ice : ')
+ IF(sn_cfctl%l_prtctl) CALL prt_ctl(tab2d_1=wndm_ice , clinfo1=' blk_ice: wndm_ice : ', mask1=tmask )
!
END SUBROUTINE blk_ice_1
@@ -1137,6 +1141,7 @@ CONTAINS
REAL(wp) :: zst, zst3, zsq, zsipt ! local variable
REAL(wp) :: zcoef_dqlw, zcoef_dqla ! - -
REAL(wp) :: zztmp, zzblk, zztmp1, z1_rLsub ! - -
+ REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: zmsk ! temporary mask for prt_ctl
REAL(wp), DIMENSION(jpi,jpj,jpl) :: z_qlw ! long wave heat flux over ice
REAL(wp), DIMENSION(jpi,jpj,jpl) :: z_qsb ! sensible heat flux over ice
REAL(wp), DIMENSION(jpi,jpj,jpl) :: z_dqlw ! long wave heat sensitivity over ice
@@ -1301,12 +1306,23 @@ CONTAINS
ENDIF
!
IF(sn_cfctl%l_prtctl) THEN
- CALL prt_ctl(tab3d_1=qla_ice , clinfo1=' blk_ice: qla_ice : ', tab3d_2=z_qsb , clinfo2=' z_qsb : ', kdim=jpl)
- CALL prt_ctl(tab3d_1=z_qlw , clinfo1=' blk_ice: z_qlw : ', tab3d_2=dqla_ice, clinfo2=' dqla_ice : ', kdim=jpl)
- CALL prt_ctl(tab3d_1=z_dqsb , clinfo1=' blk_ice: z_dqsb : ', tab3d_2=z_dqlw , clinfo2=' z_dqlw : ', kdim=jpl)
- CALL prt_ctl(tab3d_1=dqns_ice, clinfo1=' blk_ice: dqns_ice : ', tab3d_2=qsr_ice , clinfo2=' qsr_ice : ', kdim=jpl)
- CALL prt_ctl(tab3d_1=ptsu , clinfo1=' blk_ice: ptsu : ', tab3d_2=qns_ice , clinfo2=' qns_ice : ', kdim=jpl)
- CALL prt_ctl(tab2d_1=tprecip , clinfo1=' blk_ice: tprecip : ', tab2d_2=sprecip , clinfo2=' sprecip : ')
+ ALLOCATE(zmsk(jpi,jpj,jpl))
+ DO jl = 1, jpl
+ zmsk(:,:,jl) = tmask(:,:,1)
+ END DO
+ CALL prt_ctl(tab3d_1=qla_ice , clinfo1=' blk_ice: qla_ice : ', mask1=zmsk, &
+ & tab3d_2=z_qsb , clinfo2=' z_qsb : ' , mask2=zmsk, kdim=jpl)
+ CALL prt_ctl(tab3d_1=z_qlw , clinfo1=' blk_ice: z_qlw : ', mask1=zmsk, &
+ & tab3d_2=dqla_ice, clinfo2=' dqla_ice : ' , mask2=zmsk, kdim=jpl)
+ CALL prt_ctl(tab3d_1=z_dqsb , clinfo1=' blk_ice: z_dqsb : ', mask1=zmsk, &
+ & tab3d_2=z_dqlw , clinfo2=' z_dqlw : ' , mask2=zmsk, kdim=jpl)
+ CALL prt_ctl(tab3d_1=dqns_ice, clinfo1=' blk_ice: dqns_ice : ', mask1=zmsk, &
+ & tab3d_2=qsr_ice , clinfo2=' qsr_ice : ' , mask2=zmsk, kdim=jpl)
+ CALL prt_ctl(tab3d_1=ptsu , clinfo1=' blk_ice: ptsu : ', mask1=zmsk, &
+ & tab3d_2=qns_ice , clinfo2=' qns_ice : ' , mask2=zmsk, kdim=jpl)
+ CALL prt_ctl(tab2d_1=tprecip , clinfo1=' blk_ice: tprecip : ', mask1=tmask, &
+ & tab2d_2=sprecip , clinfo2=' sprecip : ' , mask2=tmask )
+ DEALLOCATE(zmsk)
ENDIF
!#LB:
diff --git a/src/OCE/SBC/sbcfwb.F90 b/src/OCE/SBC/sbcfwb.F90
index cc5ca7b32..58a17f028 100644
--- a/src/OCE/SBC/sbcfwb.F90
+++ b/src/OCE/SBC/sbcfwb.F90
@@ -104,11 +104,11 @@ CONTAINS
! isf cavities are excluded because it can feedback to the melting with generation of inhibition of plumes
! and in case of no melt, it can generate HSSW.
!
-#if ! defined key_si3 && ! defined key_cice
- snwice_mass_b(:,:) = 0.e0 ! no sea-ice model is being used : no snow+ice mass
- snwice_mass (:,:) = 0.e0
- snwice_fmass (:,:) = 0.e0
-#endif
+ IF( nn_ice == 0 ) THEN
+ snwice_mass_b(:,:) = 0.e0 ! no sea-ice model is being used : no snow+ice mass
+ snwice_mass (:,:) = 0.e0
+ snwice_fmass (:,:) = 0.e0
+ ENDIF
!
ENDIF
diff --git a/src/OCE/SBC/sbcmod.F90 b/src/OCE/SBC/sbcmod.F90
index 220693366..30ee55c7c 100644
--- a/src/OCE/SBC/sbcmod.F90
+++ b/src/OCE/SBC/sbcmod.F90
@@ -39,6 +39,7 @@ MODULE sbcmod
USE sbcice_if ! surface boundary condition: ice-if sea-ice model
#if defined key_si3
USE icestp ! surface boundary condition: SI3 sea-ice model
+ USE ice
#endif
USE sbcice_cice ! surface boundary condition: CICE sea-ice model
USE sbccpl ! surface boundary condition: coupled formulation
@@ -325,8 +326,14 @@ CONTAINS
IF( ln_apr_dyn ) CALL sbc_apr_init ! Atmo Pressure Forcing initialization
!
#if defined key_si3
- IF( lk_agrif .AND. nn_ice == 0 ) THEN ! allocate ice arrays in case agrif + ice-model + no-ice in child grid
- IF( sbc_ice_alloc() /= 0 ) CALL ctl_stop('STOP', 'sbc_ice_alloc : unable to allocate arrays' )
+ IF( nn_ice == 0 ) THEN
+ ! allocate ice arrays in case agrif + ice-model + no-ice in child grid
+ jpl = 1 ; nlay_i = 1 ; nlay_s = 1
+ IF( sbc_ice_alloc() /= 0 ) CALL ctl_stop('STOP', 'sbc_ice_alloc : unable to allocate arrays' )
+#if defined key_agrif
+ CALL Agrif_Declare_Var_ice ! " " " " " Sea ice
+#endif
+
ELSEIF( nn_ice == 2 ) THEN
CALL ice_init( Kbb, Kmm, Kaa ) ! ICE initialization
ENDIF
diff --git a/src/OCE/SBC/sbcrnf.F90 b/src/OCE/SBC/sbcrnf.F90
index be3f221b8..758d4f4ba 100644
--- a/src/OCE/SBC/sbcrnf.F90
+++ b/src/OCE/SBC/sbcrnf.F90
@@ -374,9 +374,9 @@ CONTAINS
IF( .NOT. sn_dep_rnf%ln_clim ) THEN ; WRITE(rn_dep_file, '(a,"_y",i4)' ) TRIM( rn_dep_file ), nyear ! add year
IF( sn_dep_rnf%clftyp == 'monthly' ) WRITE(rn_dep_file, '(a,"m",i2)' ) TRIM( rn_dep_file ), nmonth ! add month
ENDIF
- CALL iom_open ( rn_dep_file, inum ) ! open file
- CALL iom_get ( inum, jpdom_global, sn_dep_rnf%clvar, h_rnf ) ! read the river mouth array
- CALL iom_close( inum ) ! close file
+ CALL iom_open ( rn_dep_file, inum ) ! open file
+ CALL iom_get ( inum, jpdom_global, sn_dep_rnf%clvar, h_rnf, kfill = jpfillcopy ) ! read the river mouth. no 0 on halos!
+ CALL iom_close( inum ) ! close file
!
nk_rnf(:,:) = 0 ! set the number of level over which river runoffs are applied
DO_2D( nn_hls, nn_hls, nn_hls, nn_hls )
diff --git a/src/OCE/SBC/sbcssm.F90 b/src/OCE/SBC/sbcssm.F90
index 7ff88f827..1db7fb00a 100644
--- a/src/OCE/SBC/sbcssm.F90
+++ b/src/OCE/SBC/sbcssm.F90
@@ -239,7 +239,7 @@ CONTAINS
!
IF( zf_sbc /= REAL( nn_fsbc, wp ) ) THEN ! nn_fsbc has changed between 2 runs
IF(lwp) WRITE(numout,*) ' restart with a change in the frequency of mean from ', zf_sbc, ' to ', nn_fsbc
- zcoef = REAL( nn_fsbc - 1, wp ) / zf_sbc
+ zcoef = REAL( nn_fsbc - 1, wp ) / ( zf_sbc - 1._wp ) ! zf_sbc /= 1 as it was written in the restart
ssu_m(:,:) = zcoef * ssu_m(:,:)
ssv_m(:,:) = zcoef * ssv_m(:,:)
sst_m(:,:) = zcoef * sst_m(:,:)
diff --git a/src/OCE/TRA/eosbn2.F90 b/src/OCE/TRA/eosbn2.F90
index fcefee7a7..81fda25e0 100644
--- a/src/OCE/TRA/eosbn2.F90
+++ b/src/OCE/TRA/eosbn2.F90
@@ -44,10 +44,8 @@ MODULE eosbn2
USE stopts ! Stochastic T/S fluctuations
!
USE in_out_manager ! I/O manager
- USE lib_mpp ! MPP library
- USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined)
+ USE lib_mpp ! for ctl_stop
USE prtctl ! Print control
- USE lbclnk ! ocean lateral boundary conditions
USE timing ! Timing
IMPLICIT NONE
diff --git a/src/OCE/TRA/traldf_lap_blp.F90 b/src/OCE/TRA/traldf_lap_blp.F90
index 3113c0fcb..16e5df16c 100644
--- a/src/OCE/TRA/traldf_lap_blp.F90
+++ b/src/OCE/TRA/traldf_lap_blp.F90
@@ -240,8 +240,8 @@ CONTAINS
!
IF (nn_hls==1) CALL lbc_lnk( 'traldf_lap_blp', zlap(:,:,:,:) , 'T', 1.0_wp ) ! Lateral boundary conditions (unchanged sign)
! ! Partial top/bottom cell: GRADh( zlap )
- IF( ln_isfcav .AND. ln_zps ) THEN ; CALL zps_hde_isf( kt, Kmm, kjpt, zlap, zglu, zglv, zgui, zgvi ) ! both top & bottom
- ELSEIF( ln_zps ) THEN ; CALL zps_hde ( kt, Kmm, kjpt, zlap, zglu, zglv ) ! only bottom
+ IF( ln_isfcav .AND. ln_zps ) THEN ; CALL zps_hde_isf( kt, kjpt, zlap, zglu, zglv, zgui, zgvi ) ! both top & bottom
+ ELSEIF( ln_zps ) THEN ; CALL zps_hde ( kt, kjpt, zlap, zglu, zglv ) ! only bottom
ENDIF
!
SELECT CASE ( kldf ) !== 2nd laplacian applied to zlap (output in pt_rhs) ==!
diff --git a/src/OCE/TRA/traldf_triad.F90 b/src/OCE/TRA/traldf_triad.F90
index 84ccde88b..70eddec54 100644
--- a/src/OCE/TRA/traldf_triad.F90
+++ b/src/OCE/TRA/traldf_triad.F90
@@ -21,7 +21,6 @@ MODULE traldf_triad
USE traldf_iso ! lateral diffusion (Madec operator) (tra_ldf_iso routine)
USE diaptr ! poleward transport diagnostics
USE diaar5 ! AR5 diagnostics
- USE zpshde ! partial step: hor. derivative (zps_hde routine)
!
USE in_out_manager ! I/O manager
USE iom ! I/O library
diff --git a/src/OCE/TRA/zpshde.F90 b/src/OCE/TRA/zpshde.F90
index 67d8dca6f..2b786d65e 100644
--- a/src/OCE/TRA/zpshde.F90
+++ b/src/OCE/TRA/zpshde.F90
@@ -40,11 +40,9 @@ MODULE zpshde
!!----------------------------------------------------------------------
CONTAINS
- SUBROUTINE zps_hde( kt, Kmm, kjpt, pta, pgtu, pgtv, &
- & prd, pgru, pgrv )
+ SUBROUTINE zps_hde( kt, kjpt, pta, pgtu, pgtv, prd, pgru, pgrv )
!!
INTEGER , INTENT(in ) :: kt ! ocean time-step index
- INTEGER , INTENT(in ) :: Kmm ! ocean time level index
INTEGER , INTENT(in ) :: kjpt ! number of tracers
REAL(wp), DIMENSION(:,:,:,:), INTENT(in ) :: pta ! 4D tracers fields
REAL(wp), DIMENSION(:,:,:) , INTENT( out) :: pgtu, pgtv ! hor. grad. of ptra at u- & v-pts
@@ -56,13 +54,13 @@ CONTAINS
IF( PRESENT(prd) ) THEN ; itrd = is_tile(prd) ; ELSE ; itrd = 0 ; ENDIF
IF( PRESENT(pgru) ) THEN ; itgr = is_tile(pgru) ; ELSE ; itgr = 0 ; ENDIF
- CALL zps_hde_t( kt, Kmm, kjpt, pta, is_tile(pta), pgtu, pgtv, is_tile(pgtu), &
- & prd, itrd, pgru, pgrv, itgr )
+ CALL zps_hde_t( kt, kjpt, pta, is_tile(pta), pgtu, pgtv, is_tile(pgtu), &
+ & prd, itrd, pgru, pgrv, itgr )
END SUBROUTINE zps_hde
- SUBROUTINE zps_hde_t( kt, Kmm, kjpt, pta, ktta, pgtu, pgtv, ktgt, &
- & prd, ktrd, pgru, pgrv, ktgr )
+ SUBROUTINE zps_hde_t( kt, kjpt, pta, ktta, pgtu, pgtv, ktgt, &
+ & prd, ktrd, pgru, pgrv, ktgr )
!!----------------------------------------------------------------------
!! *** ROUTINE zps_hde ***
!!
@@ -87,13 +85,13 @@ CONTAINS
!! | |____ ____| |
!! ___ | | | ___ | | |
!!
- !! case 1-> e3w(i+1,:,:,Kmm) >= e3w(i,:,:,Kmm) ( and e3w(:,j+1,:,Kmm) >= e3w(:,j,:,Kmm) ) then
- !! t~ = t(i+1,j ,k) + (e3w(i+1,j,k,Kmm) - e3w(i,j,k,Kmm)) * dk(Ti+1)/e3w(i+1,j,k,Kmm)
- !! ( t~ = t(i ,j+1,k) + (e3w(i,j+1,k,Kmm) - e3w(i,j,k,Kmm)) * dk(Tj+1)/e3w(i,j+1,k,Kmm) )
+ !! case 1-> e3w_0(i+1,:,:) >= e3w_0(i,:,:) ( and e3w_0(:,j+1,:) >= e3w_0(:,j,:) ) then
+ !! t~ = t(i+1,j ,k) + (e3w_0(i+1,j,k) - e3w_0(i,j,k)) * dk(Ti+1)/e3w_0(i+1,j,k)
+ !! ( t~ = t(i ,j+1,k) + (e3w_0(i,j+1,k) - e3w_0(i,j,k)) * dk(Tj+1)/e3w_0(i,j+1,k) )
!! or
- !! case 2-> e3w(i+1,:,:,Kmm) <= e3w(i,:,:,Kmm) ( and e3w(:,j+1,:,Kmm) <= e3w(:,j,:,Kmm) ) then
- !! t~ = t(i,j,k) + (e3w(i,j,k,Kmm) - e3w(i+1,j,k,Kmm)) * dk(Ti)/e3w(i,j,k,Kmm)
- !! ( t~ = t(i,j,k) + (e3w(i,j,k,Kmm) - e3w(i,j+1,k,Kmm)) * dk(Tj)/e3w(i,j,k,Kmm) )
+ !! case 2-> e3w_0(i+1,:,:) <= e3w_0(i,:,:) ( and e3w_0(:,j+1,:) <= e3w_0(:,j,:) ) then
+ !! t~ = t(i,j,k) + (e3w_0(i,j,k) - e3w_0(i+1,j,k)) * dk(Ti)/e3w_0(i,j,k)
+ !! ( t~ = t(i,j,k) + (e3w_0(i,j,k) - e3w_0(i,j+1,k)) * dk(Tj)/e3w_0(i,j,k) )
!! Idem for di(s) and dj(s)
!!
!! For rho, we call eos which will compute rd~(t~,s~) at the right
@@ -107,7 +105,6 @@ CONTAINS
!! - pgru, pgrv: horizontal gradient of rho (if present) at u- & v-points
!!----------------------------------------------------------------------
INTEGER , INTENT(in ) :: kt ! ocean time-step index
- INTEGER , INTENT(in ) :: Kmm ! ocean time level index
INTEGER , INTENT(in ) :: kjpt ! number of tracers
INTEGER , INTENT(in ) :: ktta, ktgt, ktrd, ktgr
REAL(wp), DIMENSION(A2D_T(ktta),JPK,KJPT), INTENT(in ) :: pta ! 4D tracers fields
@@ -132,19 +129,18 @@ CONTAINS
DO_2D( nn_hls, nn_hls-1, nn_hls, nn_hls-1 ) ! Gradient of density at the last level
iku = mbku(ji,jj) ; ikum1 = MAX( iku - 1 , 1 ) ! last and before last ocean level at u- & v-points
ikv = mbkv(ji,jj) ; ikvm1 = MAX( ikv - 1 , 1 ) ! if level first is a p-step, ik.m1=1
-!!gm BUG ? when applied to before fields, e3w(:,:,k,Kbb) should be used....
- ze3wu = e3w(ji+1,jj ,iku,Kmm) - e3w(ji,jj,iku,Kmm)
- ze3wv = e3w(ji ,jj+1,ikv,Kmm) - e3w(ji,jj,ikv,Kmm)
+ ze3wu = e3w_0(ji+1,jj ,iku) - e3w_0(ji,jj,iku)
+ ze3wv = e3w_0(ji ,jj+1,ikv) - e3w_0(ji,jj,ikv)
!
! i- direction
IF( ze3wu >= 0._wp ) THEN ! case 1
- zmaxu = ze3wu / e3w(ji+1,jj,iku,Kmm)
+ zmaxu = ze3wu / e3w_0(ji+1,jj,iku)
! interpolated values of tracers
zti (ji,jj,jn) = pta(ji+1,jj,iku,jn) + zmaxu * ( pta(ji+1,jj,ikum1,jn) - pta(ji+1,jj,iku,jn) )
! gradient of tracers
pgtu(ji,jj,jn) = umask(ji,jj,1) * ( zti(ji,jj,jn) - pta(ji,jj,iku,jn) )
ELSE ! case 2
- zmaxu = -ze3wu / e3w(ji,jj,iku,Kmm)
+ zmaxu = -ze3wu / e3w_0(ji,jj,iku)
! interpolated values of tracers
zti (ji,jj,jn) = pta(ji,jj,iku,jn) + zmaxu * ( pta(ji,jj,ikum1,jn) - pta(ji,jj,iku,jn) )
! gradient of tracers
@@ -153,13 +149,13 @@ CONTAINS
!
! j- direction
IF( ze3wv >= 0._wp ) THEN ! case 1
- zmaxv = ze3wv / e3w(ji,jj+1,ikv,Kmm)
+ zmaxv = ze3wv / e3w_0(ji,jj+1,ikv)
! interpolated values of tracers
ztj (ji,jj,jn) = pta(ji,jj+1,ikv,jn) + zmaxv * ( pta(ji,jj+1,ikvm1,jn) - pta(ji,jj+1,ikv,jn) )
! gradient of tracers
pgtv(ji,jj,jn) = vmask(ji,jj,1) * ( ztj(ji,jj,jn) - pta(ji,jj,ikv,jn) )
ELSE ! case 2
- zmaxv = -ze3wv / e3w(ji,jj,ikv,Kmm)
+ zmaxv = -ze3wv / e3w_0(ji,jj,ikv)
! interpolated values of tracers
ztj (ji,jj,jn) = pta(ji,jj,ikv,jn) + zmaxv * ( pta(ji,jj,ikvm1,jn) - pta(ji,jj,ikv,jn) )
! gradient of tracers
@@ -176,13 +172,13 @@ CONTAINS
DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 )
iku = mbku(ji,jj)
ikv = mbkv(ji,jj)
- ze3wu = e3w(ji+1,jj ,iku,Kmm) - e3w(ji,jj,iku,Kmm)
- ze3wv = e3w(ji ,jj+1,ikv,Kmm) - e3w(ji,jj,ikv,Kmm)
- IF( ze3wu >= 0._wp ) THEN ; zhi(ji,jj) = gdept(ji ,jj,iku,Kmm) ! i-direction: case 1
- ELSE ; zhi(ji,jj) = gdept(ji+1,jj,iku,Kmm) ! - - case 2
+ ze3wu = e3w_0(ji+1,jj ,iku) - e3w_0(ji,jj,iku)
+ ze3wv = e3w_0(ji ,jj+1,ikv) - e3w_0(ji,jj,ikv)
+ IF( ze3wu >= 0._wp ) THEN ; zhi(ji,jj) = gdept_0(ji ,jj,iku) ! i-direction: case 1
+ ELSE ; zhi(ji,jj) = gdept_0(ji+1,jj,iku) ! - - case 2
ENDIF
- IF( ze3wv >= 0._wp ) THEN ; zhj(ji,jj) = gdept(ji,jj ,ikv,Kmm) ! j-direction: case 1
- ELSE ; zhj(ji,jj) = gdept(ji,jj+1,ikv,Kmm) ! - - case 2
+ IF( ze3wv >= 0._wp ) THEN ; zhj(ji,jj) = gdept_0(ji,jj ,ikv) ! j-direction: case 1
+ ELSE ; zhj(ji,jj) = gdept_0(ji,jj+1,ikv) ! - - case 2
ENDIF
END_2D
!
@@ -192,8 +188,8 @@ CONTAINS
DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 ) ! Gradient of density at the last level
iku = mbku(ji,jj)
ikv = mbkv(ji,jj)
- ze3wu = e3w(ji+1,jj ,iku,Kmm) - e3w(ji,jj,iku,Kmm)
- ze3wv = e3w(ji ,jj+1,ikv,Kmm) - e3w(ji,jj,ikv,Kmm)
+ ze3wu = e3w_0(ji+1,jj ,iku) - e3w_0(ji,jj,iku)
+ ze3wv = e3w_0(ji ,jj+1,ikv) - e3w_0(ji,jj,ikv)
IF( ze3wu >= 0._wp ) THEN ; pgru(ji,jj) = umask(ji,jj,1) * ( zri(ji ,jj ) - prd(ji,jj,iku) ) ! i: 1
ELSE ; pgru(ji,jj) = umask(ji,jj,1) * ( prd(ji+1,jj,iku) - zri(ji,jj ) ) ! i: 2
ENDIF
@@ -210,11 +206,10 @@ CONTAINS
END SUBROUTINE zps_hde_t
- SUBROUTINE zps_hde_isf( kt, Kmm, kjpt, pta, pgtu, pgtv, pgtui, pgtvi, &
- & prd, pgru, pgrv, pgrui, pgrvi )
+ SUBROUTINE zps_hde_isf( kt, kjpt, pta, pgtu, pgtv, pgtui, pgtvi, &
+ & prd, pgru, pgrv, pgrui, pgrvi )
!!
INTEGER , INTENT(in ) :: kt ! ocean time-step index
- INTEGER , INTENT(in ) :: Kmm ! ocean time level index
INTEGER , INTENT(in ) :: kjpt ! number of tracers
REAL(wp), DIMENSION(:,:,:,:), INTENT(in ) :: pta ! 4D tracers fields
REAL(wp), DIMENSION(:,:,:) , INTENT( out) :: pgtu, pgtv ! hor. grad. of ptra at u- & v-pts
@@ -229,13 +224,13 @@ CONTAINS
IF( PRESENT(pgru) ) THEN ; itgr = is_tile(pgru) ; ELSE ; itgr = 0 ; ENDIF
IF( PRESENT(pgrui) ) THEN ; itgri = is_tile(pgrui) ; ELSE ; itgri = 0 ; ENDIF
- CALL zps_hde_isf_t( kt, Kmm, kjpt, pta, is_tile(pta), pgtu, pgtv, is_tile(pgtu), pgtui, pgtvi, is_tile(pgtui), &
- & prd, itrd, pgru, pgrv, itgr, pgrui, pgrvi, itgri )
+ CALL zps_hde_isf_t( kt, kjpt, pta, is_tile(pta), pgtu, pgtv, is_tile(pgtu), pgtui, pgtvi, is_tile(pgtui), &
+ & prd, itrd, pgru, pgrv, itgr, pgrui, pgrvi, itgri )
END SUBROUTINE zps_hde_isf
- SUBROUTINE zps_hde_isf_t( kt, Kmm, kjpt, pta, ktta, pgtu, pgtv, ktgt, pgtui, pgtvi, ktgti, &
- & prd, ktrd, pgru, pgrv, ktgr, pgrui, pgrvi, ktgri )
+ SUBROUTINE zps_hde_isf_t( kt, kjpt, pta, ktta, pgtu, pgtv, ktgt, pgtui, pgtvi, ktgti, &
+ & prd, ktrd, pgru, pgrv, ktgr, pgrui, pgrvi, ktgri )
!!----------------------------------------------------------------------
!! *** ROUTINE zps_hde_isf ***
!!
@@ -261,13 +256,13 @@ CONTAINS
!! | |____ ____| |
!! ___ | | | ___ | | |
!!
- !! case 1-> e3w(i+1,j,k,Kmm) >= e3w(i,j,k,Kmm) ( and e3w(i,j+1,k,Kmm) >= e3w(i,j,k,Kmm) ) then
- !! t~ = t(i+1,j ,k) + (e3w(i+1,j ,k,Kmm) - e3w(i,j,k,Kmm)) * dk(Ti+1)/e3w(i+1,j ,k,Kmm)
- !! ( t~ = t(i ,j+1,k) + (e3w(i ,j+1,k,Kmm) - e3w(i,j,k,Kmm)) * dk(Tj+1)/e3w(i ,j+1,k,Kmm) )
+ !! case 1-> e3w_0(i+1,j,k) >= e3w_0(i,j,k) ( and e3w_0(i,j+1,k) >= e3w_0(i,j,k) ) then
+ !! t~ = t(i+1,j ,k) + (e3w_0(i+1,j ,k) - e3w_0(i,j,k)) * dk(Ti+1)/e3w_0(i+1,j ,k)
+ !! ( t~ = t(i ,j+1,k) + (e3w_0(i ,j+1,k) - e3w_0(i,j,k)) * dk(Tj+1)/e3w_0(i ,j+1,k) )
!! or
- !! case 2-> e3w(i+1,j,k,Kmm) <= e3w(i,j,k,Kmm) ( and e3w(i,j+1,k,Kmm) <= e3w(i,j,k,Kmm) ) then
- !! t~ = t(i,j,k) + (e3w(i,j,k,Kmm) - e3w(i+1,j ,k,Kmm)) * dk(Ti)/e3w(i,j,k,Kmm)
- !! ( t~ = t(i,j,k) + (e3w(i,j,k,Kmm) - e3w(i ,j+1,k,Kmm)) * dk(Tj)/e3w(i,j,k,Kmm) )
+ !! case 2-> e3w_0(i+1,j,k) <= e3w_0(i,j,k) ( and e3w_0(i,j+1,k) <= e3w_0(i,j,k) ) then
+ !! t~ = t(i,j,k) + (e3w_0(i,j,k) - e3w_0(i+1,j ,k)) * dk(Ti)/e3w_0(i,j,k)
+ !! ( t~ = t(i,j,k) + (e3w_0(i,j,k) - e3w_0(i ,j+1,k)) * dk(Tj)/e3w_0(i,j,k) )
!! Idem for di(s) and dj(s)
!!
!! For rho, we call eos which will compute rd~(t~,s~) at the right
@@ -283,7 +278,6 @@ CONTAINS
!! - pgru, pgrv, pgrui, pgtvi: horizontal gradient of rho (if present) at u- & v-points
!!----------------------------------------------------------------------
INTEGER , INTENT(in ) :: kt ! ocean time-step index
- INTEGER , INTENT(in ) :: Kmm ! ocean time level index
INTEGER , INTENT(in ) :: kjpt ! number of tracers
INTEGER , INTENT(in ) :: ktta, ktgt, ktgti, ktrd, ktgr, ktgri
REAL(wp), DIMENSION(A2D_T(ktta),JPK,KJPT), INTENT(in ) :: pta ! 4D tracers fields
@@ -313,18 +307,18 @@ CONTAINS
iku = mbku(ji,jj); ikum1 = MAX( iku - 1 , 1 ) ! last and before last ocean level at u- & v-points
ikv = mbkv(ji,jj); ikvm1 = MAX( ikv - 1 , 1 ) ! if level first is a p-step, ik.m1=1
- ze3wu = gdept(ji+1,jj,iku,Kmm) - gdept(ji,jj,iku,Kmm)
- ze3wv = gdept(ji,jj+1,ikv,Kmm) - gdept(ji,jj,ikv,Kmm)
+ ze3wu = gdept_0(ji+1,jj,iku) - gdept_0(ji,jj,iku)
+ ze3wv = gdept_0(ji,jj+1,ikv) - gdept_0(ji,jj,ikv)
!
! i- direction
IF( ze3wu >= 0._wp ) THEN ! case 1
- zmaxu = ze3wu / e3w(ji+1,jj,iku,Kmm)
+ zmaxu = ze3wu / e3w_0(ji+1,jj,iku)
! interpolated values of tracers
zti (ji,jj,jn) = pta(ji+1,jj,iku,jn) + zmaxu * ( pta(ji+1,jj,ikum1,jn) - pta(ji+1,jj,iku,jn) )
! gradient of tracers
pgtu(ji,jj,jn) = ssumask(ji,jj) * ( zti(ji,jj,jn) - pta(ji,jj,iku,jn) )
ELSE ! case 2
- zmaxu = -ze3wu / e3w(ji,jj,iku,Kmm)
+ zmaxu = -ze3wu / e3w_0(ji,jj,iku)
! interpolated values of tracers
zti (ji,jj,jn) = pta(ji,jj,iku,jn) + zmaxu * ( pta(ji,jj,ikum1,jn) - pta(ji,jj,iku,jn) )
! gradient of tracers
@@ -333,13 +327,13 @@ CONTAINS
!
! j- direction
IF( ze3wv >= 0._wp ) THEN ! case 1
- zmaxv = ze3wv / e3w(ji,jj+1,ikv,Kmm)
+ zmaxv = ze3wv / e3w_0(ji,jj+1,ikv)
! interpolated values of tracers
ztj (ji,jj,jn) = pta(ji,jj+1,ikv,jn) + zmaxv * ( pta(ji,jj+1,ikvm1,jn) - pta(ji,jj+1,ikv,jn) )
! gradient of tracers
pgtv(ji,jj,jn) = ssvmask(ji,jj) * ( ztj(ji,jj,jn) - pta(ji,jj,ikv,jn) )
ELSE ! case 2
- zmaxv = -ze3wv / e3w(ji,jj,ikv,Kmm)
+ zmaxv = -ze3wv / e3w_0(ji,jj,ikv)
! interpolated values of tracers
ztj (ji,jj,jn) = pta(ji,jj,ikv,jn) + zmaxv * ( pta(ji,jj,ikvm1,jn) - pta(ji,jj,ikv,jn) )
! gradient of tracers
@@ -359,14 +353,14 @@ CONTAINS
iku = mbku(ji,jj)
ikv = mbkv(ji,jj)
- ze3wu = gdept(ji+1,jj,iku,Kmm) - gdept(ji,jj,iku,Kmm)
- ze3wv = gdept(ji,jj+1,ikv,Kmm) - gdept(ji,jj,ikv,Kmm)
+ ze3wu = gdept_0(ji+1,jj,iku) - gdept_0(ji,jj,iku)
+ ze3wv = gdept_0(ji,jj+1,ikv) - gdept_0(ji,jj,ikv)
!
- IF( ze3wu >= 0._wp ) THEN ; zhi(ji,jj) = gdept(ji ,jj,iku,Kmm) ! i-direction: case 1
- ELSE ; zhi(ji,jj) = gdept(ji+1,jj,iku,Kmm) ! - - case 2
+ IF( ze3wu >= 0._wp ) THEN ; zhi(ji,jj) = gdept_0(ji ,jj,iku) ! i-direction: case 1
+ ELSE ; zhi(ji,jj) = gdept_0(ji+1,jj,iku) ! - - case 2
ENDIF
- IF( ze3wv >= 0._wp ) THEN ; zhj(ji,jj) = gdept(ji,jj ,ikv,Kmm) ! j-direction: case 1
- ELSE ; zhj(ji,jj) = gdept(ji,jj+1,ikv,Kmm) ! - - case 2
+ IF( ze3wv >= 0._wp ) THEN ; zhj(ji,jj) = gdept_0(ji,jj ,ikv) ! j-direction: case 1
+ ELSE ; zhj(ji,jj) = gdept_0(ji,jj+1,ikv) ! - - case 2
ENDIF
END_2D
@@ -379,8 +373,8 @@ CONTAINS
DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 )
iku = mbku(ji,jj)
ikv = mbkv(ji,jj)
- ze3wu = gdept(ji+1,jj,iku,Kmm) - gdept(ji,jj,iku,Kmm)
- ze3wv = gdept(ji,jj+1,ikv,Kmm) - gdept(ji,jj,ikv,Kmm)
+ ze3wu = gdept_0(ji+1,jj,iku) - gdept_0(ji,jj,iku)
+ ze3wv = gdept_0(ji,jj+1,ikv) - gdept_0(ji,jj,ikv)
IF( ze3wu >= 0._wp ) THEN ; pgru(ji,jj) = ssumask(ji,jj) * ( zri(ji ,jj ) - prd(ji,jj,iku) ) ! i: 1
ELSE ; pgru(ji,jj) = ssumask(ji,jj) * ( prd(ji+1,jj,iku) - zri(ji,jj ) ) ! i: 2
@@ -404,20 +398,20 @@ CONTAINS
!
! (ISF) case partial step top and bottom in adjacent cell in vertical
! cannot used e3w because if 2 cell water column, we have ps at top and bottom
- ! in this case e3w(i,j,k,Kmm) - e3w(i,j+1,k,Kmm) is not the distance between Tj~ and Tj
+ ! in this case e3w_0(i,j,k) - e3w_0(i,j+1,k) is not the distance between Tj~ and Tj
! the only common depth between cells (i,j) and (i,j+1) is gdepw_0
- ze3wu = gdept(ji,jj,iku,Kmm) - gdept(ji+1,jj,iku,Kmm)
- ze3wv = gdept(ji,jj,ikv,Kmm) - gdept(ji,jj+1,ikv,Kmm)
+ ze3wu = gdept_0(ji,jj,iku) - gdept_0(ji+1,jj,iku)
+ ze3wv = gdept_0(ji,jj,ikv) - gdept_0(ji,jj+1,ikv)
! i- direction
IF( ze3wu >= 0._wp ) THEN ! case 1
- zmaxu = ze3wu / e3w(ji+1,jj,ikup1,Kmm)
+ zmaxu = ze3wu / e3w_0(ji+1,jj,ikup1)
! interpolated values of tracers
zti(ji,jj,jn) = pta(ji+1,jj,iku,jn) + zmaxu * ( pta(ji+1,jj,ikup1,jn) - pta(ji+1,jj,iku,jn) )
! gradient of tracers
pgtui(ji,jj,jn) = ssumask(ji,jj) * ( zti(ji,jj,jn) - pta(ji,jj,iku,jn) )
ELSE ! case 2
- zmaxu = - ze3wu / e3w(ji,jj,ikup1,Kmm)
+ zmaxu = - ze3wu / e3w_0(ji,jj,ikup1)
! interpolated values of tracers
zti(ji,jj,jn) = pta(ji,jj,iku,jn) + zmaxu * ( pta(ji,jj,ikup1,jn) - pta(ji,jj,iku,jn) )
! gradient of tracers
@@ -426,13 +420,13 @@ CONTAINS
!
! j- direction
IF( ze3wv >= 0._wp ) THEN ! case 1
- zmaxv = ze3wv / e3w(ji,jj+1,ikvp1,Kmm)
+ zmaxv = ze3wv / e3w_0(ji,jj+1,ikvp1)
! interpolated values of tracers
ztj(ji,jj,jn) = pta(ji,jj+1,ikv,jn) + zmaxv * ( pta(ji,jj+1,ikvp1,jn) - pta(ji,jj+1,ikv,jn) )
! gradient of tracers
pgtvi(ji,jj,jn) = ssvmask(ji,jj) * ( ztj(ji,jj,jn) - pta(ji,jj,ikv,jn) )
ELSE ! case 2
- zmaxv = - ze3wv / e3w(ji,jj,ikvp1,Kmm)
+ zmaxv = - ze3wv / e3w_0(ji,jj,ikvp1)
! interpolated values of tracers
ztj(ji,jj,jn) = pta(ji,jj,ikv,jn) + zmaxv * ( pta(ji,jj,ikvp1,jn) - pta(ji,jj,ikv,jn) )
! gradient of tracers
@@ -451,15 +445,15 @@ CONTAINS
iku = miku(ji,jj)
ikv = mikv(ji,jj)
- ze3wu = gdept(ji,jj,iku,Kmm) - gdept(ji+1,jj,iku,Kmm)
- ze3wv = gdept(ji,jj,ikv,Kmm) - gdept(ji,jj+1,ikv,Kmm)
+ ze3wu = gdept_0(ji,jj,iku) - gdept_0(ji+1,jj,iku)
+ ze3wv = gdept_0(ji,jj,ikv) - gdept_0(ji,jj+1,ikv)
!
- IF( ze3wu >= 0._wp ) THEN ; zhi(ji,jj) = gdept(ji ,jj,iku,Kmm) ! i-direction: case 1
- ELSE ; zhi(ji,jj) = gdept(ji+1,jj,iku,Kmm) ! - - case 2
+ IF( ze3wu >= 0._wp ) THEN ; zhi(ji,jj) = gdept_0(ji ,jj,iku) ! i-direction: case 1
+ ELSE ; zhi(ji,jj) = gdept_0(ji+1,jj,iku) ! - - case 2
ENDIF
- IF( ze3wv >= 0._wp ) THEN ; zhj(ji,jj) = gdept(ji,jj ,ikv,Kmm) ! j-direction: case 1
- ELSE ; zhj(ji,jj) = gdept(ji,jj+1,ikv,Kmm) ! - - case 2
+ IF( ze3wv >= 0._wp ) THEN ; zhj(ji,jj) = gdept_0(ji,jj ,ikv) ! j-direction: case 1
+ ELSE ; zhj(ji,jj) = gdept_0(ji,jj+1,ikv) ! - - case 2
ENDIF
END_2D
@@ -470,8 +464,8 @@ CONTAINS
DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 )
iku = miku(ji,jj)
ikv = mikv(ji,jj)
- ze3wu = gdept(ji,jj,iku,Kmm) - gdept(ji+1,jj,iku,Kmm)
- ze3wv = gdept(ji,jj,ikv,Kmm) - gdept(ji,jj+1,ikv,Kmm)
+ ze3wu = gdept_0(ji,jj,iku) - gdept_0(ji+1,jj,iku)
+ ze3wv = gdept_0(ji,jj,ikv) - gdept_0(ji,jj+1,ikv)
IF( ze3wu >= 0._wp ) THEN ; pgrui(ji,jj) = ssumask(ji,jj) * ( zri(ji ,jj ) - prd(ji,jj,iku) ) ! i: 1
ELSE ; pgrui(ji,jj) = ssumask(ji,jj) * ( prd(ji+1,jj ,iku) - zri(ji,jj ) ) ! i: 2
diff --git a/src/OCE/ZDF/zdfgls.F90 b/src/OCE/ZDF/zdfgls.F90
index e3f7c96e8..b25670799 100644
--- a/src/OCE/ZDF/zdfgls.F90
+++ b/src/OCE/ZDF/zdfgls.F90
@@ -967,8 +967,13 @@ CONTAINS
CASE( 2 ) ; WRITE(numout,*) ' ==>>> scaling with mean sea-ice thickness'
CASE( 3 ) ; WRITE(numout,*) ' ==>>> scaling with max sea-ice thickness'
CASE DEFAULT
- CALL ctl_stop( 'zdf_tke_init: wrong value for nn_mxlice, should be 0,1,2,3 ')
+ CALL ctl_stop( 'zdf_gls_init: wrong value for nn_mxlice, should be 0,1,2,3 ')
END SELECT
+ IF ( (nn_mxlice>0).AND.(nn_ice==0) ) THEN
+ CALL ctl_stop( 'zdf_gls_init: with no ice at all, nn_mxlice must be 0 ')
+ ELSEIF ( (nn_mxlice>1).AND.(nn_ice==1) ) THEN
+ CALL ctl_stop( 'zdf_gls_init: with no ice model, nn_mxlice must be 0 or 1')
+ ENDIF
WRITE(numout,*)
ENDIF
diff --git a/src/OCE/ZDF/zdftke.F90 b/src/OCE/ZDF/zdftke.F90
index 4f4ae8ac3..68dd9f6ed 100644
--- a/src/OCE/ZDF/zdftke.F90
+++ b/src/OCE/ZDF/zdftke.F90
@@ -761,6 +761,11 @@ CONTAINS
CASE DEFAULT
CALL ctl_stop( 'zdf_tke_init: wrong value for nn_mxlice, should be 0,1,2,3 or 4')
END SELECT
+ IF ( (nn_mxlice>0).AND.(nn_ice==0) ) THEN
+ CALL ctl_stop( 'zdf_tke_init: with no ice at all, nn_mxlice must be 0 ')
+ ELSEIF ( (nn_mxlice>1).AND.(nn_ice==1) ) THEN
+ CALL ctl_stop( 'zdf_tke_init: with no ice model, nn_mxlice must be 0 or 1')
+ ENDIF
ENDIF
WRITE(numout,*) ' Langmuir cells parametrization ln_lc = ', ln_lc
WRITE(numout,*) ' coef to compute vertical velocity of LC rn_lc = ', rn_lc
diff --git a/src/OCE/step.F90 b/src/OCE/step.F90
index 99e8a952d..ff2759f8c 100644
--- a/src/OCE/step.F90
+++ b/src/OCE/step.F90
@@ -189,11 +189,11 @@ CONTAINS
IF( ln_zps .OR. l_ldfslp ) CALL eos( ts(:,:,:,:,Nbb), rhd, gdept_0(:,:,:) ) ! before in situ density
IF( ln_zps .AND. .NOT. ln_isfcav) &
- & CALL zps_hde ( kstp, Nnn, jpts, ts(:,:,:,:,Nbb), gtsu, gtsv, & ! Partial steps: before horizontal gradient
+ & CALL zps_hde ( kstp, jpts, ts(:,:,:,:,Nbb), gtsu, gtsv, & ! Partial steps: before horizontal gradient
& rhd, gru , grv ) ! of t, s, rd at the last ocean level
IF( ln_zps .AND. ln_isfcav) &
- & CALL zps_hde_isf( kstp, Nnn, jpts, ts(:,:,:,:,Nbb), gtsu, gtsv, gtui, gtvi, & ! Partial steps for top cell (ISF)
+ & CALL zps_hde_isf( kstp, jpts, ts(:,:,:,:,Nbb), gtsu, gtsv, gtui, gtvi, & ! Partial steps for top cell (ISF)
& rhd, gru , grv , grui, grvi ) ! of t, s, rd at the first ocean level
IF( l_ldfslp ) THEN ! slope of lateral mixing
diff --git a/src/OCE/stpmlf.F90 b/src/OCE/stpmlf.F90
index 407e38686..bae696704 100644
--- a/src/OCE/stpmlf.F90
+++ b/src/OCE/stpmlf.F90
@@ -193,11 +193,11 @@ CONTAINS
IF( ln_zps .OR. l_ldfslp ) CALL eos( ts(:,:,:,:,Nbb), rhd, gdept_0(:,:,:) ) ! before in situ density
IF( ln_zps .AND. .NOT. ln_isfcav) &
- & CALL zps_hde ( kstp, Nnn, jpts, ts(:,:,:,:,Nbb), gtsu, gtsv, & ! Partial steps: before horizontal gradient
+ & CALL zps_hde ( kstp, jpts, ts(:,:,:,:,Nbb), gtsu, gtsv, & ! Partial steps: before horizontal gradient
& rhd, gru , grv ) ! of t, s, rd at the last ocean level
IF( ln_zps .AND. ln_isfcav) &
- & CALL zps_hde_isf( kstp, Nnn, jpts, ts(:,:,:,:,Nbb), gtsu, gtsv, gtui, gtvi, & ! Partial steps for top cell (ISF)
+ & CALL zps_hde_isf( kstp, jpts, ts(:,:,:,:,Nbb), gtsu, gtsv, gtui, gtvi, & ! Partial steps for top cell (ISF)
& rhd, gru , grv , grui, grvi ) ! of t, s, rd at the first ocean level
IF( l_ldfslp ) THEN ! slope of lateral mixing
diff --git a/src/OCE/stprk3.F90 b/src/OCE/stprk3.F90
index 5b1ca4f34..337b55f62 100644
--- a/src/OCE/stprk3.F90
+++ b/src/OCE/stprk3.F90
@@ -154,20 +154,20 @@ CONTAINS
! VERTICAL PHYSICS
!!st CALL zdf_phy( kstp, Nbb, Nnn, Nrhs ) ! vertical physics update (top/bot drag, avt, avs, avm + MLD)
CALL zdf_phy( kstp, Nbb, Nbb, Nrhs ) ! vertical physics update (top/bot drag, avt, avs, avm + MLD)
-!!gm
+!!gm gdep
! LATERAL PHYSICS
!
- IF( l_ldfslp ) THEN ! slope of lateral mixing
-!!gm gdep
- CALL eos( ts(:,:,:,:,Nbb), rhd, gdept_0(:,:,:) ) ! before in situ density
+ IF( ln_zps .OR. l_ldfslp ) CALL eos( ts(:,:,:,:,Nbb), rhd, gdept_0(:,:,:) ) ! before in situ density
- IF( ln_zps .AND. .NOT. ln_isfcav) &
- & CALL zps_hde ( kstp, Nbb, jpts, ts(:,:,:,:,Nbb), gtsu, gtsv, & ! Partial steps: before horizontal gradient
+ IF( ln_zps .AND. .NOT. ln_isfcav) &
+ & CALL zps_hde ( kstp, jpts, ts(:,:,:,:,Nbb), gtsu, gtsv, & ! Partial steps: before horizontal gradient
& rhd, gru , grv ) ! of t, s, rd at the last ocean level
- IF( ln_zps .AND. ln_isfcav) &
- & CALL zps_hde_isf( kstp, Nbb, jpts, ts(:,:,:,:,Nbb), gtsu, gtsv, gtui, gtvi, & ! Partial steps for top cell (ISF)
+ IF( ln_zps .AND. ln_isfcav) &
+ & CALL zps_hde_isf( kstp, jpts, ts(:,:,:,:,Nbb), gtsu, gtsv, gtui, gtvi, & ! Partial steps for top cell (ISF)
& rhd, gru , grv , grui, grvi ) ! of t, s, rd at the first ocean level
+
+ IF( l_ldfslp ) THEN ! slope of lateral mixing
IF( ln_traldf_triad ) THEN
CALL ldf_slp_triad( kstp, Nbb, Nbb ) ! before slope for triad operator
ELSE
diff --git a/src/OFF/dtadyn.F90 b/src/OFF/dtadyn.F90
index c55f56994..253d699db 100644
--- a/src/OFF/dtadyn.F90
+++ b/src/OFF/dtadyn.F90
@@ -675,10 +675,10 @@ CONTAINS
! Partial steps: before Horizontal DErivative
IF( ln_zps .AND. .NOT. ln_isfcav) &
- & CALL zps_hde ( kt, Kmm, jpts, pts, gtsu, gtsv, & ! Partial steps: before horizontal gradient
+ & CALL zps_hde ( kt, jpts, pts, gtsu, gtsv, & ! Partial steps: before horizontal gradient
& rhd, gru , grv ) ! of t, s, rd at the last ocean level
IF( ln_zps .AND. ln_isfcav) &
- & CALL zps_hde_isf( kt, Kmm, jpts, pts, gtsu, gtsv, gtui, gtvi, & ! Partial steps for top cell (ISF)
+ & CALL zps_hde_isf( kt, jpts, pts, gtsu, gtsv, gtui, gtvi, & ! Partial steps for top cell (ISF)
& rhd, gru , grv , grui, grvi ) ! of t, s, rd at the first ocean level
rn2b(:,:,:) = rn2(:,:,:) ! needed for zdfmxl
diff --git a/src/TOP/TRP/trctrp.F90 b/src/TOP/TRP/trctrp.F90
index 03b07bdc8..f164444ed 100644
--- a/src/TOP/TRP/trctrp.F90
+++ b/src/TOP/TRP/trctrp.F90
@@ -70,8 +70,8 @@ CONTAINS
!
! ! Partial top/bottom cell: GRADh( trb )
IF( ln_zps ) THEN
- IF( ln_isfcav ) THEN ; CALL zps_hde_isf( kt, Kmm, jptra, tr(:,:,:,:,Kbb), pgtu=gtru, pgtv=gtrv, pgtui=gtrui, pgtvi=gtrvi ) ! both top & bottom
- ELSE ; CALL zps_hde ( kt, Kmm, jptra, tr(:,:,:,:,Kbb), gtru, gtrv ) ! only bottom
+ IF( ln_isfcav ) THEN ; CALL zps_hde_isf( kt, jptra, tr(:,:,:,:,Kbb), pgtu=gtru, pgtv=gtrv, pgtui=gtrui, pgtvi=gtrvi ) ! both top & bottom
+ ELSE ; CALL zps_hde ( kt, jptra, tr(:,:,:,:,Kbb), gtru, gtrv ) ! only bottom
ENDIF
ENDIF
!
diff --git a/src/TOP/trcini.F90 b/src/TOP/trcini.F90
index 486f5208c..394f1cd16 100644
--- a/src/TOP/trcini.F90
+++ b/src/TOP/trcini.F90
@@ -234,7 +234,6 @@ CONTAINS
!! *** ROUTINE trc_ini_state ***
!! ** Purpose : Initialisation of passive tracer concentration
!!----------------------------------------------------------------------
- USE zpshde ! partial step: hor. derivative (zps_hde routine)
USE trcrst ! passive tracers restart
USE trcdta ! initialisation from files
!
--
GitLab