Newer
Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
MODULE zdfric
!!======================================================================
!! *** MODULE zdfric ***
!! Ocean physics: vertical mixing coefficient compute from the local
!! Richardson number dependent formulation
!!======================================================================
!! History : OPA ! 1987-09 (P. Andrich) Original code
!! 4.0 ! 1991-11 (G. Madec)
!! 7.0 ! 1996-01 (G. Madec) complete rewriting of multitasking suppression of common work arrays
!! 8.0 ! 1997-06 (G. Madec) complete rewriting of zdfmix
!! NEMO 1.0 ! 2002-06 (G. Madec) F90: Free form and module
!! 3.3 ! 2010-10 (C. Ethe, G. Madec) reorganisation of initialisation phase
!! 3.3.1! 2011-09 (P. Oddo) Mixed layer depth parameterization
!! 4.0 ! 2017-04 (G. Madec) remove CPP ddm key & avm at t-point only
!!----------------------------------------------------------------------
!!----------------------------------------------------------------------
!! zdf_ric_init : initialization, namelist read, & parameters control
!! zdf_ric : update momentum and tracer Kz from the Richardson number
!! ric_rst : read/write RIC restart in ocean restart file
!!----------------------------------------------------------------------
USE oce ! ocean dynamics and tracers variables
USE dom_oce ! ocean space and time domain variables
USE zdf_oce ! vertical physics: variables
USE phycst ! physical constants
USE sbc_oce, ONLY : taum
!
USE in_out_manager ! I/O manager
USE iom ! I/O manager library
USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined)
IMPLICIT NONE
PRIVATE
PUBLIC zdf_ric ! called by zdfphy.F90
PUBLIC ric_rst ! called by zdfphy.F90
PUBLIC zdf_ric_init ! called by nemogcm.F90
! !!* Namelist namzdf_ric : Richardson number dependent Kz *
INTEGER :: nn_ric ! coefficient of the parameterization
REAL(wp) :: rn_avmri ! maximum value of the vertical eddy viscosity
REAL(wp) :: rn_alp ! coefficient of the parameterization
REAL(wp) :: rn_ekmfc ! Ekman Factor Coeff
REAL(wp) :: rn_mldmin ! minimum mixed layer (ML) depth
REAL(wp) :: rn_mldmax ! maximum mixed layer depth
REAL(wp) :: rn_wtmix ! Vertical eddy Diff. in the ML
REAL(wp) :: rn_wvmix ! Vertical eddy Visc. in the ML
LOGICAL :: ln_mldw ! Use or not the MLD parameters
!! * Substitutions
# include "do_loop_substitute.h90"
# include "domzgr_substitute.h90"
!!----------------------------------------------------------------------
!! NEMO/OCE 4.0 , NEMO Consortium (2018)
!! $Id: zdfric.F90 15277 2021-09-22 13:19:40Z clem $
!! Software governed by the CeCILL license (see ./LICENSE)
!!----------------------------------------------------------------------
CONTAINS
SUBROUTINE zdf_ric_init
!!----------------------------------------------------------------------
!! *** ROUTINE zdf_ric_init ***
!!
!! ** Purpose : Initialization of the vertical eddy diffusivity and
!! viscosity coef. for the Richardson number dependent formulation.
!!
!! ** Method : Read the namzdf_ric namelist and check the parameter values
!!
!! ** input : Namelist namzdf_ric
!!
!! ** Action : increase by 1 the nstop flag is setting problem encounter
!!----------------------------------------------------------------------
INTEGER :: ji, jj, jk ! dummy loop indices
INTEGER :: ios ! Local integer output status for namelist read
!!
NAMELIST/namzdf_ric/ rn_avmri, rn_alp , nn_ric , rn_ekmfc, &
& rn_mldmin, rn_mldmax, rn_wtmix, rn_wvmix, ln_mldw
!!----------------------------------------------------------------------
!
READ ( numnam_ref, namzdf_ric, IOSTAT = ios, ERR = 901)
901 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namzdf_ric in reference namelist' )
READ ( numnam_cfg, namzdf_ric, IOSTAT = ios, ERR = 902 )
902 IF( ios > 0 ) CALL ctl_nam ( ios , 'namzdf_ric in configuration namelist' )
IF(lwm) WRITE ( numond, namzdf_ric )
!
IF(lwp) THEN ! Control print
WRITE(numout,*)
WRITE(numout,*) 'zdf_ric_init : Ri depend vertical mixing scheme'
WRITE(numout,*) '~~~~~~~~~~~~'
WRITE(numout,*) ' Namelist namzdf_ric : set Kz=F(Ri) parameters'
WRITE(numout,*) ' maximum vertical viscosity rn_avmri = ', rn_avmri
WRITE(numout,*) ' coefficient rn_alp = ', rn_alp
WRITE(numout,*) ' exponent nn_ric = ', nn_ric
WRITE(numout,*) ' Ekman layer enhanced mixing ln_mldw = ', ln_mldw
WRITE(numout,*) ' Ekman Factor Coeff rn_ekmfc = ', rn_ekmfc
WRITE(numout,*) ' minimum mixed layer depth rn_mldmin = ', rn_mldmin
WRITE(numout,*) ' maximum mixed layer depth rn_mldmax = ', rn_mldmax
WRITE(numout,*) ' Vertical eddy Diff. in the ML rn_wtmix = ', rn_wtmix
WRITE(numout,*) ' Vertical eddy Visc. in the ML rn_wvmix = ', rn_wvmix
ENDIF
!
CALL ric_rst( nit000, 'READ' ) !* read or initialize all required files
!
END SUBROUTINE zdf_ric_init
SUBROUTINE zdf_ric( kt, Kbb, Kmm, p_avm, p_avt )
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
!!----------------------------------------------------------------------
!! *** ROUTINE zdfric ***
!!
!! ** Purpose : Compute the before eddy viscosity and diffusivity as
!! a function of the local richardson number.
!!
!! ** Method : Local richardson number dependent formulation of the
!! vertical eddy viscosity and diffusivity coefficients.
!! The eddy coefficients are given by:
!! avm = avm0 + avmb
!! avt = avm0 / (1 + rn_alp*ri)
!! with ri = N^2 / dz(u)**2
!! = e3w**2 * rn2/[ mi( dk(uu(:,:,:,Kbb)) )+mj( dk(vv(:,:,:,Kbb)) ) ]
!! avm0= rn_avmri / (1 + rn_alp*Ri)**nn_ric
!! where ri is the before local Richardson number,
!! rn_avmri is the maximum value reaches by avm and avt
!! and rn_alp, nn_ric are adjustable parameters.
!! Typical values : rn_alp=5. and nn_ric=2.
!!
!! As second step compute Ekman depth from wind stress forcing
!! and apply namelist provided vertical coeff within this depth.
!! The Ekman depth is:
!! Ustar = SQRT(Taum/rho0)
!! ekd= rn_ekmfc * Ustar / f0
!! Large et al. (1994, eq.24) suggest rn_ekmfc=0.7; however, the derivation
!! of the above equation indicates the value is somewhat arbitrary; therefore
!! we allow the freedom to increase or decrease this value, if the
!! Ekman depth estimate appears too shallow or too deep, respectively.
!! Ekd is then limited by rn_mldmin and rn_mldmax provided in the
!! namelist
!! N.B. the mask are required for implicit scheme, and surface
!! and bottom value already set in zdfphy.F90
!!
!! ** Action : avm, avt mixing coeff (inner domain values only)
!!
!! References : Pacanowski & Philander 1981, JPO, 1441-1451.
!! PFJ Lermusiaux 2001.
!!----------------------------------------------------------------------
INTEGER , INTENT(in ) :: kt ! ocean time-step
INTEGER , INTENT(in ) :: Kmm , Kbb ! ocean time level index
REAL(wp), DIMENSION(:,:,:) , INTENT(inout) :: p_avm, p_avt ! momentum and tracer Kz (w-points)
!!
INTEGER :: ji, jj, jk ! dummy loop indices
REAL(wp) :: zcfRi, zav, zustar, zhek, zdku, zdkv, zwx ! local scalars
REAL(wp), DIMENSION(A2D(nn_hls)) :: zh_ekm ! 2D workspace
!!----------------------------------------------------------------------
!
! !== avm and avt = F(Richardson number) ==!
DO_3D_OVR( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 2, jpkm1 ) ! coefficient = F(richardson number) (avm-weighted Ri)
zdku = 0.5 / e3uw(ji,jj,jk,Kbb) * ( uu(ji-1,jj,jk-1,Kbb) + uu(ji,jj,jk-1,Kbb) &
& - uu(ji-1,jj,jk,Kbb) - uu(ji,jj,jk,Kbb) ) * wumask(ji,jj,jk)
zdkv = 0.5 / e3vw(ji,jj,jk,Kbb) * ( vv(ji,jj-1,jk-1,Kbb) + vv(ji,jj,jk-1,Kbb) &
& - vv(ji,jj-1,jk,Kbb) - vv(ji,jj,jk,Kbb) ) * wvmask(ji,jj,jk)
zwx = zdku*zdku + zdkv*zdkv
zcfRi = 1._wp / ( 1._wp + rn_alp * MAX( 0._wp , rn2(ji,jj,jk) / ( zwx + 1.e-20 ) ) )
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
zav = rn_avmri * zcfRi**nn_ric
! ! avm and avt coefficients
p_avm(ji,jj,jk) = MAX( zav , avmb(jk) ) * wmask(ji,jj,jk)
p_avt(ji,jj,jk) = MAX( zav * zcfRi , avtb(jk) ) * wmask(ji,jj,jk)
END_3D
!
!!gm BUG <<<<==== This param can't work at low latitude
!!gm it provides there much to thick mixed layer ( summer 150m in GYRE configuration !!! )
!
IF( ln_mldw ) THEN !== set a minimum value in the Ekman layer ==!
!
DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 )
zustar = SQRT( taum(ji,jj) * r1_rho0 )
zhek = rn_ekmfc * zustar / ( ABS( ff_t(ji,jj) ) + rsmall ) ! Ekman depth
zh_ekm(ji,jj) = MAX( rn_mldmin , MIN( zhek , rn_mldmax ) ) ! set allowed range
END_2D
DO_3D_OVR( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 2, jpkm1 ) !* minimum mixing coeff. within the Ekman layer
IF( gdept(ji,jj,jk,Kmm) < zh_ekm(ji,jj) ) THEN
p_avm(ji,jj,jk) = MAX( p_avm(ji,jj,jk), rn_wvmix ) * wmask(ji,jj,jk)
p_avt(ji,jj,jk) = MAX( p_avt(ji,jj,jk), rn_wtmix ) * wmask(ji,jj,jk)
ENDIF
END_3D
ENDIF
!
END SUBROUTINE zdf_ric
SUBROUTINE ric_rst( kt, cdrw )
!!---------------------------------------------------------------------
!! *** ROUTINE ric_rst ***
!!
!! ** Purpose : Read or write TKE file (en) in restart file
!!
!! ** Method : use of IOM library
!! if the restart does not contain TKE, en is either
!! set to rn_emin or recomputed
!!----------------------------------------------------------------------
INTEGER , INTENT(in) :: kt ! ocean time-step
CHARACTER(len=*), INTENT(in) :: cdrw ! "READ"/"WRITE" flag
!
INTEGER :: jit, jk ! dummy loop indices
INTEGER :: id1, id2 ! local integers
!!----------------------------------------------------------------------
!
IF( TRIM(cdrw) == 'READ' ) THEN ! Read/initialise
! ! ---------------
! !* Read the restart file
IF( ln_rstart ) THEN
id1 = iom_varid( numror, 'avt_k', ldstop = .FALSE. )
id2 = iom_varid( numror, 'avm_k', ldstop = .FALSE. )
!
IF( MIN( id1, id2 ) > 0 ) THEN ! restart exists => read it
CALL iom_get( numror, jpdom_auto, 'avt_k', avt_k )
CALL iom_get( numror, jpdom_auto, 'avm_k', avm_k )
ENDIF
ENDIF
! !* otherwise Kz already set to the background value in zdf_phy_init
!
ELSEIF( TRIM(cdrw) == 'WRITE' ) THEN ! Create restart file
! ! -------------------
IF(lwp) WRITE(numout,*) '---- ric-rst ----'
CALL iom_rstput( kt, nitrst, numrow, 'avt_k', avt_k )
CALL iom_rstput( kt, nitrst, numrow, 'avm_k', avm_k)
!
ENDIF
!
END SUBROUTINE ric_rst
!!======================================================================
END MODULE zdfric