Changeset 11949 for NEMO/branches/2019/dev_r11943_MERGE_2019/src/OCE/LDF
- Timestamp:
- 2019-11-22T15:29:17+01:00 (4 years ago)
- Location:
- NEMO/branches/2019/dev_r11943_MERGE_2019/src
- Files:
-
- 4 edited
-
. (modified) (1 prop)
-
OCE/LDF/ldfdyn.F90 (modified) (10 diffs)
-
OCE/LDF/ldfslp.F90 (modified) (22 diffs)
-
OCE/LDF/ldftra.F90 (modified) (21 diffs)
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2019/dev_r11943_MERGE_2019/src
- Property svn:mergeinfo deleted
-
NEMO/branches/2019/dev_r11943_MERGE_2019/src/OCE/LDF/ldfdyn.F90
r11653 r11949 339 339 340 340 341 SUBROUTINE ldf_dyn( kt )341 SUBROUTINE ldf_dyn( kt, Kbb ) 342 342 !!---------------------------------------------------------------------- 343 343 !! *** ROUTINE ldf_dyn *** … … 357 357 !!---------------------------------------------------------------------- 358 358 INTEGER, INTENT(in) :: kt ! time step index 359 INTEGER, INTENT(in) :: Kbb ! ocean time level indices 359 360 ! 360 361 INTEGER :: ji, jj, jk ! dummy loop indices … … 373 374 DO jj = 2, jpjm1 374 375 DO ji = fs_2, fs_jpim1 375 zu2pv2_ij = u b(ji ,jj ,jk) * ub(ji ,jj ,jk) + vb(ji ,jj ,jk) * vb(ji ,jj ,jk)376 zu2pv2_ij_m1 = u b(ji-1,jj ,jk) * ub(ji-1,jj ,jk) + vb(ji ,jj-1,jk) * vb(ji ,jj-1,jk)376 zu2pv2_ij = uu(ji ,jj ,jk,Kbb) * uu(ji ,jj ,jk,Kbb) + vv(ji ,jj ,jk,Kbb) * vv(ji ,jj ,jk,Kbb) 377 zu2pv2_ij_m1 = uu(ji-1,jj ,jk,Kbb) * uu(ji-1,jj ,jk,Kbb) + vv(ji ,jj-1,jk,Kbb) * vv(ji ,jj-1,jk,Kbb) 377 378 zemax = MAX( e1t(ji,jj) , e2t(ji,jj) ) 378 379 ahmt(ji,jj,jk) = SQRT( (zu2pv2_ij + zu2pv2_ij_m1) * r1_288 ) * zemax * tmask(ji,jj,jk) ! 288= 12*12 * 2 … … 381 382 DO jj = 1, jpjm1 382 383 DO ji = 1, fs_jpim1 383 zu2pv2_ij_p1 = u b(ji ,jj+1,jk) * ub(ji ,jj+1,jk) + vb(ji+1,jj ,jk) * vb(ji+1,jj ,jk)384 zu2pv2_ij = u b(ji ,jj ,jk) * ub(ji ,jj ,jk) + vb(ji ,jj ,jk) * vb(ji ,jj ,jk)384 zu2pv2_ij_p1 = uu(ji ,jj+1,jk, Kbb) * uu(ji ,jj+1,jk, Kbb) + vv(ji+1,jj ,jk, Kbb) * vv(ji+1,jj ,jk, Kbb) 385 zu2pv2_ij = uu(ji ,jj ,jk, Kbb) * uu(ji ,jj ,jk, Kbb) + vv(ji ,jj ,jk, Kbb) * vv(ji ,jj ,jk, Kbb) 385 386 zemax = MAX( e1f(ji,jj) , e2f(ji,jj) ) 386 387 ahmf(ji,jj,jk) = SQRT( (zu2pv2_ij + zu2pv2_ij_p1) * r1_288 ) * zemax * fmask(ji,jj,jk) ! 288= 12*12 * 2 … … 392 393 DO jj = 2, jpjm1 393 394 DO ji = fs_2, fs_jpim1 394 zu2pv2_ij = u b(ji ,jj ,jk) * ub(ji ,jj ,jk) + vb(ji ,jj ,jk) * vb(ji ,jj ,jk)395 zu2pv2_ij_m1 = u b(ji-1,jj ,jk) * ub(ji-1,jj ,jk) + vb(ji ,jj-1,jk) * vb(ji ,jj-1,jk)395 zu2pv2_ij = uu(ji ,jj ,jk,Kbb) * uu(ji ,jj ,jk,Kbb) + vv(ji ,jj ,jk,Kbb) * vv(ji ,jj ,jk,Kbb) 396 zu2pv2_ij_m1 = uu(ji-1,jj ,jk,Kbb) * uu(ji-1,jj ,jk,Kbb) + vv(ji ,jj-1,jk,Kbb) * vv(ji ,jj-1,jk,Kbb) 396 397 zemax = MAX( e1t(ji,jj) , e2t(ji,jj) ) 397 398 ahmt(ji,jj,jk) = SQRT( SQRT( (zu2pv2_ij + zu2pv2_ij_m1) * r1_288 ) * zemax ) * zemax * tmask(ji,jj,jk) … … 400 401 DO jj = 1, jpjm1 401 402 DO ji = 1, fs_jpim1 402 zu2pv2_ij_p1 = u b(ji ,jj+1,jk) * ub(ji ,jj+1,jk) + vb(ji+1,jj ,jk) * vb(ji+1,jj ,jk)403 zu2pv2_ij = u b(ji ,jj ,jk) * ub(ji ,jj ,jk) + vb(ji ,jj ,jk) * vb(ji ,jj ,jk)403 zu2pv2_ij_p1 = uu(ji ,jj+1,jk, Kbb) * uu(ji ,jj+1,jk, Kbb) + vv(ji+1,jj ,jk, Kbb) * vv(ji+1,jj ,jk, Kbb) 404 zu2pv2_ij = uu(ji ,jj ,jk, Kbb) * uu(ji ,jj ,jk, Kbb) + vv(ji ,jj ,jk, Kbb) * vv(ji ,jj ,jk, Kbb) 404 405 zemax = MAX( e1f(ji,jj) , e2f(ji,jj) ) 405 406 ahmf(ji,jj,jk) = SQRT( SQRT( (zu2pv2_ij + zu2pv2_ij_p1) * r1_288 ) * zemax ) * zemax * fmask(ji,jj,jk) … … 425 426 DO jj = 2, jpjm1 426 427 DO ji = 2, jpim1 427 zdb = ( ub(ji,jj,jk) * r1_e2u(ji,jj) - ub(ji-1,jj,jk) * r1_e2u(ji-1,jj) ) * r1_e1t(ji,jj) * e2t(ji,jj) & 428 & - ( vb(ji,jj,jk) * r1_e1v(ji,jj) - vb(ji,jj-1,jk) * r1_e1v(ji,jj-1) ) * r1_e2t(ji,jj) * e1t(ji,jj) 428 zdb = ( uu(ji,jj,jk,Kbb) * r1_e2u(ji,jj) - uu(ji-1,jj,jk,Kbb) * r1_e2u(ji-1,jj) ) & 429 & * r1_e1t(ji,jj) * e2t(ji,jj) & 430 & - ( vv(ji,jj,jk,Kbb) * r1_e1v(ji,jj) - vv(ji,jj-1,jk,Kbb) * r1_e1v(ji,jj-1) ) & 431 & * r1_e2t(ji,jj) * e1t(ji,jj) 429 432 dtensq(ji,jj,jk) = zdb * zdb * tmask(ji,jj,jk) 430 433 END DO … … 433 436 DO jj = 1, jpjm1 434 437 DO ji = 1, jpim1 435 zdb = ( ub(ji,jj+1,jk) * r1_e1u(ji,jj+1) - ub(ji,jj,jk) * r1_e1u(ji,jj) ) * r1_e2f(ji,jj) * e1f(ji,jj) & 436 & + ( vb(ji+1,jj,jk) * r1_e2v(ji+1,jj) - vb(ji,jj,jk) * r1_e2v(ji,jj) ) * r1_e1f(ji,jj) * e2f(ji,jj) 438 zdb = ( uu(ji,jj+1,jk,Kbb) * r1_e1u(ji,jj+1) - uu(ji,jj,jk,Kbb) * r1_e1u(ji,jj) ) & 439 & * r1_e2f(ji,jj) * e1f(ji,jj) & 440 & + ( vv(ji+1,jj,jk,Kbb) * r1_e2v(ji+1,jj) - vv(ji,jj,jk,Kbb) * r1_e2v(ji,jj) ) & 441 & * r1_e1f(ji,jj) * e2f(ji,jj) 437 442 dshesq(ji,jj,jk) = zdb * zdb * fmask(ji,jj,jk) 438 443 END DO … … 448 453 DO ji = fs_2, fs_jpim1 449 454 ! 450 zu2pv2_ij = u b(ji ,jj ,jk) * ub(ji ,jj ,jk) + vb(ji ,jj ,jk) * vb(ji ,jj ,jk)451 zu2pv2_ij_m1 = u b(ji-1,jj ,jk) * ub(ji-1,jj ,jk) + vb(ji ,jj-1,jk) * vb(ji ,jj-1,jk)455 zu2pv2_ij = uu(ji ,jj ,jk,Kbb) * uu(ji ,jj ,jk,Kbb) + vv(ji ,jj ,jk,Kbb) * vv(ji ,jj ,jk,Kbb) 456 zu2pv2_ij_m1 = uu(ji-1,jj ,jk,Kbb) * uu(ji-1,jj ,jk,Kbb) + vv(ji ,jj-1,jk,Kbb) * vv(ji ,jj-1,jk,Kbb) 452 457 ! 453 458 zdelta = zcmsmag * esqt(ji,jj) ! L^2 * (C_smag/pi)^2 … … 464 469 DO ji = 1, fs_jpim1 465 470 ! 466 zu2pv2_ij_p1 = u b(ji ,jj+1,jk) * ub(ji ,jj+1,jk) + vb(ji+1,jj ,jk) * vb(ji+1,jj ,jk)467 zu2pv2_ij = u b(ji ,jj ,jk) * ub(ji ,jj ,jk) + vb(ji ,jj ,jk) * vb(ji ,jj ,jk)471 zu2pv2_ij_p1 = uu(ji ,jj+1,jk, kbb) * uu(ji ,jj+1,jk, kbb) + vv(ji+1,jj ,jk, kbb) * vv(ji+1,jj ,jk, kbb) 472 zu2pv2_ij = uu(ji ,jj ,jk, kbb) * uu(ji ,jj ,jk, kbb) + vv(ji ,jj ,jk, kbb) * vv(ji ,jj ,jk, kbb) 468 473 ! 469 474 zdelta = zcmsmag * esqf(ji,jj) ! L^2 * (C_smag/pi)^2 -
NEMO/branches/2019/dev_r11943_MERGE_2019/src/OCE/LDF/ldfslp.F90
r10425 r11949 81 81 CONTAINS 82 82 83 SUBROUTINE ldf_slp( kt, prd, pn2 )83 SUBROUTINE ldf_slp( kt, prd, pn2, Kbb, Kmm ) 84 84 !!---------------------------------------------------------------------- 85 85 !! *** ROUTINE ldf_slp *** … … 107 107 !!---------------------------------------------------------------------- 108 108 INTEGER , INTENT(in) :: kt ! ocean time-step index 109 INTEGER , INTENT(in) :: Kbb, Kmm ! ocean time level indices 109 110 REAL(wp), INTENT(in), DIMENSION(:,:,:) :: prd ! in situ density 110 111 REAL(wp), INTENT(in), DIMENSION(:,:,:) :: pn2 ! Brunt-Vaisala frequency (locally ref.) … … 171 172 ! 172 173 ! !== Slopes just below the mixed layer ==! 173 CALL ldf_slp_mxl( prd, pn2, zgru, zgrv, zdzr ) ! output: uslpml, vslpml, wslpiml, wslpjml174 CALL ldf_slp_mxl( prd, pn2, zgru, zgrv, zdzr, Kmm ) ! output: uslpml, vslpml, wslpiml, wslpjml 174 175 175 176 … … 205 206 ! ! bound the slopes: abs(zw.)<= 1/100 and zb..<0 206 207 ! ! + kxz max= ah slope max =< e1 e3 /(pi**2 2 dt) 207 zbu = MIN( zbu, - z1_slpmax * ABS( zau ) , -7.e+3_wp/e3u _n(ji,jj,jk)* ABS( zau ) )208 zbv = MIN( zbv, - z1_slpmax * ABS( zav ) , -7.e+3_wp/e3v _n(ji,jj,jk)* ABS( zav ) )208 zbu = MIN( zbu, - z1_slpmax * ABS( zau ) , -7.e+3_wp/e3u(ji,jj,jk,Kmm)* ABS( zau ) ) 209 zbv = MIN( zbv, - z1_slpmax * ABS( zav ) , -7.e+3_wp/e3v(ji,jj,jk,Kmm)* ABS( zav ) ) 209 210 ! ! uslp and vslp output in zwz and zww, resp. 210 211 zfi = MAX( omlmask(ji,jj,jk), omlmask(ji+1,jj,jk) ) 211 212 zfj = MAX( omlmask(ji,jj,jk), omlmask(ji,jj+1,jk) ) 212 213 ! thickness of water column between surface and level k at u/v point 213 zdepu = 0.5_wp * ( ( gdept _n (ji,jj,jk) + gdept_n (ji+1,jj,jk) ) &214 - 2 * MAX( risfdep(ji,jj), risfdep(ji+1,jj) ) - e3u _n(ji,jj,miku(ji,jj)) )215 zdepv = 0.5_wp * ( ( gdept _n (ji,jj,jk) + gdept_n (ji,jj+1,jk) ) &216 - 2 * MAX( risfdep(ji,jj), risfdep(ji,jj+1) ) - e3v _n(ji,jj,mikv(ji,jj)) )214 zdepu = 0.5_wp * ( ( gdept (ji,jj,jk,Kmm) + gdept (ji+1,jj,jk,Kmm) ) & 215 - 2 * MAX( risfdep(ji,jj), risfdep(ji+1,jj) ) - e3u(ji,jj,miku(ji,jj),Kmm) ) 216 zdepv = 0.5_wp * ( ( gdept (ji,jj,jk,Kmm) + gdept (ji,jj+1,jk,Kmm) ) & 217 - 2 * MAX( risfdep(ji,jj), risfdep(ji,jj+1) ) - e3v(ji,jj,mikv(ji,jj),Kmm) ) 217 218 ! 218 219 zwz(ji,jj,jk) = ( ( 1._wp - zfi) * zau / ( zbu - zeps ) & … … 224 225 ! zfi = REAL( 1 - 1/(1 + jk / MAX( nmln(ji+1,jj), nmln(ji,jj) ) ), wp ) 225 226 ! zfj = REAL( 1 - 1/(1 + jk / MAX( nmln(ji,jj+1), nmln(ji,jj) ) ), wp ) 226 ! zci = 0.5 * ( gdept _n(ji+1,jj,jk)+gdept_n(ji,jj,jk) ) / MAX( hmlpt(ji,jj), hmlpt(ji+1,jj), 10. ) )227 ! zcj = 0.5 * ( gdept _n(ji,jj+1,jk)+gdept_n(ji,jj,jk) ) / MAX( hmlpt(ji,jj), hmlpt(ji,jj+1), 10. ) )227 ! zci = 0.5 * ( gdept(ji+1,jj,jk,Kmm)+gdept(ji,jj,jk,Kmm) ) / MAX( hmlpt(ji,jj), hmlpt(ji+1,jj), 10. ) ) 228 ! zcj = 0.5 * ( gdept(ji,jj+1,jk,Kmm)+gdept(ji,jj,jk,Kmm) ) / MAX( hmlpt(ji,jj), hmlpt(ji,jj+1), 10. ) ) 228 229 ! zwz(ji,jj,jk) = ( zfi * zai / ( zbi - zeps ) + ( 1._wp - zfi ) * wslpiml(ji,jj) * zci ) * tmask(ji,jj,jk) 229 230 ! zww(ji,jj,jk) = ( zfj * zaj / ( zbj - zeps ) + ( 1._wp - zfj ) * wslpjml(ji,jj) * zcj ) * tmask(ji,jj,jk) … … 296 297 ! ! bound the slopes: abs(zw.)<= 1/100 and zb..<0. 297 298 ! ! + kxz max= ah slope max =< e1 e3 /(pi**2 2 dt) 298 zbi = MIN( zbw ,- 100._wp* ABS( zai ) , -7.e+3_wp/e3w _n(ji,jj,jk)* ABS( zai ) )299 zbj = MIN( zbw , -100._wp* ABS( zaj ) , -7.e+3_wp/e3w _n(ji,jj,jk)* ABS( zaj ) )299 zbi = MIN( zbw ,- 100._wp* ABS( zai ) , -7.e+3_wp/e3w(ji,jj,jk,Kmm)* ABS( zai ) ) 300 zbj = MIN( zbw , -100._wp* ABS( zaj ) , -7.e+3_wp/e3w(ji,jj,jk,Kmm)* ABS( zaj ) ) 300 301 ! ! wslpi and wslpj with ML flattening (output in zwz and zww, resp.) 301 302 zfk = MAX( omlmask(ji,jj,jk), omlmask(ji,jj,jk-1) ) ! zfk=1 in the ML otherwise zfk=0 302 zck = ( gdepw _n(ji,jj,jk) - gdepw_n(ji,jj,mikt(ji,jj) ) ) / MAX( hmlp(ji,jj) - gdepw_n(ji,jj,mikt(ji,jj)), 10._wp )303 zck = ( gdepw(ji,jj,jk,Kmm) - gdepw(ji,jj,mikt(ji,jj),Kmm) ) / MAX( hmlp(ji,jj) - gdepw(ji,jj,mikt(ji,jj),Kmm), 10._wp ) 303 304 zwz(ji,jj,jk) = ( zai / ( zbi - zeps ) * ( 1._wp - zfk ) + zck * wslpiml(ji,jj) * zfk ) * wmask(ji,jj,jk) 304 305 zww(ji,jj,jk) = ( zaj / ( zbj - zeps ) * ( 1._wp - zfk ) + zck * wslpjml(ji,jj) * zfk ) * wmask(ji,jj,jk) … … 375 376 376 377 377 SUBROUTINE ldf_slp_triad ( kt )378 SUBROUTINE ldf_slp_triad ( kt, Kbb, Kmm ) 378 379 !!---------------------------------------------------------------------- 379 380 !! *** ROUTINE ldf_slp_triad *** … … 390 391 !!---------------------------------------------------------------------- 391 392 INTEGER, INTENT( in ) :: kt ! ocean time-step index 393 INTEGER , INTENT(in) :: Kbb, Kmm ! ocean time level indices 392 394 !! 393 395 INTEGER :: ji, jj, jk, jl, ip, jp, kp ! dummy loop indices … … 419 421 DO jj = 1, jpjm1 ! NB: not masked ==> a minimum value is set 420 422 DO ji = 1, fs_jpim1 ! vector opt. 421 zdit = ( ts b(ji+1,jj,jk,jp_tem) - tsb(ji,jj,jk,jp_tem) ) ! i-gradient of T & S at u-point422 zdis = ( ts b(ji+1,jj,jk,jp_sal) - tsb(ji,jj,jk,jp_sal) )423 zdjt = ( ts b(ji,jj+1,jk,jp_tem) - tsb(ji,jj,jk,jp_tem) ) ! j-gradient of T & S at v-point424 zdjs = ( ts b(ji,jj+1,jk,jp_sal) - tsb(ji,jj,jk,jp_sal) )423 zdit = ( ts(ji+1,jj,jk,jp_tem,Kbb) - ts(ji,jj,jk,jp_tem,Kbb) ) ! i-gradient of T & S at u-point 424 zdis = ( ts(ji+1,jj,jk,jp_sal,Kbb) - ts(ji,jj,jk,jp_sal,Kbb) ) 425 zdjt = ( ts(ji,jj+1,jk,jp_tem,Kbb) - ts(ji,jj,jk,jp_tem,Kbb) ) ! j-gradient of T & S at v-point 426 zdjs = ( ts(ji,jj+1,jk,jp_sal,Kbb) - ts(ji,jj,jk,jp_sal,Kbb) ) 425 427 zdxrho_raw = ( - rab_b(ji+ip,jj ,jk,jp_tem) * zdit + rab_b(ji+ip,jj ,jk,jp_sal) * zdis ) * r1_e1u(ji,jj) 426 428 zdyrho_raw = ( - rab_b(ji ,jj+jp,jk,jp_tem) * zdjt + rab_b(ji ,jj+jp,jk,jp_sal) * zdjs ) * r1_e2v(ji,jj) … … 452 454 DO ji = 1, jpi ! vector opt. 453 455 IF( jk+kp > 1 ) THEN ! k-gradient of T & S a jk+kp 454 zdkt = ( ts b(ji,jj,jk+kp-1,jp_tem) - tsb(ji,jj,jk+kp,jp_tem) )455 zdks = ( ts b(ji,jj,jk+kp-1,jp_sal) - tsb(ji,jj,jk+kp,jp_sal) )456 zdkt = ( ts(ji,jj,jk+kp-1,jp_tem,Kbb) - ts(ji,jj,jk+kp,jp_tem,Kbb) ) 457 zdks = ( ts(ji,jj,jk+kp-1,jp_sal,Kbb) - ts(ji,jj,jk+kp,jp_sal,Kbb) ) 456 458 ELSE 457 459 zdkt = 0._wp ! 1st level gradient set to zero … … 460 462 zdzrho_raw = ( - rab_b(ji,jj,jk+kp,jp_tem) * zdkt & 461 463 & + rab_b(ji,jj,jk+kp,jp_sal) * zdks & 462 & ) / e3w _n(ji,jj,jk+kp)464 & ) / e3w(ji,jj,jk+kp,Kmm) 463 465 zdzrho(ji,jj,jk,kp) = - MIN( - repsln , zdzrho_raw ) ! force zdzrho >= repsln 464 466 END DO … … 470 472 DO ji = 1, jpi 471 473 jk = MIN( nmln(ji,jj), mbkt(ji,jj) ) + 1 ! MIN in case ML depth is the ocean depth 472 z1_mlbw(ji,jj) = 1._wp / gdepw _n(ji,jj,jk)474 z1_mlbw(ji,jj) = 1._wp / gdepw(ji,jj,jk,Kmm) 473 475 END DO 474 476 END DO … … 499 501 ! Add s-coordinate slope at t-points (do this by *subtracting* gradient of depth) 500 502 zti_g_raw = ( zdxrho(ji+ip,jj,jk-kp,1-ip) / zdzrho(ji+ip,jj,jk-kp,kp) & 501 & - ( gdept _n(ji+1,jj,jk-kp) - gdept_n(ji,jj,jk-kp) ) * r1_e1u(ji,jj) ) * umask(ji,jj,jk)502 ze3_e1 = e3w _n(ji+ip,jj,jk-kp) * r1_e1u(ji,jj)503 & - ( gdept(ji+1,jj,jk-kp,Kmm) - gdept(ji,jj,jk-kp,Kmm) ) * r1_e1u(ji,jj) ) * umask(ji,jj,jk) 504 ze3_e1 = e3w(ji+ip,jj,jk-kp,Kmm) * r1_e1u(ji,jj) 503 505 zti_mlb(ji+ip,jj ,1-ip,kp) = SIGN( MIN( rn_slpmax, 5.0_wp * ze3_e1 , ABS( zti_g_raw ) ), zti_g_raw ) 504 506 ENDIF … … 509 511 ELSE 510 512 ztj_g_raw = ( zdyrho(ji,jj+jp,jk-kp,1-jp) / zdzrho(ji,jj+jp,jk-kp,kp) & 511 & - ( gdept _n(ji,jj+1,jk-kp) - gdept_n(ji,jj,jk-kp) ) / e2v(ji,jj) ) * vmask(ji,jj,jk)512 ze3_e2 = e3w _n(ji,jj+jp,jk-kp) / e2v(ji,jj)513 & - ( gdept(ji,jj+1,jk-kp,Kmm) - gdept(ji,jj,jk-kp,Kmm) ) / e2v(ji,jj) ) * vmask(ji,jj,jk) 514 ze3_e2 = e3w(ji,jj+jp,jk-kp,Kmm) / e2v(ji,jj) 513 515 ztj_mlb(ji ,jj+jp,1-jp,kp) = SIGN( MIN( rn_slpmax, 5.0_wp * ze3_e2 , ABS( ztj_g_raw ) ), ztj_g_raw ) 514 516 ENDIF … … 542 544 ! 543 545 ! Must mask contribution to slope for triad jk=1,kp=0 that poke up though ocean surface 544 zti_coord = znot_thru_surface * ( gdept _n(ji+1,jj ,jk) - gdept_n(ji,jj,jk) ) * r1_e1u(ji,jj)545 ztj_coord = znot_thru_surface * ( gdept _n(ji ,jj+1,jk) - gdept_n(ji,jj,jk) ) * r1_e2v(ji,jj) ! unmasked546 zti_coord = znot_thru_surface * ( gdept(ji+1,jj ,jk,Kmm) - gdept(ji,jj,jk,Kmm) ) * r1_e1u(ji,jj) 547 ztj_coord = znot_thru_surface * ( gdept(ji ,jj+1,jk,Kmm) - gdept(ji,jj,jk,Kmm) ) * r1_e2v(ji,jj) ! unmasked 546 548 zti_g_raw = zti_raw - zti_coord ! ref to geopot surfaces 547 549 ztj_g_raw = ztj_raw - ztj_coord 548 550 ! additional limit required in bilaplacian case 549 ze3_e1 = e3w _n(ji+ip,jj ,jk+kp) * r1_e1u(ji,jj)550 ze3_e2 = e3w _n(ji ,jj+jp,jk+kp) * r1_e2v(ji,jj)551 ze3_e1 = e3w(ji+ip,jj ,jk+kp,Kmm) * r1_e1u(ji,jj) 552 ze3_e2 = e3w(ji ,jj+jp,jk+kp,Kmm) * r1_e2v(ji,jj) 551 553 ! NB: hard coded factor 5 (can be a namelist parameter...) 552 554 zti_g_lim = SIGN( MIN( rn_slpmax, 5.0_wp * ze3_e1, ABS( zti_g_raw ) ), zti_g_raw ) … … 561 563 zti_g_lim = ( zfacti * zti_g_lim & 562 564 & + ( 1._wp - zfacti ) * zti_mlb(ji+ip,jj,1-ip,kp) & 563 & * gdepw _n(ji+ip,jj,jk+kp) * z1_mlbw(ji+ip,jj) ) * umask(ji,jj,jk+kp)565 & * gdepw(ji+ip,jj,jk+kp,Kmm) * z1_mlbw(ji+ip,jj) ) * umask(ji,jj,jk+kp) 564 566 ztj_g_lim = ( zfactj * ztj_g_lim & 565 567 & + ( 1._wp - zfactj ) * ztj_mlb(ji,jj+jp,1-jp,kp) & 566 & * gdepw _n(ji,jj+jp,jk+kp) * z1_mlbw(ji,jj+jp) ) * vmask(ji,jj,jk+kp)568 & * gdepw(ji,jj+jp,jk+kp,Kmm) * z1_mlbw(ji,jj+jp) ) * vmask(ji,jj,jk+kp) 567 569 ! 568 570 triadi_g(ji+ip,jj ,jk,1-ip,kp) = zti_g_lim … … 596 598 triadj(ji ,jj+jp,jk,1-jp,kp) = ztj_lim * zjsw 597 599 ! 598 zbu = e1e2u(ji ,jj ) * e3u _n(ji ,jj ,jk)599 zbv = e1e2v(ji ,jj ) * e3v _n(ji ,jj ,jk)600 zbti = e1e2t(ji+ip,jj ) * e3w _n(ji+ip,jj ,jk+kp)601 zbtj = e1e2t(ji ,jj+jp) * e3w _n(ji ,jj+jp,jk+kp)600 zbu = e1e2u(ji ,jj ) * e3u(ji ,jj ,jk ,Kmm) 601 zbv = e1e2v(ji ,jj ) * e3v(ji ,jj ,jk ,Kmm) 602 zbti = e1e2t(ji+ip,jj ) * e3w(ji+ip,jj ,jk+kp,Kmm) 603 zbtj = e1e2t(ji ,jj+jp) * e3w(ji ,jj+jp,jk+kp,Kmm) 602 604 ! 603 605 wslp2(ji+ip,jj,jk+kp) = wslp2(ji+ip,jj,jk+kp) + 0.25_wp * zbu / zbti * zti_g_lim*zti_g_lim ! masked … … 618 620 619 621 620 SUBROUTINE ldf_slp_mxl( prd, pn2, p_gru, p_grv, p_dzr )622 SUBROUTINE ldf_slp_mxl( prd, pn2, p_gru, p_grv, p_dzr, Kmm ) 621 623 !!---------------------------------------------------------------------- 622 624 !! *** ROUTINE ldf_slp_mxl *** … … 638 640 REAL(wp), DIMENSION(:,:,:), INTENT(in) :: p_gru, p_grv ! i- & j-gradient of density (u- & v-pts) 639 641 REAL(wp), DIMENSION(:,:,:), INTENT(in) :: p_dzr ! z-gradient of density (T-point) 642 INTEGER , INTENT(in) :: Kmm ! ocean time level indices 640 643 !! 641 644 INTEGER :: ji , jj , jk ! dummy loop indices … … 694 697 ! !- bound the slopes: abs(zw.)<= 1/100 and zb..<0 695 698 ! kxz max= ah slope max =< e1 e3 /(pi**2 2 dt) 696 zbu = MIN( zbu , - z1_slpmax * ABS( zau ) , -7.e+3_wp/e3u _n(ji,jj,iku)* ABS( zau ) )697 zbv = MIN( zbv , - z1_slpmax * ABS( zav ) , -7.e+3_wp/e3v _n(ji,jj,ikv)* ABS( zav ) )699 zbu = MIN( zbu , - z1_slpmax * ABS( zau ) , -7.e+3_wp/e3u(ji,jj,iku,Kmm)* ABS( zau ) ) 700 zbv = MIN( zbv , - z1_slpmax * ABS( zav ) , -7.e+3_wp/e3v(ji,jj,ikv,Kmm)* ABS( zav ) ) 698 701 ! !- Slope at u- & v-points (uslpml, vslpml) 699 702 uslpml(ji,jj) = zau / ( zbu - zeps ) * umask(ji,jj,iku) … … 717 720 ! !- bound the slopes: abs(zw.)<= 1/100 and zb..<0. 718 721 ! kxz max= ah slope max =< e1 e3 /(pi**2 2 dt) 719 zbi = MIN( zbw , -100._wp* ABS( zai ) , -7.e+3_wp/e3w _n(ji,jj,ik)* ABS( zai ) )720 zbj = MIN( zbw , -100._wp* ABS( zaj ) , -7.e+3_wp/e3w _n(ji,jj,ik)* ABS( zaj ) )722 zbi = MIN( zbw , -100._wp* ABS( zai ) , -7.e+3_wp/e3w(ji,jj,ik,Kmm)* ABS( zai ) ) 723 zbj = MIN( zbw , -100._wp* ABS( zaj ) , -7.e+3_wp/e3w(ji,jj,ik,Kmm)* ABS( zaj ) ) 721 724 ! !- i- & j-slope at w-points (wslpiml, wslpjml) 722 725 wslpiml(ji,jj) = zai / ( zbi - zeps ) * tmask (ji,jj,ik) … … 786 789 ! DO jj = 2, jpjm1 787 790 ! DO ji = fs_2, fs_jpim1 ! vector opt. 788 ! uslp (ji,jj,jk) = - ( gdept _n(ji+1,jj,jk) - gdept_n(ji ,jj ,jk) ) * r1_e1u(ji,jj) * umask(ji,jj,jk)789 ! vslp (ji,jj,jk) = - ( gdept _n(ji,jj+1,jk) - gdept_n(ji ,jj ,jk) ) * r1_e2v(ji,jj) * vmask(ji,jj,jk)790 ! wslpi(ji,jj,jk) = - ( gdepw _n(ji+1,jj,jk) - gdepw_n(ji-1,jj,jk) ) * r1_e1t(ji,jj) * wmask(ji,jj,jk) * 0.5791 ! wslpj(ji,jj,jk) = - ( gdepw _n(ji,jj+1,jk) - gdepw_n(ji,jj-1,jk) ) * r1_e2t(ji,jj) * wmask(ji,jj,jk) * 0.5791 ! uslp (ji,jj,jk) = - ( gdept(ji+1,jj,jk,Kmm) - gdept(ji ,jj ,jk,Kmm) ) * r1_e1u(ji,jj) * umask(ji,jj,jk) 792 ! vslp (ji,jj,jk) = - ( gdept(ji,jj+1,jk,Kmm) - gdept(ji ,jj ,jk,Kmm) ) * r1_e2v(ji,jj) * vmask(ji,jj,jk) 793 ! wslpi(ji,jj,jk) = - ( gdepw(ji+1,jj,jk,Kmm) - gdepw(ji-1,jj,jk,Kmm) ) * r1_e1t(ji,jj) * wmask(ji,jj,jk) * 0.5 794 ! wslpj(ji,jj,jk) = - ( gdepw(ji,jj+1,jk,Kmm) - gdepw(ji,jj-1,jk,Kmm) ) * r1_e2t(ji,jj) * wmask(ji,jj,jk) * 0.5 792 795 ! END DO 793 796 ! END DO -
NEMO/branches/2019/dev_r11943_MERGE_2019/src/OCE/LDF/ldftra.F90
r11536 r11949 382 382 383 383 384 SUBROUTINE ldf_tra( kt )384 SUBROUTINE ldf_tra( kt, Kbb, Kmm ) 385 385 !!---------------------------------------------------------------------- 386 386 !! *** ROUTINE ldf_tra *** … … 403 403 !!---------------------------------------------------------------------- 404 404 INTEGER, INTENT(in) :: kt ! time step 405 INTEGER, INTENT(in) :: Kbb, Kmm ! ocean time level indices 405 406 ! 406 407 INTEGER :: ji, jj, jk ! dummy loop indices … … 411 412 ! ! =F(growth rate of baroclinic instability) 412 413 ! ! max value aeiv_0 ; decreased to 0 within 20N-20S 413 CALL ldf_eiv( kt, aei0, aeiu, aeiv )414 CALL ldf_eiv( kt, aei0, aeiu, aeiv, Kmm ) 414 415 ENDIF 415 416 ! … … 424 425 ahtv(:,:,1) = aeiv(:,:,1) 425 426 ELSE ! compute aht. 426 CALL ldf_eiv( kt, aht0, ahtu, ahtv )427 CALL ldf_eiv( kt, aht0, ahtu, ahtv, Kmm ) 427 428 ENDIF 428 429 ! … … 448 449 IF( ln_traldf_lap ) THEN ! laplacian operator |u| e /12 449 450 DO jk = 1, jpkm1 450 ahtu(:,:,jk) = ABS( u b(:,:,jk) ) * e1u(:,:) * r1_12 ! n.b. ub,vbare masked451 ahtv(:,:,jk) = ABS( v b(:,:,jk) ) * e2v(:,:) * r1_12451 ahtu(:,:,jk) = ABS( uu(:,:,jk,Kbb) ) * e1u(:,:) * r1_12 ! n.b. uu,vv are masked 452 ahtv(:,:,jk) = ABS( vv(:,:,jk,Kbb) ) * e2v(:,:) * r1_12 452 453 END DO 453 454 ELSEIF( ln_traldf_blp ) THEN ! bilaplacian operator sqrt( |u| e^3 /12 ) = sqrt( |u| e /12 ) * e 454 455 DO jk = 1, jpkm1 455 ahtu(:,:,jk) = SQRT( ABS( u b(:,:,jk) ) * e1u(:,:) * r1_12 ) * e1u(:,:)456 ahtv(:,:,jk) = SQRT( ABS( v b(:,:,jk) ) * e2v(:,:) * r1_12 ) * e2v(:,:)456 ahtu(:,:,jk) = SQRT( ABS( uu(:,:,jk,Kbb) ) * e1u(:,:) * r1_12 ) * e1u(:,:) 457 ahtv(:,:,jk) = SQRT( ABS( vv(:,:,jk,Kbb) ) * e2v(:,:) * r1_12 ) * e2v(:,:) 457 458 END DO 458 459 ENDIF … … 625 626 626 627 627 SUBROUTINE ldf_eiv( kt, paei0, paeiu, paeiv )628 SUBROUTINE ldf_eiv( kt, paei0, paeiu, paeiv, Kmm ) 628 629 !!---------------------------------------------------------------------- 629 630 !! *** ROUTINE ldf_eiv *** … … 637 638 !!---------------------------------------------------------------------- 638 639 INTEGER , INTENT(in ) :: kt ! ocean time-step index 640 INTEGER , INTENT(in ) :: Kmm ! ocean time level indices 639 641 REAL(wp) , INTENT(inout) :: paei0 ! max value [m2/s] 640 642 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) :: paeiu, paeiv ! eiv coefficient [m2/s] … … 658 660 ! internal Rossby radius Ro = .5 * sum_jpk(N) / f 659 661 zn2 = MAX( rn2b(ji,jj,jk), 0._wp ) 660 zn(ji,jj) = zn(ji,jj) + SQRT( zn2 ) * e3w _n(ji,jj,jk)662 zn(ji,jj) = zn(ji,jj) + SQRT( zn2 ) * e3w(ji,jj,jk,Kmm) 661 663 ! Compute elements required for the inverse time scale of baroclinic 662 664 ! eddies using the isopycnal slopes calculated in ldfslp.F : 663 665 ! T^-1 = sqrt(m_jpk(N^2*(r1^2+r2^2)*e3w)) 664 ze3w = e3w _n(ji,jj,jk) * tmask(ji,jj,jk)666 ze3w = e3w(ji,jj,jk,Kmm) * tmask(ji,jj,jk) 665 667 zah(ji,jj) = zah(ji,jj) + zn2 * wslp2(ji,jj,jk) * ze3w 666 668 zhw(ji,jj) = zhw(ji,jj) + ze3w … … 676 678 ! internal Rossby radius Ro = .5 * sum_jpk(N) / f 677 679 zn2 = MAX( rn2b(ji,jj,jk), 0._wp ) 678 zn(ji,jj) = zn(ji,jj) + SQRT( zn2 ) * e3w _n(ji,jj,jk)680 zn(ji,jj) = zn(ji,jj) + SQRT( zn2 ) * e3w(ji,jj,jk,Kmm) 679 681 ! Compute elements required for the inverse time scale of baroclinic 680 682 ! eddies using the isopycnal slopes calculated in ldfslp.F : 681 683 ! T^-1 = sqrt(m_jpk(N^2*(r1^2+r2^2)*e3w)) 682 ze3w = e3w _n(ji,jj,jk) * tmask(ji,jj,jk)684 ze3w = e3w(ji,jj,jk,Kmm) * tmask(ji,jj,jk) 683 685 zah(ji,jj) = zah(ji,jj) + zn2 * ( wslpi(ji,jj,jk) * wslpi(ji,jj,jk) & 684 686 & + wslpj(ji,jj,jk) * wslpj(ji,jj,jk) ) * ze3w … … 725 727 726 728 727 SUBROUTINE ldf_eiv_trp( kt, kit000, pu n, pvn, pwn, cdtype)729 SUBROUTINE ldf_eiv_trp( kt, kit000, pu, pv, pw, cdtype, Kmm, Krhs ) 728 730 !!---------------------------------------------------------------------- 729 731 !! *** ROUTINE ldf_eiv_trp *** … … 741 743 !! velocity and heat transport (call ldf_eiv_dia) 742 744 !! 743 !! ** Action : pun, pvn increased by the eiv transport 744 !!---------------------------------------------------------------------- 745 INTEGER , INTENT(in ) :: kt ! ocean time-step index 746 INTEGER , INTENT(in ) :: kit000 ! first time step index 747 CHARACTER(len=3) , INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator) 748 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) :: pun ! in : 3 ocean transport components [m3/s] 749 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) :: pvn ! out: 3 ocean transport components [m3/s] 750 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) :: pwn ! increased by the eiv [m3/s] 745 !! ** Action : pu, pv increased by the eiv transport 746 !!---------------------------------------------------------------------- 747 INTEGER , INTENT(in ) :: kt ! ocean time-step index 748 INTEGER , INTENT(in ) :: kit000 ! first time step index 749 INTEGER , INTENT(in ) :: Kmm, Krhs ! ocean time level indices 750 CHARACTER(len=3) , INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator) 751 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) :: pu ! in : 3 ocean transport components [m3/s] 752 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) :: pv ! out: 3 ocean transport components [m3/s] 753 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) :: pw ! increased by the eiv [m3/s] 751 754 !! 752 755 INTEGER :: ji, jj, jk ! dummy loop indices … … 780 783 DO jj = 1, jpjm1 781 784 DO ji = 1, fs_jpim1 ! vector opt. 782 pu n(ji,jj,jk) = pun(ji,jj,jk) - ( zpsi_uw(ji,jj,jk) - zpsi_uw(ji,jj,jk+1) )783 pv n(ji,jj,jk) = pvn(ji,jj,jk) - ( zpsi_vw(ji,jj,jk) - zpsi_vw(ji,jj,jk+1) )785 pu(ji,jj,jk) = pu(ji,jj,jk) - ( zpsi_uw(ji,jj,jk) - zpsi_uw(ji,jj,jk+1) ) 786 pv(ji,jj,jk) = pv(ji,jj,jk) - ( zpsi_vw(ji,jj,jk) - zpsi_vw(ji,jj,jk+1) ) 784 787 END DO 785 788 END DO … … 788 791 DO jj = 2, jpjm1 789 792 DO ji = fs_2, fs_jpim1 ! vector opt. 790 pw n(ji,jj,jk) = pwn(ji,jj,jk) + ( zpsi_uw(ji,jj,jk) - zpsi_uw(ji-1,jj ,jk) &793 pw(ji,jj,jk) = pw(ji,jj,jk) + ( zpsi_uw(ji,jj,jk) - zpsi_uw(ji-1,jj ,jk) & 791 794 & + zpsi_vw(ji,jj,jk) - zpsi_vw(ji ,jj-1,jk) ) 792 795 END DO … … 795 798 ! 796 799 ! ! diagnose the eddy induced velocity and associated heat transport 797 IF( ln_ldfeiv_dia .AND. cdtype == 'TRA' ) CALL ldf_eiv_dia( zpsi_uw, zpsi_vw )800 IF( ln_ldfeiv_dia .AND. cdtype == 'TRA' ) CALL ldf_eiv_dia( zpsi_uw, zpsi_vw, Kmm ) 798 801 ! 799 802 END SUBROUTINE ldf_eiv_trp 800 803 801 804 802 SUBROUTINE ldf_eiv_dia( psi_uw, psi_vw )805 SUBROUTINE ldf_eiv_dia( psi_uw, psi_vw, Kmm ) 803 806 !!---------------------------------------------------------------------- 804 807 !! *** ROUTINE ldf_eiv_dia *** … … 811 814 !!---------------------------------------------------------------------- 812 815 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) :: psi_uw, psi_vw ! streamfunction [m3/s] 816 INTEGER , INTENT(in ) :: Kmm ! ocean time level indices 813 817 ! 814 818 INTEGER :: ji, jj, jk ! dummy loop indices … … 831 835 ! 832 836 DO jk = 1, jpkm1 ! e2u e3u u_eiv = -dk[psi_uw] 833 zw3d(:,:,jk) = ( psi_uw(:,:,jk+1) - psi_uw(:,:,jk) ) / ( e2u(:,:) * e3u _n(:,:,jk) )837 zw3d(:,:,jk) = ( psi_uw(:,:,jk+1) - psi_uw(:,:,jk) ) / ( e2u(:,:) * e3u(:,:,jk,Kmm) ) 834 838 END DO 835 839 CALL iom_put( "uoce_eiv", zw3d ) 836 840 ! 837 841 DO jk = 1, jpkm1 ! e1v e3v v_eiv = -dk[psi_vw] 838 zw3d(:,:,jk) = ( psi_vw(:,:,jk+1) - psi_vw(:,:,jk) ) / ( e1v(:,:) * e3v _n(:,:,jk) )842 zw3d(:,:,jk) = ( psi_vw(:,:,jk+1) - psi_vw(:,:,jk) ) / ( e1v(:,:) * e3v(:,:,jk,Kmm) ) 839 843 END DO 840 844 CALL iom_put( "voce_eiv", zw3d ) … … 859 863 DO jj = 2, jpjm1 860 864 DO ji = fs_2, fs_jpim1 ! vector opt. 861 zw3d(ji,jj,jk) = zw3d(ji,jj,jk) + ( psi_uw(ji,jj,jk+1) - psi_uw(ji,jj,jk)) &862 & * ( ts n (ji,jj,jk,jp_tem) + tsn (ji+1,jj,jk,jp_tem) )865 zw3d(ji,jj,jk) = zw3d(ji,jj,jk) + ( psi_uw(ji,jj,jk+1) - psi_uw(ji ,jj,jk) ) & 866 & * ( ts (ji,jj,jk,jp_tem,Kmm) + ts (ji+1,jj,jk,jp_tem,Kmm) ) 863 867 zw2d(ji,jj) = zw2d(ji,jj) + zw3d(ji,jj,jk) 864 868 END DO … … 875 879 DO jj = 2, jpjm1 876 880 DO ji = fs_2, fs_jpim1 ! vector opt. 877 zw3d(ji,jj,jk) = zw3d(ji,jj,jk) + ( psi_vw(ji,jj,jk+1) - psi_vw(ji,jj,jk)) &878 & * ( ts n (ji,jj,jk,jp_tem) + tsn (ji,jj+1,jk,jp_tem) )881 zw3d(ji,jj,jk) = zw3d(ji,jj,jk) + ( psi_vw(ji,jj,jk+1) - psi_vw(ji,jj ,jk) ) & 882 & * ( ts (ji,jj,jk,jp_tem,Kmm) + ts (ji,jj+1,jk,jp_tem,Kmm) ) 879 883 zw2d(ji,jj) = zw2d(ji,jj) + zw3d(ji,jj,jk) 880 884 END DO … … 894 898 DO jj = 2, jpjm1 895 899 DO ji = fs_2, fs_jpim1 ! vector opt. 896 zw3d(ji,jj,jk) = zw3d(ji,jj,jk) * ( psi_uw(ji,jj,jk+1) - psi_uw(ji,jj,jk)) &897 & * ( ts n (ji,jj,jk,jp_sal) + tsn (ji+1,jj,jk,jp_sal) )900 zw3d(ji,jj,jk) = zw3d(ji,jj,jk) * ( psi_uw(ji,jj,jk+1) - psi_uw(ji ,jj,jk) ) & 901 & * ( ts (ji,jj,jk,jp_sal,Kmm) + ts (ji+1,jj,jk,jp_sal,Kmm) ) 898 902 zw2d(ji,jj) = zw2d(ji,jj) + zw3d(ji,jj,jk) 899 903 END DO … … 903 907 CALL lbc_lnk( 'ldftra', zw3d, 'U', -1. ) 904 908 CALL iom_put( "ueiv_salttr", zztmp * zw2d ) ! salt transport in i-direction 905 CALL iom_put( "ueiv_salttr3d", zztmp * zw3d ) 909 CALL iom_put( "ueiv_salttr3d", zztmp * zw3d ) ! salt transport in i-direction 906 910 ENDIF 907 911 zw2d(:,:) = 0._wp … … 910 914 DO jj = 2, jpjm1 911 915 DO ji = fs_2, fs_jpim1 ! vector opt. 912 zw3d(ji,jj,jk) = zw3d(ji,jj,jk) + ( psi_vw(ji,jj,jk+1) - psi_vw(ji,jj,jk)) &913 & * ( ts n (ji,jj,jk,jp_sal) + tsn (ji,jj+1,jk,jp_sal) )916 zw3d(ji,jj,jk) = zw3d(ji,jj,jk) + ( psi_vw(ji,jj,jk+1) - psi_vw(ji,jj ,jk) ) & 917 & * ( ts (ji,jj,jk,jp_sal,Kmm) + ts (ji,jj+1,jk,jp_sal,Kmm) ) 914 918 zw2d(ji,jj) = zw2d(ji,jj) + zw3d(ji,jj,jk) 915 919 END DO
Note: See TracChangeset
for help on using the changeset viewer.
