- Timestamp:
- 2020-09-29T12:41:06+02:00 (4 years ago)
- Location:
- NEMO/branches/2020/r12377_ticket2386
- Files:
-
- 2 edited
-
. (modified) (1 prop)
-
src/OCE/ZDF/zdftke.F90 (modified) (26 diffs)
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2020/r12377_ticket2386
- Property svn:externals
-
old new 3 3 ^/utils/build/mk@HEAD mk 4 4 ^/utils/tools@HEAD tools 5 ^/vendors/AGRIF/dev @HEADext/AGRIF5 ^/vendors/AGRIF/dev_r12970_AGRIF_CMEMS ext/AGRIF 6 6 ^/vendors/FCM@HEAD ext/FCM 7 7 ^/vendors/IOIPSL@HEAD ext/IOIPSL 8 8 9 9 # SETTE 10 ^/utils/CI/sette@ HEADsette10 ^/utils/CI/sette@13507 sette
-
- Property svn:externals
-
NEMO/branches/2020/r12377_ticket2386/src/OCE/ZDF/zdftke.F90
r12511 r13540 28 28 !! 3.6 ! 2014-11 (P. Mathiot) add ice shelf capability 29 29 !! 4.0 ! 2017-04 (G. Madec) remove CPP ddm key & avm at t-point only 30 !! - ! 2017-05 (G. Madec) add top/bottom friction as boundary condition (ln_drg)30 !! - ! 2017-05 (G. Madec) add top/bottom friction as boundary condition 31 31 !!---------------------------------------------------------------------- 32 32 … … 45 45 USE zdfdrg ! vertical physics: top/bottom drag coef. 46 46 USE zdfmxl ! vertical physics: mixed layer 47 #if defined key_si3 48 USE ice, ONLY: hm_i, h_i 49 #endif 50 #if defined key_cice 51 USE sbc_ice, ONLY: h_i 52 #endif 47 53 ! 48 54 USE in_out_manager ! I/O manager … … 62 68 ! !!** Namelist namzdf_tke ** 63 69 LOGICAL :: ln_mxl0 ! mixing length scale surface value as function of wind stress or not 70 INTEGER :: nn_mxlice ! type of scaling under sea-ice (=0/1/2/3) 71 REAL(wp) :: rn_mxlice ! ice thickness value when scaling under sea-ice 64 72 INTEGER :: nn_mxl ! type of mixing length (=0/1/2/3) 65 73 REAL(wp) :: rn_mxl0 ! surface min value of mixing length (kappa*z_o=0.4*0.1 m) [m] … … 71 79 REAL(wp) :: rn_emin0 ! surface minimum value of tke [m2/s2] 72 80 REAL(wp) :: rn_bshear ! background shear (>0) currently a numerical threshold (do not change it) 73 LOGICAL :: ln_drg ! top/bottom friction forcing flag74 81 INTEGER :: nn_etau ! type of depth penetration of surface tke (=0/1/2/3) 75 82 INTEGER :: nn_htau ! type of tke profile of penetration (=0/1) 76 83 REAL(wp) :: rn_efr ! fraction of TKE surface value which penetrates in the ocean 77 REAL(wp) :: rn_eice ! =0 ON below sea-ice, =4 OFF when ice fraction > 1/478 84 LOGICAL :: ln_lc ! Langmuir cells (LC) as a source term of TKE or not 79 85 REAL(wp) :: rn_lc ! coef to compute vertical velocity of Langmuir cells 86 INTEGER :: nn_eice ! attenutaion of langmuir & surface wave breaking under ice (=0/1/2/3) 80 87 81 88 REAL(wp) :: ri_cri ! critic Richardson number (deduced from rn_ediff and rn_ediss values) … … 90 97 !! * Substitutions 91 98 # include "do_loop_substitute.h90" 99 # include "domzgr_substitute.h90" 92 100 !!---------------------------------------------------------------------- 93 101 !! NEMO/OCE 4.0 , NEMO Consortium (2018) … … 191 199 REAL(wp), DIMENSION(:,:,:) , INTENT(in ) :: p_avm, p_avt ! vertical eddy viscosity & diffusivity (w-points) 192 200 ! 193 INTEGER :: ji, jj, jk ! dummy loop arguments201 INTEGER :: ji, jj, jk ! dummy loop arguments 194 202 REAL(wp) :: zetop, zebot, zmsku, zmskv ! local scalars 195 203 REAL(wp) :: zrhoa = 1.22 ! Air density kg/m3 196 204 REAL(wp) :: zcdrag = 1.5e-3 ! drag coefficient 197 REAL(wp) :: zbbrau, z ri! local scalars198 REAL(wp) :: zfact1, zfact2, zfact3 ! - 199 REAL(wp) :: ztx2 , zty2 , zcof ! - 200 REAL(wp) :: ztau , zdif ! - 201 REAL(wp) :: zus , zwlc , zind ! - 202 REAL(wp) :: zzd_up, zzd_lw ! - 205 REAL(wp) :: zbbrau, zbbirau, zri ! local scalars 206 REAL(wp) :: zfact1, zfact2, zfact3 ! - - 207 REAL(wp) :: ztx2 , zty2 , zcof ! - - 208 REAL(wp) :: ztau , zdif ! - - 209 REAL(wp) :: zus , zwlc , zind ! - - 210 REAL(wp) :: zzd_up, zzd_lw ! - - 203 211 INTEGER , DIMENSION(jpi,jpj) :: imlc 204 REAL(wp), DIMENSION(jpi,jpj) :: z hlc, zfr_i212 REAL(wp), DIMENSION(jpi,jpj) :: zice_fra, zhlc, zus3 205 213 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zpelc, zdiag, zd_up, zd_lw 206 214 !!-------------------------------------------------------------------- 207 215 ! 208 zbbrau = rn_ebb / rho0 ! Local constant initialisation 209 zfact1 = -.5_wp * rn_Dt 210 zfact2 = 1.5_wp * rn_Dt * rn_ediss 211 zfact3 = 0.5_wp * rn_ediss 216 zbbrau = rn_ebb / rho0 ! Local constant initialisation 217 zbbirau = 3.75_wp / rho0 218 zfact1 = -.5_wp * rn_Dt 219 zfact2 = 1.5_wp * rn_Dt * rn_ediss 220 zfact3 = 0.5_wp * rn_ediss 221 ! 222 ! ice fraction considered for attenuation of langmuir & wave breaking 223 SELECT CASE ( nn_eice ) 224 CASE( 0 ) ; zice_fra(:,:) = 0._wp 225 CASE( 1 ) ; zice_fra(:,:) = TANH( fr_i(:,:) * 10._wp ) 226 CASE( 2 ) ; zice_fra(:,:) = fr_i(:,:) 227 CASE( 3 ) ; zice_fra(:,:) = MIN( 4._wp * fr_i(:,:) , 1._wp ) 228 END SELECT 212 229 ! 213 230 ! !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 214 231 ! ! Surface/top/bottom boundary condition on tke 215 232 ! !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 216 217 DO_2D_00_00 233 ! 234 DO_2D( 0, 0, 0, 0 ) ! en(1) = rn_ebb taum / rau0 (min value rn_emin0) 235 !! clem: this should be the right formulation but it makes the model unstable unless drags are calculated implicitly 236 !! one way around would be to increase zbbirau 237 !! en(ji,jj,1) = MAX( rn_emin0, ( ( 1._wp - fr_i(ji,jj) ) * zbbrau + & 238 !! & fr_i(ji,jj) * zbbirau ) * taum(ji,jj) ) * tmask(ji,jj,1) 218 239 en(ji,jj,1) = MAX( rn_emin0, zbbrau * taum(ji,jj) ) * tmask(ji,jj,1) 219 240 END_2D 220 IF ( ln_isfcav ) THEN221 DO_2D_00_00222 en(ji,jj,mikt(ji,jj)) = rn_emin * tmask(ji,jj,1)223 END_2D224 ENDIF225 241 ! 226 242 ! !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< … … 232 248 ! Note that stress averaged is done using an wet-only calculation of u and v at t-point like in zdfsh2 233 249 ! 234 IF( ln_drg ) THEN!== friction used as top/bottom boundary condition on TKE235 ! 236 DO_2D _00_00250 IF( .NOT.ln_drg_OFF ) THEN !== friction used as top/bottom boundary condition on TKE 251 ! 252 DO_2D( 0, 0, 0, 0 ) ! bottom friction 237 253 zmsku = ( 2. - umask(ji-1,jj,mbkt(ji,jj)) * umask(ji,jj,mbkt(ji,jj)) ) 238 254 zmskv = ( 2. - vmask(ji,jj-1,mbkt(ji,jj)) * vmask(ji,jj,mbkt(ji,jj)) ) … … 242 258 en(ji,jj,mbkt(ji,jj)+1) = MAX( zebot, rn_emin ) * ssmask(ji,jj) 243 259 END_2D 244 IF( ln_isfcav ) THEN ! top friction245 DO_2D _00_00260 IF( ln_isfcav ) THEN 261 DO_2D( 0, 0, 0, 0 ) ! top friction 246 262 zmsku = ( 2. - umask(ji-1,jj,mikt(ji,jj)) * umask(ji,jj,mikt(ji,jj)) ) 247 263 zmskv = ( 2. - vmask(ji,jj-1,mikt(ji,jj)) * vmask(ji,jj,mikt(ji,jj)) ) … … 249 265 zetop = - 0.001875_wp * rCdU_top(ji,jj) * SQRT( ( zmsku*( uu(ji,jj,mikt(ji,jj),Kbb)+uu(ji-1,jj,mikt(ji,jj),Kbb) ) )**2 & 250 266 & + ( zmskv*( vv(ji,jj,mikt(ji,jj),Kbb)+vv(ji,jj-1,mikt(ji,jj),Kbb) ) )**2 ) 251 en(ji,jj,mikt(ji,jj)) = MAX( zetop, rn_emin ) * (1._wp - tmask(ji,jj,1)) ! masked at ocean surface 267 ! (1._wp - tmask(ji,jj,1)) * ssmask(ji,jj) = 1 where ice shelves are present 268 en(ji,jj,mikt(ji,jj)) = en(ji,jj,1) * tmask(ji,jj,1) & 269 & + MAX( zetop, rn_emin ) * (1._wp - tmask(ji,jj,1)) * ssmask(ji,jj) 252 270 END_2D 253 271 ENDIF … … 262 280 zpelc(:,:,1) = MAX( rn2b(:,:,1), 0._wp ) * gdepw(:,:,1,Kmm) * e3w(:,:,1,Kmm) 263 281 DO jk = 2, jpk 264 zpelc(:,:,jk) = zpelc(:,:,jk-1) + MAX( rn2b(:,:,jk), 0._wp ) * gdepw(:,:,jk,Kmm) * e3w(:,:,jk,Kmm) 282 zpelc(:,:,jk) = zpelc(:,:,jk-1) + & 283 & MAX( rn2b(:,:,jk), 0._wp ) * gdepw(:,:,jk,Kmm) * e3w(:,:,jk,Kmm) 265 284 END DO 266 285 ! !* finite Langmuir Circulation depth 267 286 zcof = 0.5 * 0.016 * 0.016 / ( zrhoa * zcdrag ) 268 287 imlc(:,:) = mbkt(:,:) + 1 ! Initialization to the number of w ocean point (=2 over land) 269 DO_3DS _11_11( jpkm1, 2, -1 )270 zus = zcof * taum(ji,jj)288 DO_3DS( 1, 1, 1, 1, jpkm1, 2, -1 ) ! Last w-level at which zpelc>=0.5*us*us 289 zus = zcof * taum(ji,jj) ! with us=0.016*wind(starting from jpk-1) 271 290 IF( zpelc(ji,jj,jk) > zus ) imlc(ji,jj) = jk 272 291 END_3D 273 292 ! ! finite LC depth 274 DO_2D _11_11293 DO_2D( 1, 1, 1, 1 ) 275 294 zhlc(ji,jj) = gdepw(ji,jj,imlc(ji,jj),Kmm) 276 295 END_2D 277 296 zcof = 0.016 / SQRT( zrhoa * zcdrag ) 278 DO_2D _00_00297 DO_2D( 0, 0, 0, 0 ) 279 298 zus = zcof * SQRT( taum(ji,jj) ) ! Stokes drift 280 zfr_i(ji,jj) = ( 1._wp - 4._wp * fr_i(ji,jj) ) * zus * zus * zus * tmask(ji,jj,1) ! zus > 0. ok 281 IF (zfr_i(ji,jj) < 0. ) zfr_i(ji,jj) = 0. 299 zus3(ji,jj) = MAX( 0._wp, 1._wp - zice_fra(ji,jj) ) * zus * zus * zus * tmask(ji,jj,1) ! zus > 0. ok 282 300 END_2D 283 DO_3D_00_00( 2, jpkm1 ) 284 IF ( zfr_i(ji,jj) /= 0. ) THEN 285 ! vertical velocity due to LC 301 DO_3D( 0, 0, 0, 0, 2, jpkm1 ) !* TKE Langmuir circulation source term added to en 302 IF ( zus3(ji,jj) /= 0._wp ) THEN 286 303 IF ( gdepw(ji,jj,jk,Kmm) - zhlc(ji,jj) < 0 .AND. wmask(ji,jj,jk) /= 0. ) THEN 287 304 ! ! vertical velocity due to LC 288 zwlc = rn_lc * SIN( rpi * gdepw(ji,jj,jk,Kmm) / zhlc(ji,jj) ) ! warning: optimization: zus^3 is in zfr_i305 zwlc = rn_lc * SIN( rpi * gdepw(ji,jj,jk,Kmm) / zhlc(ji,jj) ) 289 306 ! ! TKE Langmuir circulation source term 290 en(ji,jj,jk) = en(ji,jj,jk) + rn_Dt * z fr_i(ji,jj) * ( zwlc * zwlc * zwlc ) / zhlc(ji,jj)307 en(ji,jj,jk) = en(ji,jj,jk) + rn_Dt * zus3(ji,jj) * ( zwlc * zwlc * zwlc ) / zhlc(ji,jj) 291 308 ENDIF 292 309 ENDIF … … 302 319 ! ! zdiag : diagonal zd_up : upper diagonal zd_lw : lower diagonal 303 320 ! 304 IF( nn_pdl == 1 ) THEN !* Prandtl number = F( Ri )305 DO_3D _00_00(2, jpkm1 )321 IF( nn_pdl == 1 ) THEN !* Prandtl number = F( Ri ) 322 DO_3D( 0, 0, 0, 0, 2, jpkm1 ) 306 323 ! ! local Richardson number 307 zri = MAX( rn2b(ji,jj,jk), 0._wp ) * p_avm(ji,jj,jk) / ( p_sh2(ji,jj,jk) + rn_bshear ) 324 IF (rn2b(ji,jj,jk) <= 0.0_wp) then 325 zri = 0.0_wp 326 ELSE 327 zri = rn2b(ji,jj,jk) * p_avm(ji,jj,jk) / ( p_sh2(ji,jj,jk) + rn_bshear ) 328 ENDIF 308 329 ! ! inverse of Prandtl number 309 330 apdlr(ji,jj,jk) = MAX( 0.1_wp, ri_cri / MAX( ri_cri , zri ) ) … … 311 332 ENDIF 312 333 ! 313 DO_3D _00_00( 2, jpkm1 )334 DO_3D( 0, 0, 0, 0, 2, jpkm1 ) !* Matrix and right hand side in en 314 335 zcof = zfact1 * tmask(ji,jj,jk) 315 336 ! ! A minimum of 2.e-5 m2/s is imposed on TKE vertical … … 331 352 END_3D 332 353 ! !* Matrix inversion from level 2 (tke prescribed at level 1) 333 DO_3D _00_00( 3, jpkm1 )354 DO_3D( 0, 0, 0, 0, 3, jpkm1 ) ! First recurrence : Dk = Dk - Lk * Uk-1 / Dk-1 334 355 zdiag(ji,jj,jk) = zdiag(ji,jj,jk) - zd_lw(ji,jj,jk) * zd_up(ji,jj,jk-1) / zdiag(ji,jj,jk-1) 335 356 END_3D 336 DO_2D _00_00357 DO_2D( 0, 0, 0, 0 ) ! Second recurrence : Lk = RHSk - Lk / Dk-1 * Lk-1 337 358 zd_lw(ji,jj,2) = en(ji,jj,2) - zd_lw(ji,jj,2) * en(ji,jj,1) ! Surface boudary conditions on tke 338 359 END_2D 339 DO_3D _00_00(3, jpkm1 )360 DO_3D( 0, 0, 0, 0, 3, jpkm1 ) 340 361 zd_lw(ji,jj,jk) = en(ji,jj,jk) - zd_lw(ji,jj,jk) / zdiag(ji,jj,jk-1) *zd_lw(ji,jj,jk-1) 341 362 END_3D 342 DO_2D _00_00363 DO_2D( 0, 0, 0, 0 ) ! thrid recurrence : Ek = ( Lk - Uk * Ek+1 ) / Dk 343 364 en(ji,jj,jpkm1) = zd_lw(ji,jj,jpkm1) / zdiag(ji,jj,jpkm1) 344 365 END_2D 345 DO_3DS _00_00(jpk-2, 2, -1 )366 DO_3DS( 0, 0, 0, 0, jpk-2, 2, -1 ) 346 367 en(ji,jj,jk) = ( zd_lw(ji,jj,jk) - zd_up(ji,jj,jk) * en(ji,jj,jk+1) ) / zdiag(ji,jj,jk) 347 368 END_3D 348 DO_3D _00_00( 2, jpkm1 )369 DO_3D( 0, 0, 0, 0, 2, jpkm1 ) ! set the minimum value of tke 349 370 en(ji,jj,jk) = MAX( en(ji,jj,jk), rn_emin ) * wmask(ji,jj,jk) 350 371 END_3D … … 355 376 !!gm BUG : in the exp remove the depth of ssh !!! 356 377 !!gm i.e. use gde3w in argument (gdepw(:,:,:,Kmm)) 357 358 378 ! 379 ! penetration is partly switched off below sea-ice if nn_eice/=0 380 ! 359 381 IF( nn_etau == 1 ) THEN !* penetration below the mixed layer (rn_efr fraction) 360 DO_3D _00_00( 2, jpkm1 )382 DO_3D( 0, 0, 0, 0, 2, jpkm1 ) 361 383 en(ji,jj,jk) = en(ji,jj,jk) + rn_efr * en(ji,jj,1) * EXP( -gdepw(ji,jj,jk,Kmm) / htau(ji,jj) ) & 362 & * MAX( 0.,1._wp - rn_eice *fr_i(ji,jj) )* wmask(ji,jj,jk) * tmask(ji,jj,1)384 & * MAX( 0._wp, 1._wp - zice_fra(ji,jj) ) * wmask(ji,jj,jk) * tmask(ji,jj,1) 363 385 END_3D 364 386 ELSEIF( nn_etau == 2 ) THEN !* act only at the base of the mixed layer (jk=nmln) (rn_efr fraction) 365 DO_2D _00_00387 DO_2D( 0, 0, 0, 0 ) 366 388 jk = nmln(ji,jj) 367 389 en(ji,jj,jk) = en(ji,jj,jk) + rn_efr * en(ji,jj,1) * EXP( -gdepw(ji,jj,jk,Kmm) / htau(ji,jj) ) & 368 & * MAX( 0.,1._wp - rn_eice *fr_i(ji,jj) )* wmask(ji,jj,jk) * tmask(ji,jj,1)390 & * MAX( 0._wp, 1._wp - zice_fra(ji,jj) ) * wmask(ji,jj,jk) * tmask(ji,jj,1) 369 391 END_2D 370 392 ELSEIF( nn_etau == 3 ) THEN !* penetration belox the mixed layer (HF variability) 371 DO_3D _00_00(2, jpkm1 )393 DO_3D( 0, 0, 0, 0, 2, jpkm1 ) 372 394 ztx2 = utau(ji-1,jj ) + utau(ji,jj) 373 395 zty2 = vtau(ji ,jj-1) + vtau(ji,jj) … … 376 398 zdif = rhftau_scl * MAX( 0._wp, zdif + rhftau_add ) ! apply some modifications... 377 399 en(ji,jj,jk) = en(ji,jj,jk) + zbbrau * zdif * EXP( -gdepw(ji,jj,jk,Kmm) / htau(ji,jj) ) & 378 & * MAX( 0.,1._wp - rn_eice *fr_i(ji,jj) ) * wmask(ji,jj,jk) * tmask(ji,jj,1)400 & * MAX( 0._wp, 1._wp - zice_fra(ji,jj) ) * wmask(ji,jj,jk) * tmask(ji,jj,1) 379 401 END_3D 380 402 ENDIF … … 425 447 REAL(wp) :: zrn2, zraug, zcoef, zav ! local scalars 426 448 REAL(wp) :: zdku, zdkv, zsqen ! - - 427 REAL(wp) :: zemxl, zemlm, zemlp ! - -449 REAL(wp) :: zemxl, zemlm, zemlp, zmaxice ! - - 428 450 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zmxlm, zmxld ! 3D workspace 429 451 !!-------------------------------------------------------------------- … … 438 460 zmxlm(:,:,:) = rmxl_min 439 461 zmxld(:,:,:) = rmxl_min 440 ! 441 IF( ln_mxl0 ) THEN ! surface mixing length = F(stress) : l=vkarmn*2.e5*taum/(rho0*g) 462 ! 463 IF( ln_mxl0 ) THEN ! surface mixing length = F(stress) : l=vkarmn*2.e5*taum/(rho0*g) 464 ! 442 465 zraug = vkarmn * 2.e5_wp / ( rho0 * grav ) 443 DO_2D_00_00 444 zmxlm(ji,jj,1) = MAX( rn_mxl0, zraug * taum(ji,jj) * tmask(ji,jj,1) ) 466 #if ! defined key_si3 && ! defined key_cice 467 DO_2D( 0, 0, 0, 0 ) ! No sea-ice 468 zmxlm(ji,jj,1) = zraug * taum(ji,jj) * tmask(ji,jj,1) 445 469 END_2D 446 ELSE 470 #else 471 SELECT CASE( nn_mxlice ) ! Type of scaling under sea-ice 472 ! 473 CASE( 0 ) ! No scaling under sea-ice 474 DO_2D( 0, 0, 0, 0 ) 475 zmxlm(ji,jj,1) = zraug * taum(ji,jj) * tmask(ji,jj,1) 476 END_2D 477 ! 478 CASE( 1 ) ! scaling with constant sea-ice thickness 479 DO_2D( 0, 0, 0, 0 ) 480 zmxlm(ji,jj,1) = ( ( 1._wp - fr_i(ji,jj) ) * zraug * taum(ji,jj) + & 481 & fr_i(ji,jj) * rn_mxlice ) * tmask(ji,jj,1) 482 END_2D 483 ! 484 CASE( 2 ) ! scaling with mean sea-ice thickness 485 DO_2D( 0, 0, 0, 0 ) 486 #if defined key_si3 487 zmxlm(ji,jj,1) = ( ( 1._wp - fr_i(ji,jj) ) * zraug * taum(ji,jj) + & 488 & fr_i(ji,jj) * hm_i(ji,jj) * 2._wp ) * tmask(ji,jj,1) 489 #elif defined key_cice 490 zmaxice = MAXVAL( h_i(ji,jj,:) ) 491 zmxlm(ji,jj,1) = ( ( 1._wp - fr_i(ji,jj) ) * zraug * taum(ji,jj) + & 492 & fr_i(ji,jj) * zmaxice ) * tmask(ji,jj,1) 493 #endif 494 END_2D 495 ! 496 CASE( 3 ) ! scaling with max sea-ice thickness 497 DO_2D( 0, 0, 0, 0 ) 498 zmaxice = MAXVAL( h_i(ji,jj,:) ) 499 zmxlm(ji,jj,1) = ( ( 1._wp - fr_i(ji,jj) ) * zraug * taum(ji,jj) + & 500 & fr_i(ji,jj) * zmaxice ) * tmask(ji,jj,1) 501 END_2D 502 ! 503 END SELECT 504 #endif 505 ! 506 DO_2D( 0, 0, 0, 0 ) 507 zmxlm(ji,jj,1) = MAX( rn_mxl0, zmxlm(ji,jj,1) ) 508 END_2D 509 ! 510 ELSE 447 511 zmxlm(:,:,1) = rn_mxl0 448 512 ENDIF 449 ! 450 DO_3D_00_00( 2, jpkm1 ) 513 514 ! 515 DO_3D( 0, 0, 0, 0, 2, jpkm1 ) 451 516 zrn2 = MAX( rn2(ji,jj,jk), rsmall ) 452 517 zmxlm(ji,jj,jk) = MAX( rmxl_min, SQRT( 2._wp * en(ji,jj,jk) / zrn2 ) ) … … 463 528 ! where wmask = 0 set zmxlm == e3w(:,:,:,Kmm) 464 529 CASE ( 0 ) ! bounded by the distance to surface and bottom 465 DO_3D _00_00(2, jpkm1 )530 DO_3D( 0, 0, 0, 0, 2, jpkm1 ) 466 531 zemxl = MIN( gdepw(ji,jj,jk,Kmm) - gdepw(ji,jj,mikt(ji,jj),Kmm), zmxlm(ji,jj,jk), & 467 532 & gdepw(ji,jj,mbkt(ji,jj)+1,Kmm) - gdepw(ji,jj,jk,Kmm) ) 468 533 ! wmask prevent zmxlm = 0 if jk = mikt(ji,jj) 469 zmxlm(ji,jj,jk) = zemxl * wmask(ji,jj,jk) + MIN( zmxlm(ji,jj,jk) , e3w(ji,jj,jk,Kmm) ) * (1 - wmask(ji,jj,jk)) 470 zmxld(ji,jj,jk) = zemxl * wmask(ji,jj,jk) + MIN( zmxlm(ji,jj,jk) , e3w(ji,jj,jk,Kmm) ) * (1 - wmask(ji,jj,jk)) 534 zmxlm(ji,jj,jk) = zemxl * wmask(ji,jj,jk) & 535 & + MIN( zmxlm(ji,jj,jk) , e3w(ji,jj,jk,Kmm) ) * (1 - wmask(ji,jj,jk)) 536 zmxld(ji,jj,jk) = zemxl * wmask(ji,jj,jk) & 537 & + MIN( zmxlm(ji,jj,jk) , e3w(ji,jj,jk,Kmm) ) * (1 - wmask(ji,jj,jk)) 471 538 END_3D 472 539 ! 473 540 CASE ( 1 ) ! bounded by the vertical scale factor 474 DO_3D _00_00(2, jpkm1 )541 DO_3D( 0, 0, 0, 0, 2, jpkm1 ) 475 542 zemxl = MIN( e3w(ji,jj,jk,Kmm), zmxlm(ji,jj,jk) ) 476 543 zmxlm(ji,jj,jk) = zemxl … … 479 546 ! 480 547 CASE ( 2 ) ! |dk[xml]| bounded by e3t : 481 DO_3D_00_00( 2, jpkm1 ) 482 zmxlm(ji,jj,jk) = MIN( zmxlm(ji,jj,jk-1) + e3t(ji,jj,jk-1,Kmm), zmxlm(ji,jj,jk) ) 483 END_3D 484 DO_3DS_00_00( jpkm1, 2, -1 ) 548 DO_3D( 0, 0, 0, 0, 2, jpkm1 ) ! from the surface to the bottom : 549 zmxlm(ji,jj,jk) = & 550 & MIN( zmxlm(ji,jj,jk-1) + e3t(ji,jj,jk-1,Kmm), zmxlm(ji,jj,jk) ) 551 END_3D 552 DO_3DS( 0, 0, 0, 0, jpkm1, 2, -1 ) ! from the bottom to the surface : 485 553 zemxl = MIN( zmxlm(ji,jj,jk+1) + e3t(ji,jj,jk+1,Kmm), zmxlm(ji,jj,jk) ) 486 554 zmxlm(ji,jj,jk) = zemxl … … 489 557 ! 490 558 CASE ( 3 ) ! lup and ldown, |dk[xml]| bounded by e3t : 491 DO_3D_00_00( 2, jpkm1 ) 492 zmxld(ji,jj,jk) = MIN( zmxld(ji,jj,jk-1) + e3t(ji,jj,jk-1,Kmm), zmxlm(ji,jj,jk) ) 493 END_3D 494 DO_3DS_00_00( jpkm1, 2, -1 ) 495 zmxlm(ji,jj,jk) = MIN( zmxlm(ji,jj,jk+1) + e3t(ji,jj,jk+1,Kmm), zmxlm(ji,jj,jk) ) 496 END_3D 497 DO_3D_00_00( 2, jpkm1 ) 559 DO_3D( 0, 0, 0, 0, 2, jpkm1 ) ! from the surface to the bottom : lup 560 zmxld(ji,jj,jk) = & 561 & MIN( zmxld(ji,jj,jk-1) + e3t(ji,jj,jk-1,Kmm), zmxlm(ji,jj,jk) ) 562 END_3D 563 DO_3DS( 0, 0, 0, 0, jpkm1, 2, -1 ) ! from the bottom to the surface : ldown 564 zmxlm(ji,jj,jk) = & 565 & MIN( zmxlm(ji,jj,jk+1) + e3t(ji,jj,jk+1,Kmm), zmxlm(ji,jj,jk) ) 566 END_3D 567 DO_3D( 0, 0, 0, 0, 2, jpkm1 ) 498 568 zemlm = MIN ( zmxld(ji,jj,jk), zmxlm(ji,jj,jk) ) 499 569 zemlp = SQRT( zmxld(ji,jj,jk) * zmxlm(ji,jj,jk) ) … … 507 577 ! ! Vertical eddy viscosity and diffusivity (avm and avt) 508 578 ! !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 509 DO_3D _00_00( 1, jpkm1 )579 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) !* vertical eddy viscosity & diffivity at w-points 510 580 zsqen = SQRT( en(ji,jj,jk) ) 511 581 zav = rn_ediff * zmxlm(ji,jj,jk) * zsqen … … 516 586 ! 517 587 ! 518 IF( nn_pdl == 1 ) THEN !* Prandtl number case: update avt519 DO_3D _00_00(2, jpkm1 )520 p_avt(ji,jj,jk) = MAX( apdlr(ji,jj,jk) * p_avt(ji,jj,jk), avtb_2d(ji,jj) * avtb(jk) ) * tmask(ji,jj,jk)588 IF( nn_pdl == 1 ) THEN !* Prandtl number case: update avt 589 DO_3D( 0, 0, 0, 0, 2, jpkm1 ) 590 p_avt(ji,jj,jk) = MAX( apdlr(ji,jj,jk) * p_avt(ji,jj,jk), avtb_2d(ji,jj) * avtb(jk) ) * wmask(ji,jj,jk) 521 591 END_3D 522 592 ENDIF … … 550 620 INTEGER :: ios 551 621 !! 552 NAMELIST/namzdf_tke/ rn_ediff, rn_ediss , rn_ebb , rn_emin , & 553 & rn_emin0, rn_bshear, nn_mxl , ln_mxl0 , & 554 & rn_mxl0 , nn_pdl , ln_drg , ln_lc , rn_lc, & 555 & nn_etau , nn_htau , rn_efr , rn_eice 622 NAMELIST/namzdf_tke/ rn_ediff, rn_ediss , rn_ebb , rn_emin , & 623 & rn_emin0, rn_bshear, nn_mxl , ln_mxl0 , & 624 & rn_mxl0 , nn_mxlice, rn_mxlice, & 625 & nn_pdl , ln_lc , rn_lc , & 626 & nn_etau , nn_htau , rn_efr , nn_eice 556 627 !!---------------------------------------------------------------------- 557 628 ! … … 580 651 WRITE(numout,*) ' surface mixing length = F(stress) or not ln_mxl0 = ', ln_mxl0 581 652 WRITE(numout,*) ' surface mixing length minimum value rn_mxl0 = ', rn_mxl0 582 WRITE(numout,*) ' top/bottom friction forcing flag ln_drg = ', ln_drg 653 IF( ln_mxl0 ) THEN 654 WRITE(numout,*) ' type of scaling under sea-ice nn_mxlice = ', nn_mxlice 655 IF( nn_mxlice == 1 ) & 656 WRITE(numout,*) ' ice thickness when scaling under sea-ice rn_mxlice = ', rn_mxlice 657 SELECT CASE( nn_mxlice ) ! Type of scaling under sea-ice 658 CASE( 0 ) ; WRITE(numout,*) ' ==>>> No scaling under sea-ice' 659 CASE( 1 ) ; WRITE(numout,*) ' ==>>> scaling with constant sea-ice thickness' 660 CASE( 2 ) ; WRITE(numout,*) ' ==>>> scaling with mean sea-ice thickness' 661 CASE( 3 ) ; WRITE(numout,*) ' ==>>> scaling with max sea-ice thickness' 662 CASE DEFAULT 663 CALL ctl_stop( 'zdf_tke_init: wrong value for nn_mxlice, should be 0,1,2,3 or 4') 664 END SELECT 665 ENDIF 583 666 WRITE(numout,*) ' Langmuir cells parametrization ln_lc = ', ln_lc 584 667 WRITE(numout,*) ' coef to compute vertical velocity of LC rn_lc = ', rn_lc … … 586 669 WRITE(numout,*) ' type of tke penetration profile nn_htau = ', nn_htau 587 670 WRITE(numout,*) ' fraction of TKE that penetrates rn_efr = ', rn_efr 588 WRITE(numout,*) ' below sea-ice: =0 ON rn_eice = ', rn_eice 589 WRITE(numout,*) ' =4 OFF when ice fraction > 1/4 ' 590 IF( ln_drg ) THEN 671 WRITE(numout,*) ' langmuir & surface wave breaking under ice nn_eice = ', nn_eice 672 SELECT CASE( nn_eice ) 673 CASE( 0 ) ; WRITE(numout,*) ' ==>>> no impact of ice cover on langmuir & surface wave breaking' 674 CASE( 1 ) ; WRITE(numout,*) ' ==>>> weigthed by 1-TANH( fr_i(:,:) * 10 )' 675 CASE( 2 ) ; WRITE(numout,*) ' ==>>> weighted by 1-fr_i(:,:)' 676 CASE( 3 ) ; WRITE(numout,*) ' ==>>> weighted by 1-MIN( 1, 4 * fr_i(:,:) )' 677 CASE DEFAULT 678 CALL ctl_stop( 'zdf_tke_init: wrong value for nn_eice, should be 0,1,2, or 3') 679 END SELECT 680 IF( .NOT.ln_drg_OFF ) THEN 591 681 WRITE(numout,*) 592 682 WRITE(numout,*) ' Namelist namdrg_top/_bot: used values:' … … 674 764 ! 675 765 IF( MIN( id1, id2, id3, id4 ) > 0 ) THEN ! fields exist 676 CALL iom_get( numror, jpdom_auto glo, 'en' , en , ldxios = lrxios )677 CALL iom_get( numror, jpdom_auto glo, 'avt_k', avt_k, ldxios = lrxios )678 CALL iom_get( numror, jpdom_auto glo, 'avm_k', avm_k, ldxios = lrxios )679 CALL iom_get( numror, jpdom_auto glo, 'dissl', dissl, ldxios = lrxios )766 CALL iom_get( numror, jpdom_auto, 'en' , en , ldxios = lrxios ) 767 CALL iom_get( numror, jpdom_auto, 'avt_k', avt_k, ldxios = lrxios ) 768 CALL iom_get( numror, jpdom_auto, 'avm_k', avm_k, ldxios = lrxios ) 769 CALL iom_get( numror, jpdom_auto, 'dissl', dissl, ldxios = lrxios ) 680 770 ELSE ! start TKE from rest 681 771 IF(lwp) WRITE(numout,*)
Note: See TracChangeset
for help on using the changeset viewer.
