Changeset 14909
- Timestamp:
- 2021-05-26T18:04:24+02:00 (3 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2021/dev_r14122_HPC-08_Mueller_OSMOSIS_streamlining/src/OCE/ZDF/zdfosm.F90
r14901 r14909 34 34 !! 16/05/2019 (19) Fox-Kemper parametrization of restratification through mixed layer eddies added to revised code. 35 35 !! 23/05/19 (20) Old code where key_osmldpth1` is *not* set removed, together with the key key_osmldpth1 36 !! 4.2 ! 2021-05 (S. Mueller) Efficiency and source-code clarity enhancements, adaptation for tiling 36 37 !!---------------------------------------------------------------------- 37 38 … … 282 283 zdf_osm_alloc = zdf_osm_alloc + ierr 283 284 ! 284 ALLOCATE( l_conv(A2D((nn_hls-1))), l_shear(A2D((nn_hls-1))), l_coup(A2D((nn_hls-1))), l_pyc(A2D((nn_hls-1))), l_flux(A2D((nn_hls-1))), l_mle(A2D((nn_hls-1))), STAT=ierr ) 285 ALLOCATE( l_conv(A2D((nn_hls-1))), l_shear(A2D((nn_hls-1))), l_coup(A2D((nn_hls-1))), l_pyc(A2D((nn_hls-1))), & 286 & l_flux(A2D((nn_hls-1))), l_mle(A2D((nn_hls-1))), STAT=ierr ) 285 287 zdf_osm_alloc = zdf_osm_alloc + ierr 286 288 ! 287 ALLOCATE( swth0(A2D((nn_hls-1))), sws0(A2D((nn_hls-1))), swb0(A2D((nn_hls-1))), suw0(A2D((nn_hls-1))), sustar(A2D((nn_hls-1))), scos_wind(A2D((nn_hls-1))), & 288 & ssin_wind(A2D((nn_hls-1))), swthav(A2D((nn_hls-1))), swsav(A2D((nn_hls-1))), swbav(A2D((nn_hls-1))), sustke(A2D((nn_hls-1))), dstokes(A2D((nn_hls-1))), & 289 & swstrl(A2D((nn_hls-1))), swstrc(A2D((nn_hls-1))), sla(A2D((nn_hls-1))), svstr(A2D((nn_hls-1))), shol(A2D((nn_hls-1))), STAT=ierr ) 289 ALLOCATE( swth0(A2D((nn_hls-1))), sws0(A2D((nn_hls-1))), swb0(A2D((nn_hls-1))), suw0(A2D((nn_hls-1))), & 290 & sustar(A2D((nn_hls-1))), scos_wind(A2D((nn_hls-1))), ssin_wind(A2D((nn_hls-1))), swthav(A2D((nn_hls-1))), & 291 & swsav(A2D((nn_hls-1))), swbav(A2D((nn_hls-1))), sustke(A2D((nn_hls-1))), dstokes(A2D((nn_hls-1))), & 292 & swstrl(A2D((nn_hls-1))), swstrc(A2D((nn_hls-1))), sla(A2D((nn_hls-1))), svstr(A2D((nn_hls-1))), & 293 & shol(A2D((nn_hls-1))), STAT=ierr ) 290 294 zdf_osm_alloc = zdf_osm_alloc + ierr 291 295 ! 292 ALLOCATE( av_t_bl(A2D((nn_hls-1))), av_s_bl(A2D((nn_hls-1))), av_u_bl(A2D((nn_hls-1))), av_v_bl(A2D((nn_hls-1))), av_b_bl(A2D((nn_hls-1))), STAT=ierr) 296 ALLOCATE( av_t_bl(A2D((nn_hls-1))), av_s_bl(A2D((nn_hls-1))), av_u_bl(A2D((nn_hls-1))), av_v_bl(A2D((nn_hls-1))), & 297 & av_b_bl(A2D((nn_hls-1))), STAT=ierr) 293 298 zdf_osm_alloc = zdf_osm_alloc + ierr 294 299 ! 295 ALLOCATE( av_dt_bl(A2D((nn_hls-1))), av_ds_bl(A2D((nn_hls-1))), av_du_bl(A2D((nn_hls-1))), av_dv_bl(A2D((nn_hls-1))), av_db_bl(A2D((nn_hls-1))), STAT=ierr) 300 ALLOCATE( av_dt_bl(A2D((nn_hls-1))), av_ds_bl(A2D((nn_hls-1))), av_du_bl(A2D((nn_hls-1))), av_dv_bl(A2D((nn_hls-1))), & 301 & av_db_bl(A2D((nn_hls-1))), STAT=ierr) 296 302 zdf_osm_alloc = zdf_osm_alloc + ierr 297 303 ! 298 ALLOCATE( av_t_ml(A2D((nn_hls-1))), av_s_ml(A2D((nn_hls-1))), av_u_ml(A2D((nn_hls-1))), av_v_ml(A2D((nn_hls-1))), av_b_ml(A2D((nn_hls-1))), STAT=ierr) 304 ALLOCATE( av_t_ml(A2D((nn_hls-1))), av_s_ml(A2D((nn_hls-1))), av_u_ml(A2D((nn_hls-1))), av_v_ml(A2D((nn_hls-1))), & 305 & av_b_ml(A2D((nn_hls-1))), STAT=ierr) 299 306 zdf_osm_alloc = zdf_osm_alloc + ierr 300 307 ! 301 ALLOCATE( av_dt_ml(A2D((nn_hls-1))), av_ds_ml(A2D((nn_hls-1))), av_du_ml(A2D((nn_hls-1))), av_dv_ml(A2D((nn_hls-1))), av_db_ml(A2D((nn_hls-1))), STAT=ierr) 308 ALLOCATE( av_dt_ml(A2D((nn_hls-1))), av_ds_ml(A2D((nn_hls-1))), av_du_ml(A2D((nn_hls-1))), av_dv_ml(A2D((nn_hls-1))), & 309 & av_db_ml(A2D((nn_hls-1))), STAT=ierr) 302 310 zdf_osm_alloc = zdf_osm_alloc + ierr 303 311 ! 304 ALLOCATE( av_t_mle(A2D((nn_hls-1))), av_s_mle(A2D((nn_hls-1))), av_u_mle(A2D((nn_hls-1))), av_v_mle(A2D((nn_hls-1))), av_b_mle(A2D((nn_hls-1))), STAT=ierr) 312 ALLOCATE( av_t_mle(A2D((nn_hls-1))), av_s_mle(A2D((nn_hls-1))), av_u_mle(A2D((nn_hls-1))), av_v_mle(A2D((nn_hls-1))), & 313 & av_b_mle(A2D((nn_hls-1))), STAT=ierr) 305 314 zdf_osm_alloc = zdf_osm_alloc + ierr 306 315 ! … … 364 373 REAL(wp), DIMENSION(A2D((nn_hls-1))) :: zrad0 ! Surface solar temperature flux (deg m/s) 365 374 REAL(wp), DIMENSION(A2D((nn_hls-1))) :: zradh ! Radiative flux at bl base (Buoyancy units) 366 REAL(wp) :: zradav ! Radiative flux, bl average (Buoyancy Units)367 REAL(wp) :: zvw0 ! Surface v-momentum flux375 REAL(wp) :: zradav ! Radiative flux, bl average (Buoyancy Units) 376 REAL(wp) :: zvw0 ! Surface v-momentum flux 368 377 REAL(wp), DIMENSION(A2D((nn_hls-1))) :: zwb0tot ! Total surface buoyancy flux including insolation 369 378 REAL(wp), DIMENSION(A2D((nn_hls-1))) :: zwb_ent ! Buoyancy entrainment flux … … 380 389 REAL(wp), DIMENSION(A2D((nn_hls-1))) :: zhml ! ML depth - grid 381 390 ! 382 REAL(wp), DIMENSION(A2D((nn_hls-1))) :: zhmle ! MLE depth - grid 383 REAL(wp), DIMENSION(A2D(nn_hls)) :: zmld ! ML depth on grid 384 ! 385 REAL(wp), DIMENSION(A2D((nn_hls-1))) :: zdh ! Pycnocline depth - grid 386 REAL(wp), DIMENSION(A2D((nn_hls-1))) :: zdhdt ! BL depth tendency 387 REAL(wp), DIMENSION(A2D((nn_hls-1))) :: zdtdz_bl_ext, zdsdz_bl_ext, zdbdz_bl_ext ! External temperature/salinity and buoyancy gradients 388 REAL(wp), DIMENSION(A2D(nn_hls)) :: zdtdx, zdtdy, zdsdx, zdsdy ! Horizontal gradients for Fox-Kemper parametrization 391 REAL(wp), DIMENSION(A2D((nn_hls-1))) :: zhmle ! MLE depth - grid 392 REAL(wp), DIMENSION(A2D(nn_hls)) :: zmld ! ML depth on grid 393 ! 394 REAL(wp), DIMENSION(A2D((nn_hls-1))) :: zdh ! Pycnocline depth - grid 395 REAL(wp), DIMENSION(A2D((nn_hls-1))) :: zdhdt ! BL depth tendency 396 REAL(wp), DIMENSION(A2D((nn_hls-1))) :: zdtdz_bl_ext, zdsdz_bl_ext ! External temperature/salinity gradients 397 REAL(wp), DIMENSION(A2D((nn_hls-1))) :: zdbdz_bl_ext ! External buoyancy gradients 398 REAL(wp), DIMENSION(A2D(nn_hls)) :: zdtdx, zdtdy, zdsdx, zdsdy ! Horizontal gradients for Fox-Kemper parametrization 389 399 ! 390 400 REAL(wp), DIMENSION(A2D((nn_hls-1))) :: zdbds_mle ! Magnitude of horizontal buoyancy gradient … … 400 410 ! 401 411 ! Temporary variables 402 REAL(wp) :: znd! Temporary non-dimensional depth403 REAL(wp) :: zz0, zz1, zfac404 REAL(wp) :: zus_x, zus_y! Temporary Stokes drift405 REAL(wp), DIMENSION(A2D((nn_hls-1)),jpk) :: zviscos ! Viscosity406 REAL(wp), DIMENSION(A2D((nn_hls-1)),jpk) :: zdiffut ! t-diffusivity407 REAL(wp) :: zabsstke408 REAL(wp) :: zsqrtpi, z_two_thirds, zthickness409 REAL(wp) :: z2k_times_thickness, zsqrt_depth, zexp_depth, zf, zexperfc412 REAL(wp) :: znd ! Temporary non-dimensional depth 413 REAL(wp) :: zz0, zz1, zfac 414 REAL(wp) :: zus_x, zus_y ! Temporary Stokes drift 415 REAL(wp), DIMENSION(A2D((nn_hls-1)),jpk) :: zviscos ! Viscosity 416 REAL(wp), DIMENSION(A2D((nn_hls-1)),jpk) :: zdiffut ! t-diffusivity 417 REAL(wp) :: zabsstke 418 REAL(wp) :: zsqrtpi, z_two_thirds, zthickness 419 REAL(wp) :: z2k_times_thickness, zsqrt_depth, zexp_depth, zf, zexperfc 410 420 ! 411 421 ! For debugging … … 421 431 ! Mixed layer 422 432 ! No initialization of zhbl or zhml (or zdh?) 423 zhbl(:,:) = pp_large ; zhml(:,:) = pp_large ; zdh(:,:) = pp_large 433 zhbl(:,:) = pp_large 434 zhml(:,:) = pp_large 435 zdh(:,:) = pp_large 424 436 ! 425 437 IF ( ln_osm_mle ) THEN ! Only initialise arrays if needed … … 521 533 IF ( hsw(ji,jj) > 1e-4_wp ) THEN 522 534 ! Use wave fields 523 zabsstke = SQRT( ut0sd(ji,jj)**2 + vt0sd(ji,jj)**2 )535 zabsstke = SQRT( ut0sd(ji,jj)**2 + vt0sd(ji,jj)**2 ) 524 536 sustke(ji,jj) = MAX( ( scos_wind(ji,jj) * ut0sd(ji,jj) + ssin_wind(ji,jj) * vt0sd(ji,jj) ), 1e-8_wp ) 525 537 dstokes(ji,jj) = MAX( zfac * hsw(ji,jj) * hsw(ji,jj) / ( MAX( zabsstke * wmp(ji,jj), 1e-7 ) ), 5e-1_wp ) … … 578 590 zexperfc = zsqrtpi * zsqrt_depth * ERFC(zsqrt_depth) * EXP( z2k_times_thickness ) 579 591 ELSE 580 ! Asymptotic expansion of sqrt(pi)*zsqrt_depth*EXP(z2k_times_thickness)*ERFC(zsqrt_depth) for large z2k_times_thickness 592 ! Asymptotic expansion of sqrt(pi)*zsqrt_depth*EXP(z2k_times_thickness)*ERFC(zsqrt_depth) for large 593 ! z2k_times_thickness 581 594 ! See Abramowitz and Stegun, Eq. 7.1.23 582 595 ! zexperfc = 1._wp - (1/2)/(z2k_times_thickness) + (3/4)/(z2k_times_thickness**2) - (15/8)/(z2k_times_thickness**3) … … 599 612 IF ( sla(ji,jj) > 0.45_wp ) dstokes(ji,jj) = MIN( dstokes(ji,jj), 0.5_wp * hbl(ji,jj) ) 600 613 ! Velocity scale that tends to sustar for large Langmuir numbers 601 svstr(ji,jj) = ( swstrl(ji,jj)**3 + ( 1.0_wp - EXP( -0.5_wp * sla(ji,jj)**2 ) ) * sustar(ji,jj) * sustar(ji,jj) * &602 & sustar(ji,jj) )**pthird614 svstr(ji,jj) = ( swstrl(ji,jj)**3 + ( 1.0_wp - EXP( -0.5_wp * sla(ji,jj)**2 ) ) * sustar(ji,jj) * sustar(ji,jj) * & 615 & sustar(ji,jj) )**pthird 603 616 ! 604 617 ! Limit maximum value of Langmuir number as approximate treatment for shear turbulence … … 677 690 END_2D 678 691 DO_3D_OVR( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 5, jpkm1 ) 679 IF ( hmle(ji,jj) >= gdepw(ji,jj,jk,Kmm) ) mld_prof(ji,jj) = MIN( mbkt(ji,jj), jk)692 IF ( hmle(ji,jj) >= gdepw(ji,jj,jk,Kmm) ) mld_prof(ji,jj) = MIN( mbkt(ji,jj), jk) 680 693 END_3D 681 694 CALL zdf_osm_vertical_average( Kbb, Kmm, & … … 984 997 CALL zdf_osm_iomput( "zhbl", tmask(A2D(0),1) * zhbl(A2D(0)) ) ! BL depth internal to zdf_osm routine 985 998 CALL zdf_osm_iomput( "zhml", tmask(A2D(0),1) * zhml(A2D(0)) ) ! ML depth internal to zdf_osm routine 986 CALL zdf_osm_iomput( "imld", tmask(A2D(0),1) * nmld(A2D(0)) ) ! Index for ML depth internal to zdf_osm routine 987 CALL zdf_osm_iomput( "jp_ext", tmask(A2D(0),1) * jk_ext(A2D(0)) ) ! =1 if pycnocline resolved internal to zdf_osm routine 988 CALL zdf_osm_iomput( "j_ddh", tmask(A2D(0),1) * n_ddh(A2D(0)) ) ! Index forpyc thicknessh internal to zdf_osm routine 989 CALL zdf_osm_iomput( "zshear", tmask(A2D(0),1) * zshear(A2D(0)) ) ! Shear production of TKE internal to zdf_osm routine 990 CALL zdf_osm_iomput( "zdh", tmask(A2D(0),1) * zdh(A2D(0)) ) ! Pyc thicknessh internal to zdf_osm routine 999 CALL zdf_osm_iomput( "imld", tmask(A2D(0),1) * nmld(A2D(0)) ) ! Index for ML depth internal to zdf_osm 1000 ! ! routine 1001 CALL zdf_osm_iomput( "jp_ext", tmask(A2D(0),1) * jk_ext(A2D(0)) ) ! =1 if pycnocline resolved internal to 1002 ! ! zdf_osm routine 1003 CALL zdf_osm_iomput( "j_ddh", tmask(A2D(0),1) * n_ddh(A2D(0)) ) ! Index forpyc thicknessh internal to 1004 ! ! zdf_osm routine 1005 CALL zdf_osm_iomput( "zshear", tmask(A2D(0),1) * zshear(A2D(0)) ) ! Shear production of TKE internal to 1006 ! ! zdf_osm routine 1007 CALL zdf_osm_iomput( "zdh", tmask(A2D(0),1) * zdh(A2D(0)) ) ! Pyc thicknessh internal to zdf_osm 1008 ! ! routine 991 1009 CALL zdf_osm_iomput( "zhol", tmask(A2D(0),1) * shol(A2D(0)) ) ! ML depth internal to zdf_osm routine 992 1010 CALL zdf_osm_iomput( "zwb_ent", tmask(A2D(0),1) * zwb_ent(A2D(0)) ) ! Upward turb buoyancy entrainment flux … … 1008 1026 ! Lateral boundary conditions on ghamu and ghamv, currently on W-grid (sign unchanged), needed to caclulate gham[uv] on u and 1009 1027 ! v grids 1010 IF ( .NOT. l_istiled .OR. ntile == nijtile ) THEN ! Finalise ghamu, ghamv, hbl, and hmle only after full domain has been processed 1028 IF ( .NOT. l_istiled .OR. ntile == nijtile ) THEN ! Finalise ghamu, ghamv, hbl, and hmle only after full domain has been 1029 ! ! processed 1011 1030 IF ( nn_hls == 1 ) CALL lbc_lnk( 'zdfosm', ghamu, 'W', 1.0_wp, & 1012 1031 & ghamv, 'W', 1.0_wp ) … … 1014 1033 DO jj = Njs0, Nje0 1015 1034 DO ji = Nis0, Nie0 1016 ghamu(ji,jj,jk) = ( ghamu(ji,jj,jk) + ghamu(ji+1,jj,jk) ) / MAX( 1.0_wp, tmask(ji,jj,jk) + tmask (ji+1,jj,jk) ) *&1017 & umask(ji,jj,jk)1018 ghamv(ji,jj,jk) = ( ghamv(ji,jj,jk) + ghamv(ji,jj+1,jk) ) / MAX( 1.0_wp, tmask(ji,jj,jk) + tmask (ji,jj+1,jk) ) *&1019 & vmask(ji,jj,jk)1035 ghamu(ji,jj,jk) = ( ghamu(ji,jj,jk) + ghamu(ji+1,jj,jk) ) / & 1036 & MAX( 1.0_wp, tmask(ji,jj,jk) + tmask (ji+1,jj,jk) ) * umask(ji,jj,jk) 1037 ghamv(ji,jj,jk) = ( ghamv(ji,jj,jk) + ghamv(ji,jj+1,jk) ) / & 1038 & MAX( 1.0_wp, tmask(ji,jj,jk) + tmask (ji,jj+1,jk) ) * vmask(ji,jj,jk) 1020 1039 ghamt(ji,jj,jk) = ghamt(ji,jj,jk) * tmask(ji,jj,jk) 1021 1040 ghams(ji,jj,jk) = ghams(ji,jj,jk) * tmask(ji,jj,jk) … … 1027 1046 & hmle, 'T', 1.0_wp ) 1028 1047 ! 1029 CALL zdf_osm_iomput( "ghamt", tmask *ghamt )! <Tw_NL>1030 CALL zdf_osm_iomput( "ghams", tmask *ghams )! <Sw_NL>1031 CALL zdf_osm_iomput( "ghamu", umask *ghamu )! <uw_NL>1032 CALL zdf_osm_iomput( "ghamv", vmask *ghamv )! <vw_NL>1033 CALL zdf_osm_iomput( "hbl", tmask(:,:,1) * hbl ) ! Boundary-layer depth1034 CALL zdf_osm_iomput( "hmle", tmask(:,:,1) *hmle )! FK layer depth1048 CALL zdf_osm_iomput( "ghamt", tmask * ghamt ) ! <Tw_NL> 1049 CALL zdf_osm_iomput( "ghams", tmask * ghams ) ! <Sw_NL> 1050 CALL zdf_osm_iomput( "ghamu", umask * ghamu ) ! <uw_NL> 1051 CALL zdf_osm_iomput( "ghamv", vmask * ghamv ) ! <vw_NL> 1052 CALL zdf_osm_iomput( "hbl", tmask(:,:,1) * hbl ) ! Boundary-layer depth 1053 CALL zdf_osm_iomput( "hmle", tmask(:,:,1) * hmle ) ! FK layer depth 1035 1054 END IF 1036 1055 ! … … 1051 1070 !! knlev+kp_ext 1052 1071 !!---------------------------------------------------------------------- 1053 INTEGER, INTENT(in ) :: Kbb, Kmm ! Ocean time-level indices1072 INTEGER, INTENT(in ) :: Kbb, Kmm ! Ocean time-level indices 1054 1073 INTEGER, DIMENSION(A2D((nn_hls-1))), INTENT(in ) :: knlev ! Number of levels to average over. 1055 1074 REAL(wp), DIMENSION(A2D((nn_hls-1))), INTENT( out) :: pt, ps ! Average temperature and salinity … … 1057 1076 REAL(wp), DIMENSION(A2D((nn_hls-1))), INTENT( out) :: pu, pv ! Average current components 1058 1077 INTEGER, DIMENSION(A2D((nn_hls-1))), INTENT(in ), OPTIONAL :: kp_ext ! External-level offsets 1059 REAL(wp), DIMENSION(A2D((nn_hls-1))), INTENT( out), OPTIONAL :: pdt, pds ! Difference between average temperature, salinity, 1060 REAL(wp), DIMENSION(A2D((nn_hls-1))), INTENT( out), OPTIONAL :: pdb ! buoyancy, 1061 REAL(wp), DIMENSION(A2D((nn_hls-1))), INTENT( out), OPTIONAL :: pdu, pdv ! velocity components and the OSBL 1062 ! 1063 INTEGER :: jk, jkflt, jkmax, ji, jj ! Loop indices 1064 INTEGER :: ibld_ext ! External-layer index 1078 REAL(wp), DIMENSION(A2D((nn_hls-1))), INTENT( out), OPTIONAL :: pdt ! Difference between average temperature, 1079 REAL(wp), DIMENSION(A2D((nn_hls-1))), INTENT( out), OPTIONAL :: pds ! salinity, 1080 REAL(wp), DIMENSION(A2D((nn_hls-1))), INTENT( out), OPTIONAL :: pdb ! buoyancy, and 1081 REAL(wp), DIMENSION(A2D((nn_hls-1))), INTENT( out), OPTIONAL :: pdu, pdv ! velocity components and the OSBL 1082 ! 1083 INTEGER :: jk, jkflt, jkmax, ji, jj ! Loop indices 1084 INTEGER :: ibld_ext ! External-layer index 1065 1085 REAL(wp), DIMENSION(A2D((nn_hls-1))) :: zthick ! Layer thickness 1066 REAL(wp) :: zthermal, zbeta ! Thermal/haline expansion/contraction coefficients 1086 REAL(wp) :: zthermal ! Thermal expansion coefficient 1087 REAL(wp) :: zbeta ! Haline contraction coefficient 1067 1088 !!---------------------------------------------------------------------- 1068 1089 ! … … 1154 1175 !! 1155 1176 !!---------------------------------------------------------------------- 1156 REAL(wp), INTENT(inout), DIMENSION(A2D((nn_hls-1))) :: pu, pv 1157 LOGICAL, OPTIONAL, INTENT(in ) :: fwd! Forward (default) or reverse rotation1177 REAL(wp), INTENT(inout), DIMENSION(A2D((nn_hls-1))) :: pu, pv ! Components of current 1178 LOGICAL, OPTIONAL, INTENT(in ) :: fwd ! Forward (default) or reverse rotation 1158 1179 ! 1159 1180 INTEGER :: ji, jj ! Loop indices … … 1184 1205 !! 1185 1206 !!---------------------------------------------------------------------- 1186 REAL(wp), INTENT(inout), DIMENSION(A2D((nn_hls-1)),jpk) :: pu, pv 1187 LOGICAL, OPTIONAL, INTENT(in ) :: fwd! Forward (default) or reverse rotation1188 INTEGER, OPTIONAL, INTENT(in ) :: ktop! Minimum depth index1189 INTEGER, OPTIONAL, INTENT(in ), DIMENSION(A2D((nn_hls-1))) :: knlev 1207 REAL(wp), INTENT(inout), DIMENSION(A2D((nn_hls-1)),jpk) :: pu, pv ! Components of current 1208 LOGICAL, OPTIONAL, INTENT(in ) :: fwd ! Forward (default) or reverse rotation 1209 INTEGER, OPTIONAL, INTENT(in ) :: ktop ! Minimum depth index 1210 INTEGER, OPTIONAL, INTENT(in ), DIMENSION(A2D((nn_hls-1))) :: knlev ! Array of maximum depth indices 1190 1211 ! 1191 1212 INTEGER :: ji, jj, jk, jktop, jkmax ! Loop indices … … 1230 1251 !! 1231 1252 !!---------------------------------------------------------------------- 1232 INTEGER, INTENT(in ) :: Kmm ! Ocean time-level index1253 INTEGER, INTENT(in ) :: Kmm ! Ocean time-level index 1233 1254 REAL(wp), DIMENSION(A2D((nn_hls-1))), INTENT( out) :: pwb_ent ! Buoyancy fluxes at base 1234 1255 REAL(wp), DIMENSION(A2D((nn_hls-1))), INTENT( out) :: pwb_min ! of well-mixed layer … … 1243 1264 REAL(wp), DIMENSION(A2D((nn_hls-1))) :: zekman 1244 1265 REAL(wp), DIMENSION(A2D((nn_hls-1))) :: zri_p, zri_b ! Richardson numbers 1245 REAL(wp) :: zshear_u, zshear_v, zwb_shr1246 REAL(wp) :: zwcor, zrf_conv, zrf_shear, zrf_langmuir, zr_stokes1266 REAL(wp) :: zshear_u, zshear_v, zwb_shr 1267 REAL(wp) :: zwcor, zrf_conv, zrf_shear, zrf_langmuir, zr_stokes 1247 1268 ! 1248 1269 REAL(wp), PARAMETER :: pp_a_shr = 0.4_wp, pp_b_shr = 6.5_wp, pp_a_wb_s = 0.8_wp … … 1278 1299 pshear(ji,jj) = 0.0_wp 1279 1300 END_2D 1280 zekman(:,:) = EXP( -1.0_wp * pp_ek * ABS( ff_t(A2D((nn_hls-1))) ) * phbl(A2D((nn_hls-1))) / MAX( sustar(A2D((nn_hls-1))), 1.e-8 ) ) 1301 zekman(:,:) = EXP( -1.0_wp * pp_ek * ABS( ff_t(A2D((nn_hls-1))) ) * phbl(A2D((nn_hls-1))) / & 1302 & MAX( sustar(A2D((nn_hls-1))), 1.e-8 ) ) 1281 1303 ! 1282 1304 DO_2D_OVR( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 ) 1283 1305 IF ( l_conv(ji,jj) ) THEN 1284 1306 IF ( av_db_bl(ji,jj) > 0.0_wp ) THEN 1285 zri_p(ji,jj) = MAX ( SQRT( av_db_bl(ji,jj) * pdh(ji,jj) / MAX( av_du_bl(ji,jj)**2 + av_dv_bl(ji,jj)**2, 1307 zri_p(ji,jj) = MAX ( SQRT( av_db_bl(ji,jj) * pdh(ji,jj) / MAX( av_du_bl(ji,jj)**2 + av_dv_bl(ji,jj)**2, & 1286 1308 & 1e-8_wp ) ) * ( phbl(ji,jj) / pdh(ji,jj) ) * & 1287 1309 & ( svstr(ji,jj) / MAX( sustar(ji,jj), 1e-6_wp ) )**2 / & … … 1305 1327 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 1306 1328 IF ( pshear(ji,jj) > 1e-10 ) THEN 1307 IF ( zri_p(ji,jj) < pp_ri_p_thresh .AND. MIN( hu(ji,jj,Kmm), hu(ji-1,jj,Kmm), hv(ji,jj,Kmm), hv(ji,jj-1,Kmm) ) > 100.0_wp ) THEN 1329 IF ( zri_p(ji,jj) < pp_ri_p_thresh .AND. & 1330 & MIN( hu(ji,jj,Kmm), hu(ji-1,jj,Kmm), hv(ji,jj,Kmm), hv(ji,jj-1,Kmm) ) > 100.0_wp ) THEN 1308 1331 ! Growing shear layer 1309 1332 n_ddh(ji,jj) = 0 … … 1361 1384 ! Developing shear layer, additional shear production possible. 1362 1385 1363 ! 1364 ! 1365 ! 1386 ! pshear_u = MAX( zustar(ji,jj)**2 * MAX( av_du_ml(ji,jj), 0._wp ) / phbl(ji,jj), 0._wp ) 1387 ! pshear(ji,jj) = pshear(ji,jj) + pshear_u * ( 1.0 - MIN( zri_p(ji,jj) / pp_ri_p_thresh, 1.d0 )**2 ) 1388 ! pshear(ji,jj) = MIN( pshear(ji,jj), pshear_u ) 1366 1389 1367 ! 1368 ! 1390 ! zwb_shr = zwb_shr - 0.25 * MAX ( pshear_u, 0._wp) * ( 1.0 - MIN( zri_p(ji,jj) / pp_ri_p_thresh, 1._wp )**2 ) 1391 ! zwb_shr = MAX( zwb_shr, -0.25 * pshear_u ) 1369 1392 ENDIF 1370 1393 pwb_ent(ji,jj) = pwb_ent(ji,jj) + zwb_shr … … 1391 1414 !! 1392 1415 !!---------------------------------------------------------------------- 1393 INTEGER, INTENT(in ) :: Kmm ! Ocean time-level index 1394 INTEGER, DIMENSION(A2D((nn_hls-1))), INTENT(in ) :: kbase ! OSBL base layer index 1395 REAL(wp), DIMENSION(A2D((nn_hls-1))), INTENT( out) :: pdtdz, pdsdz, pdbdz ! External gradients of temperature, salinity and buoyancy 1416 INTEGER, INTENT(in ) :: Kmm ! Ocean time-level index 1417 INTEGER, DIMENSION(A2D((nn_hls-1))), INTENT(in ) :: kbase ! OSBL base layer index 1418 REAL(wp), DIMENSION(A2D((nn_hls-1))), INTENT( out) :: pdtdz, pdsdz ! External gradients of temperature, salinity 1419 REAL(wp), DIMENSION(A2D((nn_hls-1))), INTENT( out) :: pdbdz ! and buoyancy 1396 1420 ! 1397 1421 ! Local variables … … 1476 1500 ! ! entrainment > restratification 1477 1501 IF ( av_db_bl(ji,jj) > 1e-15_wp ) THEN 1478 zgamma_b_nd = MAX( pdbdz_bl_ext(ji,jj), 0.0_wp ) * pdh(ji,jj) / ( zvel_max + MAX( av_db_bl(ji,jj), 1e-15_wp ) ) 1479 zpsi = ( 1.0_wp - 0.5_wp * pdh(ji,jj) / phbl(ji,jj) ) * & 1480 & ( swb0(ji,jj) - MIN( ( pwb_min(ji,jj) + 2.0_wp * pwb_fk_b(ji,jj) ), 0.0_wp ) ) * pdh(ji,jj) / phbl(ji,jj) 1502 zgamma_b_nd = MAX( pdbdz_bl_ext(ji,jj), 0.0_wp ) * pdh(ji,jj) / & 1503 & ( zvel_max + MAX( av_db_bl(ji,jj), 1e-15_wp ) ) 1504 zpsi = ( 1.0_wp - 0.5_wp * pdh(ji,jj) / phbl(ji,jj) ) * & 1505 & ( swb0(ji,jj) - MIN( ( pwb_min(ji,jj) + 2.0_wp * pwb_fk_b(ji,jj) ), 0.0_wp ) ) * pdh(ji,jj) / & 1506 & phbl(ji,jj) 1481 1507 zpsi = zpsi + 1.75_wp * ( 1.0_wp - 0.5_wp * pdh(ji,jj) / phbl(ji,jj) ) * & 1482 & ( pdh(ji,jj) / phbl(ji,jj) + zgamma_b_nd ) * MIN( ( pwb_min(ji,jj) + 2.0_wp * pwb_fk_b(ji,jj) ), 0.0_wp ) 1508 & ( pdh(ji,jj) / phbl(ji,jj) + zgamma_b_nd ) * & 1509 & MIN( ( pwb_min(ji,jj) + 2.0_wp * pwb_fk_b(ji,jj) ), 0.0_wp ) 1483 1510 zpsi = pp_alpha_b * MAX( zpsi, 0.0_wp ) 1484 pdhdt(ji,jj) = -1.0_wp * ( pwb_ent(ji,jj) + 2.0_wp * pwb_fk_b(ji,jj) ) / &1511 pdhdt(ji,jj) = -1.0_wp * ( pwb_ent(ji,jj) + 2.0_wp * pwb_fk_b(ji,jj) ) / & 1485 1512 & ( zvel_max + MAX( av_db_bl(ji,jj), 1e-15_wp ) ) + & 1486 1513 & zpsi / ( zvel_max + MAX( av_db_bl(ji,jj), 1e-15_wp ) ) … … 1488 1515 IF ( ( swstrc(ji,jj) / svstr(ji,jj) )**3 <= 0.5_wp ) THEN 1489 1516 zari = MIN( 1.5_wp * av_db_bl(ji,jj) / & 1490 & ( phbl(ji,jj) * ( MAX( pdbdz_bl_ext(ji,jj), 0.0_wp ) + &1517 & ( phbl(ji,jj) * ( MAX( pdbdz_bl_ext(ji,jj), 0.0_wp ) + & 1491 1518 & av_db_bl(ji,jj)**2 / MAX( 4.5_wp * svstr(ji,jj)**2, & 1492 1519 & 1e-12_wp ) ) ), 0.2_wp ) 1493 1520 ELSE 1494 1521 zari = MIN( 1.5_wp * av_db_bl(ji,jj) / & 1495 & ( phbl(ji,jj) * ( MAX( pdbdz_bl_ext(ji,jj), 0.0_wp ) + &1522 & ( phbl(ji,jj) * ( MAX( pdbdz_bl_ext(ji,jj), 0.0_wp ) + & 1496 1523 & av_db_bl(ji,jj)**2 / MAX( 4.5_wp * swstrc(ji,jj)**2, & 1497 1524 & 1e-12_wp ) ) ), 0.2_wp ) … … 1601 1628 !! 1602 1629 !!---------------------------------------------------------------------- 1603 INTEGER, INTENT(in ) :: Kmm ! Ocean time-level index1630 INTEGER, INTENT(in ) :: Kmm ! Ocean time-level index 1604 1631 REAL(wp), DIMENSION(A2D((nn_hls-1))), INTENT(inout) :: pdhdt ! Rates of change of hbl 1605 1632 REAL(wp), DIMENSION(A2D((nn_hls-1))), INTENT(inout) :: phbl ! BL depth … … 1701 1728 !! 1702 1729 !!---------------------------------------------------------------------- 1703 INTEGER, INTENT(in ) :: Kmm ! Ocean time-level index1730 INTEGER, INTENT(in ) :: Kmm ! Ocean time-level index 1704 1731 REAL(wp), DIMENSION(A2D((nn_hls-1))), INTENT(inout) :: pdh ! Pycnocline thickness 1705 1732 REAL(wp), DIMENSION(A2D((nn_hls-1))), INTENT(inout) :: phml ! ML depth … … 1789 1816 ztmp = swstrc(ji,jj) 1790 1817 END IF 1791 zari = MIN( 1.5_wp * av_db_bl(ji,jj) / ( phbl(ji,jj) *&1792 & ( MAX( pdbdz_bl_ext(ji,jj), 0.0_wp) + &1793 & av_db_bl(ji,jj)**2 / MAX(4.5_wp * ztmp**2 , 1e-12_wp ) ) ), 0.2_wp )1818 zari = MIN( 1.5_wp * av_db_bl(ji,jj) / & 1819 & ( phbl(ji,jj) * ( MAX( pdbdz_bl_ext(ji,jj), 0.0_wp ) + & 1820 & av_db_bl(ji,jj)**2 / MAX( 4.5_wp * ztmp**2 , 1e-12_wp ) ) ), 0.2_wp ) 1794 1821 ELSE 1795 1822 zari = 0.2_wp … … 1807 1834 ztmp = swstrc(ji,jj) 1808 1835 END IF 1809 zari = MIN( 1.5_wp * av_db_bl(ji,jj) / ( phbl(ji,jj) *&1810 & 1811 & 1836 zari = MIN( 1.5_wp * av_db_bl(ji,jj) / & 1837 & ( phbl(ji,jj) * ( MAX( pdbdz_bl_ext(ji,jj), 0.0_wp ) + & 1838 & av_db_bl(ji,jj)**2 / MAX( 4.5_wp * ztmp**2 , 1e-12_wp ) ) ), 0.2_wp ) 1812 1839 ELSE 1813 1840 zari = 0.2_wp … … 1862 1889 !! 1863 1890 !!---------------------------------------------------------------------- 1864 INTEGER, INTENT(in ) :: Kmm ! Ocean time-level index1891 INTEGER, INTENT(in ) :: Kmm ! Ocean time-level index 1865 1892 INTEGER, DIMENSION(A2D((nn_hls-1))), INTENT(in ) :: kp_ext ! External-level offsets 1866 1893 REAL(wp), DIMENSION(A2D((nn_hls-1)),jpk), INTENT( out) :: pdbdz ! Gradients in the pycnocline … … 1946 1973 END IF ! IF (shol >=0.5) 1947 1974 END IF ! IF (av_db_bl> 0.) 1948 END IF ! IF (pdhdt >= 0) pdhdt < 0 not considered since pycnocline profile is zero and profile arrays are intialized to zero 1975 END IF ! IF (pdhdt >= 0) pdhdt < 0 not considered since pycnocline profile is zero and profile arrays are 1976 ! ! intialized to zero 1949 1977 ! 1950 1978 END IF ! IF (l_conv) … … 1972 2000 !! 1973 2001 !!---------------------------------------------------------------------- 1974 INTEGER, INTENT(in ) :: Kbb, Kmm ! Ocean time-level indices2002 INTEGER, INTENT(in ) :: Kbb, Kmm ! Ocean time-level indices 1975 2003 REAL(wp), DIMENSION(A2D((nn_hls-1)),jpk), INTENT(inout) :: pdiffut ! t-diffusivity 1976 2004 REAL(wp), DIMENSION(A2D((nn_hls-1)),jpk), INTENT(inout) :: pviscos ! Viscosity … … 2168 2196 !! 2169 2197 !!---------------------------------------------------------------------- 2170 INTEGER, INTENT(in ) :: Kmm ! Time-level index2198 INTEGER, INTENT(in ) :: Kmm ! Time-level index 2171 2199 INTEGER, DIMENSION(A2D((nn_hls-1))), INTENT(in ) :: kp_ext ! Offset for external level 2172 2200 REAL(wp), DIMENSION(A2D((nn_hls-1))), INTENT(in ) :: phbl ! BL depth … … 2185 2213 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: z3ddz_pyc_1, z3ddz_pyc_2 ! Pycnocline gradient/shear profiles 2186 2214 ! 2187 INTEGER :: ji, jj, jk, jkm_bld, jkf_mld, jkm_mld ! Loop indices2188 INTEGER :: istat ! Memory allocation status2189 REAL(wp) :: zznd_d, zznd_ml, zznd_pyc, znd ! Temporary non-dimensional depths2215 INTEGER :: ji, jj, jk, jkm_bld, jkf_mld, jkm_mld ! Loop indices 2216 INTEGER :: istat ! Memory allocation status 2217 REAL(wp) :: zznd_d, zznd_ml, zznd_pyc, znd ! Temporary non-dimensional depths 2190 2218 REAL(wp), DIMENSION(A2D((nn_hls-1))) :: zsc_wth_1,zsc_ws_1 ! Temporary scales 2191 2219 REAL(wp), DIMENSION(A2D((nn_hls-1))) :: zsc_uw_1, zsc_uw_2 ! Temporary scales 2192 2220 REAL(wp), DIMENSION(A2D((nn_hls-1))) :: zsc_vw_1, zsc_vw_2 ! Temporary scales 2193 2221 REAL(wp), DIMENSION(A2D((nn_hls-1))) :: ztau_sc_u ! Dissipation timescale at base of WML 2194 REAL(wp) :: zbuoy_pyc_sc, zdelta_pyc !2195 REAL(wp) :: zl_c,zl_l,zl_eps ! Used to calculate turbulence length scale2222 REAL(wp) :: zbuoy_pyc_sc, zdelta_pyc ! 2223 REAL(wp) :: zl_c,zl_l,zl_eps ! Used to calculate turbulence length scale 2196 2224 REAL(wp), DIMENSION(A2D((nn_hls-1))) :: za_cubic, zb_cubic ! Coefficients in cubic polynomial specifying 2197 REAL(wp), DIMENSION(A2D((nn_hls-1))) :: zc_cubic, zd_cubic ! diffusivity in pycnocline2225 REAL(wp), DIMENSION(A2D((nn_hls-1))) :: zc_cubic, zd_cubic ! diffusivity in pycnocline 2198 2226 REAL(wp), DIMENSION(A2D((nn_hls-1))) :: zwt_pyc_sc_1, zws_pyc_sc_1 ! 2199 2227 REAL(wp), DIMENSION(A2D((nn_hls-1))) :: zzeta_pyc ! 2200 REAL(wp) :: zomega, zvw_max !2228 REAL(wp) :: zomega, zvw_max ! 2201 2229 REAL(wp), DIMENSION(A2D((nn_hls-1))) :: zuw_bse,zvw_bse ! Momentum, heat, and salinity fluxes 2202 REAL(wp), DIMENSION(A2D((nn_hls-1))) :: zwth_ent,zws_ent ! at the top of the pycnocline2230 REAL(wp), DIMENSION(A2D((nn_hls-1))) :: zwth_ent,zws_ent ! at the top of the pycnocline 2203 2231 REAL(wp), DIMENSION(A2D((nn_hls-1))) :: zsc_wth_pyc, zsc_ws_pyc ! Scales for pycnocline transport term 2204 REAL(wp) :: ztmp ! 2205 REAL(wp) :: ztgrad, zsgrad, zbgrad ! Variables used to calculate pycnocline gradients 2206 REAL(wp) :: zugrad, zvgrad ! Variables for calculating pycnocline shear 2207 REAL(wp) :: zdtdz_pyc ! Parametrized gradient of temperature in pycnocline 2208 REAL(wp) :: zdsdz_pyc ! Parametrised gradient of salinity in pycnocline 2209 REAL(wp) :: zdudz_pyc ! u-shear across the pycnocline 2210 REAL(wp) :: zdvdz_pyc ! v-shear across the pycnocline 2211 !!---------------------------------------------------------------------- 2232 REAL(wp) :: ztmp ! 2233 REAL(wp) :: ztgrad, zsgrad, zbgrad ! Variables used to calculate pycnocline 2234 ! ! gradients 2235 REAL(wp) :: zugrad, zvgrad ! Variables for calculating pycnocline shear 2236 REAL(wp) :: zdtdz_pyc ! Parametrized gradient of temperature in 2237 ! ! pycnocline 2238 REAL(wp) :: zdsdz_pyc ! Parametrised gradient of salinity in 2239 ! ! pycnocline 2240 REAL(wp) :: zdudz_pyc ! u-shear across the pycnocline 2241 REAL(wp) :: zdvdz_pyc ! v-shear across the pycnocline 2212 2242 ! 2213 2243 !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> … … 2231 2261 ! ------------------------------------------------------ 2232 2262 WHERE ( l_conv(A2D((nn_hls-1))) ) 2233 zsc_wth_1(:,:) = swstrl(A2D((nn_hls-1)))**3 * swth0(A2D((nn_hls-1))) / ( svstr(A2D((nn_hls-1)))**3 + 0.5_wp * swstrc(A2D((nn_hls-1)))**3 + epsln ) 2234 zsc_ws_1(:,:) = swstrl(A2D((nn_hls-1)))**3 * sws0(A2D((nn_hls-1))) / ( svstr(A2D((nn_hls-1)))**3 + 0.5_wp * swstrc(A2D((nn_hls-1)))**3 + epsln ) 2263 zsc_wth_1(:,:) = swstrl(A2D((nn_hls-1)))**3 * swth0(A2D((nn_hls-1))) / & 2264 & ( svstr(A2D((nn_hls-1)))**3 + 0.5_wp * swstrc(A2D((nn_hls-1)))**3 + epsln ) 2265 zsc_ws_1(:,:) = swstrl(A2D((nn_hls-1)))**3 * sws0(A2D((nn_hls-1))) / & 2266 & ( svstr(A2D((nn_hls-1)))**3 + 0.5_wp * swstrc(A2D((nn_hls-1)))**3 + epsln ) 2235 2267 ELSEWHERE 2236 2268 zsc_wth_1(:,:) = 2.0_wp * swthav(A2D((nn_hls-1))) … … 2266 2298 ! --------------------------------------------------------------------- 2267 2299 WHERE ( l_conv(A2D((nn_hls-1))) ) 2268 zsc_uw_1(:,:) = ( swstrl(A2D((nn_hls-1)))**3 + 0.5_wp * swstrc(A2D((nn_hls-1)))**3 )**pthird * sustke(A2D((nn_hls-1))) / & 2269 & MAX( ( 1.0_wp - 1.0_wp * 6.5_wp * sla(A2D((nn_hls-1)))**( 8.0_wp / 3.0_wp ) ), 0.2_wp ) 2270 zsc_uw_2(:,:) = ( swstrl(A2D((nn_hls-1)))**3 + 0.5_wp * swstrc(A2D((nn_hls-1)))**3 )**pthird * sustke(A2D((nn_hls-1))) / & 2271 & MIN( sla(A2D((nn_hls-1)))**( 8.0_wp / 3.0_wp ) + epsln, 0.12_wp ) 2272 zsc_vw_1(:,:) = ff_t(A2D((nn_hls-1))) * phml(A2D((nn_hls-1))) * sustke(A2D((nn_hls-1)))**3 * MIN( sla(A2D((nn_hls-1)))**( 8.0_wp / 3.0_wp ), 0.12_wp ) / & 2300 zsc_uw_1(:,:) = ( swstrl(A2D((nn_hls-1)))**3 + & 2301 & 0.5_wp * swstrc(A2D((nn_hls-1)))**3 )**pthird * sustke(A2D((nn_hls-1))) / & 2302 & MAX( ( 1.0_wp - 1.0_wp * 6.5_wp * sla(A2D((nn_hls-1)))**( 8.0_wp / 3.0_wp ) ), 0.2_wp ) 2303 zsc_uw_2(:,:) = ( swstrl(A2D((nn_hls-1)))**3 + & 2304 & 0.5_wp * swstrc(A2D((nn_hls-1)))**3 )**pthird * sustke(A2D((nn_hls-1))) / & 2305 & MIN( sla(A2D((nn_hls-1)))**( 8.0_wp / 3.0_wp ) + epsln, 0.12_wp ) 2306 zsc_vw_1(:,:) = ff_t(A2D((nn_hls-1))) * phml(A2D((nn_hls-1))) * sustke(A2D((nn_hls-1)))**3 * & 2307 & MIN( sla(A2D((nn_hls-1)))**( 8.0_wp / 3.0_wp ), 0.12_wp ) / & 2273 2308 & ( ( svstr(A2D((nn_hls-1)))**3 + 0.5_wp * swstrc(A2D((nn_hls-1)))**3 )**( 2.0_wp / 3.0_wp ) + epsln ) 2274 2309 ELSEWHERE 2275 2310 zsc_uw_1(:,:) = sustar(A2D((nn_hls-1)))**2 2276 zsc_vw_1(:,:) = ff_t(A2D((nn_hls-1))) * phbl(A2D((nn_hls-1))) * sustke(A2D((nn_hls-1)))**3 * MIN( sla(A2D((nn_hls-1)))**( 8.0_wp / 3.0_wp ), 0.12_wp ) /&2277 & ( svstr(A2D((nn_hls-1)))**2 + epsln )2311 zsc_vw_1(:,:) = ff_t(A2D((nn_hls-1))) * phbl(A2D((nn_hls-1))) * sustke(A2D((nn_hls-1)))**3 * & 2312 & MIN( sla(A2D((nn_hls-1)))**( 8.0_wp / 3.0_wp ), 0.12_wp ) / ( svstr(A2D((nn_hls-1)))**2 + epsln ) 2278 2313 ENDWHERE 2279 2314 DO_3D_OVR( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 2, MAX( jkm_mld, jkm_bld ) ) … … 2300 2335 ! ---------------------------------------------------------------------- 2301 2336 WHERE ( l_conv(A2D((nn_hls-1))) ) 2302 zsc_wth_1(:,:) = swbav(A2D((nn_hls-1))) * swth0(A2D((nn_hls-1))) * ( 1.0_wp + EXP( 0.2_wp * shol(A2D((nn_hls-1))) ) ) * phml(A2D((nn_hls-1))) /&2303 & ( svstr(A2D((nn_hls-1)))**3 + 0.5_wp * swstrc(A2D((nn_hls-1)))**3 + epsln )2304 zsc_ws_1(:,:) = swbav(A2D((nn_hls-1))) * sws0(A2D((nn_hls-1))) * ( 1.0_wp + EXP( 0.2_wp * shol(A2D((nn_hls-1))) ) ) * phml(A2D((nn_hls-1))) /&2305 & ( svstr(A2D((nn_hls-1)))**3 + 0.5_wp * swstrc(A2D((nn_hls-1)))**3 + epsln )2337 zsc_wth_1(:,:) = swbav(A2D((nn_hls-1))) * swth0(A2D((nn_hls-1))) * ( 1.0_wp + EXP( 0.2_wp * shol(A2D((nn_hls-1))) ) ) * & 2338 & phml(A2D((nn_hls-1))) / ( svstr(A2D((nn_hls-1)))**3 + 0.5_wp * swstrc(A2D((nn_hls-1)))**3 + epsln ) 2339 zsc_ws_1(:,:) = swbav(A2D((nn_hls-1))) * sws0(A2D((nn_hls-1))) * ( 1.0_wp + EXP( 0.2_wp * shol(A2D((nn_hls-1))) ) ) * & 2340 & phml(A2D((nn_hls-1))) / ( svstr(A2D((nn_hls-1)))**3 + 0.5_wp * swstrc(A2D((nn_hls-1)))**3 + epsln ) 2306 2341 ELSEWHERE 2307 2342 zsc_wth_1(:,:) = 0.0_wp … … 2443 2478 zsc_ws_1(:,:) = sws0(A2D((nn_hls-1))) / ( 1.0_wp - 0.56_wp * EXP( shol(A2D((nn_hls-1))) ) ) 2444 2479 WHERE ( l_pyc(A2D((nn_hls-1))) ) ! Pycnocline scales 2445 zsc_wth_pyc(:,:) = -0.003_wp * swstrc(A2D((nn_hls-1))) * ( 1.0_wp - pdh(A2D((nn_hls-1))) / phbl(A2D((nn_hls-1))) ) * av_dt_ml(A2D((nn_hls-1))) 2446 zsc_ws_pyc(:,:) = -0.003_wp * swstrc(A2D((nn_hls-1))) * ( 1.0_wp - pdh(A2D((nn_hls-1))) / phbl(A2D((nn_hls-1))) ) * av_ds_ml(A2D((nn_hls-1))) 2480 zsc_wth_pyc(:,:) = -0.003_wp * swstrc(A2D((nn_hls-1))) * ( 1.0_wp - pdh(A2D((nn_hls-1))) / phbl(A2D((nn_hls-1))) ) * & 2481 & av_dt_ml(A2D((nn_hls-1))) 2482 zsc_ws_pyc(:,:) = -0.003_wp * swstrc(A2D((nn_hls-1))) * ( 1.0_wp - pdh(A2D((nn_hls-1))) / phbl(A2D((nn_hls-1))) ) * & 2483 & av_ds_ml(A2D((nn_hls-1))) 2447 2484 END WHERE 2448 2485 ELSEWHERE … … 2454 2491 IF ( ( jk > 1 ) .AND. ( jk <= nmld(ji,jj) ) ) THEN 2455 2492 zznd_ml = gdepw(ji,jj,jk,Kmm) / phml(ji,jj) 2456 ghamt(ji,jj,jk) = ghamt(ji,jj,jk) + 0.3_wp * zsc_wth_1(ji,jj) * ( -2.0_wp + 2.75_wp * ( ( 1.0_wp + 0.6_wp * zznd_ml**4 ) - EXP( -6.0_wp * zznd_ml ) ) ) * & 2457 & ( 1.0_wp - EXP( -15.0_wp * ( 1.0_wp - zznd_ml ) ) ) 2458 ghams(ji,jj,jk) = ghams(ji,jj,jk) + 0.3_wp * zsc_ws_1(ji,jj) * ( -2.0_wp + 2.75_wp * ( ( 1.0_wp + 0.6_wp * zznd_ml**4 ) - EXP( -6.0_wp * zznd_ml ) ) ) * & 2459 & ( 1.0_wp - EXP( -15.0_wp * ( 1.0_wp - zznd_ml ) ) ) 2493 ghamt(ji,jj,jk) = ghamt(ji,jj,jk) + 0.3_wp * zsc_wth_1(ji,jj) * & 2494 & ( -2.0_wp + 2.75_wp * ( ( 1.0_wp + 0.6_wp * zznd_ml**4 ) - & 2495 & EXP( -6.0_wp * zznd_ml ) ) ) * & 2496 & ( 1.0_wp - EXP( -15.0_wp * ( 1.0_wp - zznd_ml ) ) ) 2497 ghams(ji,jj,jk) = ghams(ji,jj,jk) + 0.3_wp * zsc_ws_1(ji,jj) * & 2498 & ( -2.0_wp + 2.75_wp * ( ( 1.0_wp + 0.6_wp * zznd_ml**4 ) - & 2499 & EXP( -6.0_wp * zznd_ml ) ) ) * ( 1.0_wp - EXP( -15.0_wp * ( 1.0_wp - zznd_ml ) ) ) 2460 2500 END IF 2461 2501 ! … … 2463 2503 IF ( l_pyc(ji,jj) .AND. ( jk >= nmld(ji,jj) ) .AND. ( jk <= nbld(ji,jj) ) ) THEN ! Pycnocline 2464 2504 zznd_pyc = -1.0_wp * ( gdepw(ji,jj,jk,Kmm) - phbl(ji,jj) ) / pdh(ji,jj) 2465 ghamt(ji,jj,jk) = ghamt(ji,jj,jk) + 4.0_wp * zsc_wth_pyc(ji,jj) * ( 0.48_wp - EXP( -1.5_wp * ( zznd_pyc - 0.3_wp )**2 ) ) 2466 ghams(ji,jj,jk) = ghams(ji,jj,jk) + 4.0_wp * zsc_ws_pyc(ji,jj) * ( 0.48_wp - EXP( -1.5_wp * ( zznd_pyc - 0.3_wp )**2 ) ) 2505 ghamt(ji,jj,jk) = ghamt(ji,jj,jk) + 4.0_wp * zsc_wth_pyc(ji,jj) * & 2506 & ( 0.48_wp - EXP( -1.5_wp * ( zznd_pyc - 0.3_wp )**2 ) ) 2507 ghams(ji,jj,jk) = ghams(ji,jj,jk) + 4.0_wp * zsc_ws_pyc(ji,jj) * & 2508 & ( 0.48_wp - EXP( -1.5_wp * ( zznd_pyc - 0.3_wp )**2 ) ) 2467 2509 END IF 2468 2510 ELSE … … 2517 2559 ghamv(ji,jj,jk) = ghamv(ji,jj,jk) + 0.3_wp * 0.15_wp * SIN( 3.14159_wp * ( 0.65_wp * zznd_d ) ) * & 2518 2560 & EXP( -0.25_wp * zznd_d**2 ) * zsc_vw_1(ji,jj) 2519 ghamv(ji,jj,jk) = ghamv(ji,jj,jk) + 0.3_wp * 0.15_wp * EXP( -5.0 * ( 1.0 - znd ) ) * ( 1.0 - EXP( -20.0 * ( 1.0 - znd ) ) ) * zsc_vw_2(ji,jj) 2561 ghamv(ji,jj,jk) = ghamv(ji,jj,jk) + 0.3_wp * 0.15_wp * EXP( -5.0 * ( 1.0 - znd ) ) * & 2562 & ( 1.0 - EXP( -20.0 * ( 1.0 - znd ) ) ) * zsc_vw_2(ji,jj) 2520 2563 END IF 2521 2564 END IF … … 2617 2660 ENDIF ! IF (shol >=0.5) 2618 2661 ENDIF ! IF (av_db_bl> 0.) 2619 ENDIF ! IF (zdhdt >= 0) zdhdt < 0 not considered since pycnocline profile is zero and profile arrays are intialized to zero 2662 ENDIF ! IF (zdhdt >= 0) zdhdt < 0 not considered since pycnocline profile is zero and profile arrays are 2663 ! ! intialized to zero 2620 2664 END IF 2621 2665 END IF … … 2689 2733 !! 2690 2734 !!---------------------------------------------------------------------- 2691 INTEGER, INTENT(in ) :: Kmm ! Time-level index2735 INTEGER, INTENT(in ) :: Kmm ! Time-level index 2692 2736 REAL(wp), DIMENSION(A2D(nn_hls)), INTENT( out) :: pmld ! == Estimated FK BLD used for MLE horizontal gradients == ! 2693 2737 REAL(wp), DIMENSION(A2D(nn_hls)), INTENT(inout) :: pdtdx ! Horizontal gradient for Fox-Kemper parametrization … … 2698 2742 ! 2699 2743 ! Local variables 2700 INTEGER :: ji, jj, jk ! Dummy loop indices2701 INTEGER, DIMENSION(A2D(nn_hls)) :: jk_mld_prof ! Base level of MLE layer2702 INTEGER :: ikt, ikmax ! Local integers2703 REAL(wp) :: zc2704 REAL(wp) :: zN2_c ! Local buoyancy difference from 10m value2744 INTEGER :: ji, jj, jk ! Dummy loop indices 2745 INTEGER, DIMENSION(A2D(nn_hls)) :: jk_mld_prof ! Base level of MLE layer 2746 INTEGER :: ikt, ikmax ! Local integers 2747 REAL(wp) :: zc 2748 REAL(wp) :: zN2_c ! Local buoyancy difference from 10m value 2705 2749 REAL(wp), DIMENSION(A2D(nn_hls)) :: ztm 2706 2750 REAL(wp), DIMENSION(A2D(nn_hls)) :: zsm … … 2735 2779 zsm(:,:) = 0.0_wp 2736 2780 DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, ikmax ) 2737 zc = e3t(ji,jj,jk,Kmm) * REAL( MIN( MAX( 0, jk_mld_prof(ji,jj) - jk ), 1 ), KIND=wp ) ! zc being 0 outside the ML t-points 2781 zc = e3t(ji,jj,jk,Kmm) * REAL( MIN( MAX( 0, jk_mld_prof(ji,jj) - jk ), 1 ), KIND=wp ) ! zc being 0 outside the ML 2782 ! ! t-points 2738 2783 ztm(ji,jj) = ztm(ji,jj) + zc * ts(ji,jj,jk,jp_tem,Kmm) 2739 2784 zsm(ji,jj) = zsm(ji,jj) + zc * ts(ji,jj,jk,jp_sal,Kmm) … … 2796 2841 !!---------------------------------------------------------------------- 2797 2842 ! Outputs 2798 INTEGER, INTENT(in ) :: Kmm ! Time-level index2843 INTEGER, INTENT(in ) :: Kmm ! Time-level index 2799 2844 REAL(wp), DIMENSION(A2D((nn_hls-1))), INTENT(inout) :: pwb_fk 2800 2845 REAL(wp), DIMENSION(A2D((nn_hls-1))), INTENT(in ) :: phbl ! BL depth … … 2804 2849 ! 2805 2850 ! Local variables 2806 INTEGER :: ji, jj, jk ! Dummy loop indices2851 INTEGER :: ji, jj, jk ! Dummy loop indices 2807 2852 REAL(wp), DIMENSION(A2D((nn_hls-1))) :: znd_param 2808 REAL(wp) :: zthermal, zbeta2809 REAL(wp) :: zbuoy2810 REAL(wp) :: ztmp2811 REAL(wp) :: zpe_mle_layer2812 REAL(wp) :: zpe_mle_ref2813 REAL(wp) :: zdbdz_mle_int2853 REAL(wp) :: zthermal, zbeta 2854 REAL(wp) :: zbuoy 2855 REAL(wp) :: ztmp 2856 REAL(wp) :: zpe_mle_layer 2857 REAL(wp) :: zpe_mle_ref 2858 REAL(wp) :: zdbdz_mle_int 2814 2859 ! 2815 2860 znd_param(:,:) = 0.0_wp … … 2908 2953 !! 2909 2954 !!---------------------------------------------------------------------- 2910 INTEGER, INTENT(in ) :: Kmm ! Time-level index2955 INTEGER, INTENT(in ) :: Kmm ! Time-level index 2911 2956 REAL(wp), DIMENSION(A2D((nn_hls-1))), INTENT(in ) :: pmld ! == Estimated FK BLD used for MLE horiz gradients == ! 2912 2957 REAL(wp), DIMENSION(A2D((nn_hls-1))), INTENT(inout) :: phmle ! MLE depth
Note: See TracChangeset
for help on using the changeset viewer.
