| 1 | MODULE zdftke_crs |
|---|
| 2 | !!====================================================================== |
|---|
| 3 | !! *** MODULE zdftke *** |
|---|
| 4 | !! Ocean physics: vertical mixing coefficient computed from the tke |
|---|
| 5 | !! turbulent closure parameterization |
|---|
| 6 | !!===================================================================== |
|---|
| 7 | !! History : OPA ! 1991-03 (b. blanke) Original code |
|---|
| 8 | !! 7.0 ! 1991-11 (G. Madec) bug fix |
|---|
| 9 | !! 7.1 ! 1992-10 (G. Madec) new mixing length and eav |
|---|
| 10 | !! 7.2 ! 1993-03 (M. Guyon) symetrical conditions |
|---|
| 11 | !! 7.3 ! 1994-08 (G. Madec, M. Imbard) nn_pdl flag |
|---|
| 12 | !! 7.5 ! 1996-01 (G. Madec) s-coordinates |
|---|
| 13 | !! 8.0 ! 1997-07 (G. Madec) lbc |
|---|
| 14 | !! 8.1 ! 1999-01 (E. Stretta) new option for the mixing length |
|---|
| 15 | !! NEMO 1.0 ! 2002-06 (G. Madec) add tke_init routine |
|---|
| 16 | !! - ! 2004-10 (C. Ethe ) 1D configuration |
|---|
| 17 | !! 2.0 ! 2006-07 (S. Masson) distributed restart using iom |
|---|
| 18 | !! 3.0 ! 2008-05 (C. Ethe, G.Madec) : update TKE physics: |
|---|
| 19 | !! ! - tke penetration (wind steering) |
|---|
| 20 | !! ! - suface condition for tke & mixing length |
|---|
| 21 | !! ! - Langmuir cells |
|---|
| 22 | !! - ! 2008-05 (J.-M. Molines, G. Madec) 2D form of avtb |
|---|
| 23 | !! - ! 2008-06 (G. Madec) style + DOCTOR name for namelist parameters |
|---|
| 24 | !! - ! 2008-12 (G. Reffray) stable discretization of the production term |
|---|
| 25 | !! 3.2 ! 2009-06 (G. Madec, S. Masson) TKE restart compatible with key_cpl |
|---|
| 26 | !! ! + cleaning of the parameters + bugs correction |
|---|
| 27 | !! 3.3 ! 2010-10 (C. Ethe, G. Madec) reorganisation of initialisation phase |
|---|
| 28 | !!---------------------------------------------------------------------- |
|---|
| 29 | #if defined key_zdftke || defined key_esopa |
|---|
| 30 | !!---------------------------------------------------------------------- |
|---|
| 31 | !! 'key_zdftke' TKE vertical physics |
|---|
| 32 | !!---------------------------------------------------------------------- |
|---|
| 33 | !! zdf_tke : update momentum and tracer Kz from a tke scheme |
|---|
| 34 | !! tke_tke : tke time stepping: update tke at now time step (en) |
|---|
| 35 | !! tke_avn : compute mixing length scale and deduce avm and avt |
|---|
| 36 | !! zdf_tke_init : initialization, namelist read, and parameters control |
|---|
| 37 | !! tke_rst : read/write tke restart in ocean restart file |
|---|
| 38 | !!---------------------------------------------------------------------- |
|---|
| 39 | USE crs |
|---|
| 40 | USE zdf_oce , ONLY : avtb, avmb, avtb_2d |
|---|
| 41 | USE zdftke |
|---|
| 42 | USE crslbclnk |
|---|
| 43 | !USE oce ! ocean: dynamics and active tracers variables |
|---|
| 44 | USE phycst ! physical constants |
|---|
| 45 | !USE dom_oce ! domain: ocean |
|---|
| 46 | !USE domvvl ! domain: variable volume layer |
|---|
| 47 | !USE sbc_oce ! surface boundary condition: ocean |
|---|
| 48 | !USE zdf_oce ! vertical physics: ocean variables |
|---|
| 49 | !USE zdfmxl ! vertical physics: mixed layer |
|---|
| 50 | !USE lbclnk ! ocean lateral boundary conditions (or mpp link) |
|---|
| 51 | !USE prtctl ! Print control |
|---|
| 52 | !USE in_out_manager ! I/O manager |
|---|
| 53 | !USE iom ! I/O manager library |
|---|
| 54 | !USE lib_mpp ! MPP library |
|---|
| 55 | USE wrk_nemo ! work arrays |
|---|
| 56 | USE timing ! Timing |
|---|
| 57 | !USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) |
|---|
| 58 | |
|---|
| 59 | IMPLICIT NONE |
|---|
| 60 | PRIVATE |
|---|
| 61 | |
|---|
| 62 | PUBLIC tke_avn_crs |
|---|
| 63 | PUBLIC tke_avn_ini_crs |
|---|
| 64 | |
|---|
| 65 | !LOGICAL , PUBLIC, PARAMETER :: lk_zdftke = .TRUE. !: TKE vertical mixing flag |
|---|
| 66 | |
|---|
| 67 | ! !!** Namelist namzdf_tke ** |
|---|
| 68 | !LOGICAL :: ln_mxl0 ! mixing length scale surface value as function of wind stress or not |
|---|
| 69 | !INTEGER :: nn_mxl ! type of mixing length (=0/1/2/3) |
|---|
| 70 | !REAL(wp) :: rn_mxl0 ! surface min value of mixing length (kappa*z_o=0.4*0.1 m) [m] |
|---|
| 71 | !INTEGER :: nn_pdl ! Prandtl number or not (ratio avt/avm) (=0/1) |
|---|
| 72 | !REAL(wp) :: rn_ediff ! coefficient for avt: avt=rn_ediff*mxl*sqrt(e) |
|---|
| 73 | !REAL(wp) :: rn_ediss ! coefficient of the Kolmogoroff dissipation |
|---|
| 74 | !REAL(wp) :: rn_ebb ! coefficient of the surface input of tke |
|---|
| 75 | !REAL(wp) :: rn_emin ! minimum value of tke [m2/s2] |
|---|
| 76 | !REAL(wp) :: rn_emin0 ! surface minimum value of tke [m2/s2] |
|---|
| 77 | !REAL(wp) :: rn_bshear ! background shear (>0) currently a numerical threshold (do not change it) |
|---|
| 78 | !INTEGER :: nn_etau ! type of depth penetration of surface tke (=0/1/2/3) |
|---|
| 79 | !INTEGER :: nn_htau ! type of tke profile of penetration (=0/1) |
|---|
| 80 | !REAL(wp) :: rn_efr ! fraction of TKE surface value which penetrates in the ocean |
|---|
| 81 | !LOGICAL :: ln_lc ! Langmuir cells (LC) as a source term of TKE or not |
|---|
| 82 | !REAL(wp) :: rn_lc ! coef to compute vertical velocity of Langmuir cells |
|---|
| 83 | |
|---|
| 84 | !REAL(wp) :: ri_cri ! critic Richardson number (deduced from rn_ediff and rn_ediss values) |
|---|
| 85 | !REAL(wp) :: rmxl_min ! minimum mixing length value (deduced from rn_ediff and rn_emin values) [m] |
|---|
| 86 | !REAL(wp) :: rhftau_add = 1.e-3_wp ! add offset applied to HF part of taum (nn_etau=3) |
|---|
| 87 | !REAL(wp) :: rhftau_scl = 1.0_wp ! scale factor applied to HF part of taum (nn_etau=3) |
|---|
| 88 | |
|---|
| 89 | !REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: en !: now turbulent kinetic energy [m2/s2] |
|---|
| 90 | !REAL(wp) , ALLOCATABLE, SAVE, DIMENSION(:,:) :: htau ! depth of tke penetration (nn_htau) |
|---|
| 91 | !REAL(wp) , ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: dissl ! now mixing lenght of dissipation |
|---|
| 92 | !REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: avt_k , avm_k ! not enhanced Kz |
|---|
| 93 | !REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: avmu_k, avmv_k ! not enhanced Kz |
|---|
| 94 | #if defined key_c1d |
|---|
| 95 | ! !!** 1D cfg only ** ('key_c1d') |
|---|
| 96 | !REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: e_dis, e_mix !: dissipation and mixing turbulent lengh scales |
|---|
| 97 | !REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: e_pdl, e_ric !: prandl and local Richardson numbers |
|---|
| 98 | #endif |
|---|
| 99 | |
|---|
| 100 | !! * Substitutions |
|---|
| 101 | # include "domzgr_substitute.h90" |
|---|
| 102 | # include "vectopt_loop_substitute.h90" |
|---|
| 103 | !!---------------------------------------------------------------------- |
|---|
| 104 | !! NEMO/OPA 4.0 , NEMO Consortium (2011) |
|---|
| 105 | !! $Id: zdftke.F90 4990 2014-12-15 16:42:49Z timgraham $ |
|---|
| 106 | !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) |
|---|
| 107 | !!---------------------------------------------------------------------- |
|---|
| 108 | CONTAINS |
|---|
| 109 | |
|---|
| 110 | ! INTEGER FUNCTION zdf_tke_alloc() |
|---|
| 111 | ! !!---------------------------------------------------------------------- |
|---|
| 112 | ! !! *** FUNCTION zdf_tke_alloc *** |
|---|
| 113 | ! !!---------------------------------------------------------------------- |
|---|
| 114 | ! ALLOCATE( & |
|---|
| 115 | !#if defined key_c1d |
|---|
| 116 | ! & e_dis(jpi,jpj,jpk) , e_mix(jpi,jpj,jpk) , & |
|---|
| 117 | ! & e_pdl(jpi,jpj,jpk) , e_ric(jpi,jpj,jpk) , & |
|---|
| 118 | !#endif |
|---|
| 119 | ! & en (jpi,jpj,jpk) , htau (jpi,jpj) , dissl(jpi,jpj,jpk) , & |
|---|
| 120 | ! & avt_k (jpi,jpj,jpk) , avm_k (jpi,jpj,jpk), & |
|---|
| 121 | ! & avmu_k(jpi,jpj,jpk) , avmv_k(jpi,jpj,jpk), STAT= zdf_tke_alloc ) |
|---|
| 122 | ! ! |
|---|
| 123 | ! IF( lk_mpp ) CALL mpp_sum ( zdf_tke_alloc ) |
|---|
| 124 | ! IF( zdf_tke_alloc /= 0 ) CALL ctl_warn('zdf_tke_alloc: failed to allocate arrays') |
|---|
| 125 | ! ! |
|---|
| 126 | ! END FUNCTION zdf_tke_alloc |
|---|
| 127 | |
|---|
| 128 | |
|---|
| 129 | SUBROUTINE tke_avn_crs |
|---|
| 130 | !!---------------------------------------------------------------------- |
|---|
| 131 | !! *** ROUTINE tke_avn *** |
|---|
| 132 | !! |
|---|
| 133 | !! ** Purpose : Compute the vertical eddy viscosity and diffusivity |
|---|
| 134 | !! |
|---|
| 135 | !! ** Method : At this stage, en, the now TKE, is known (computed in |
|---|
| 136 | !! the tke_tke routine). First, the now mixing lenth is |
|---|
| 137 | !! computed from en and the strafification (N^2), then the mixings |
|---|
| 138 | !! coefficients are computed. |
|---|
| 139 | !! - Mixing length : a first evaluation of the mixing lengh |
|---|
| 140 | !! scales is: |
|---|
| 141 | !! mxl = sqrt(2*en) / N |
|---|
| 142 | !! where N is the brunt-vaisala frequency, with a minimum value set |
|---|
| 143 | !! to rmxl_min (rn_mxl0) in the interior (surface) ocean. |
|---|
| 144 | !! The mixing and dissipative length scale are bound as follow : |
|---|
| 145 | !! nn_mxl=0 : mxl bounded by the distance to surface and bottom. |
|---|
| 146 | !! zmxld = zmxlm = mxl |
|---|
| 147 | !! nn_mxl=1 : mxl bounded by the e3w and zmxld = zmxlm = mxl |
|---|
| 148 | !! nn_mxl=2 : mxl bounded such that the vertical derivative of mxl is |
|---|
| 149 | !! less than 1 (|d/dz(mxl)|<1) and zmxld = zmxlm = mxl |
|---|
| 150 | !! nn_mxl=3 : mxl is bounded from the surface to the bottom usings |
|---|
| 151 | !! |d/dz(xml)|<1 to obtain lup, and from the bottom to |
|---|
| 152 | !! the surface to obtain ldown. the resulting length |
|---|
| 153 | !! scales are: |
|---|
| 154 | !! zmxld = sqrt( lup * ldown ) |
|---|
| 155 | !! zmxlm = min ( lup , ldown ) |
|---|
| 156 | !! - Vertical eddy viscosity and diffusivity: |
|---|
| 157 | !! avm = max( avtb, rn_ediff * zmxlm * en^1/2 ) |
|---|
| 158 | !! avt = max( avmb, pdlr * avm ) |
|---|
| 159 | !! with pdlr=1 if nn_pdl=0, pdlr=1/pdl=F(Ri) otherwise. |
|---|
| 160 | !! |
|---|
| 161 | !! ** Action : - avt : now vertical eddy diffusivity (w-point) |
|---|
| 162 | !! - avmu, avmv : now vertical eddy viscosity at uw- and vw-points |
|---|
| 163 | !!---------------------------------------------------------------------- |
|---|
| 164 | INTEGER :: ji, jj, jk ! dummy loop indices |
|---|
| 165 | REAL(wp) :: zrn2, zraug, zcoef, zav ! local scalars |
|---|
| 166 | REAL(wp) :: zdku, zdkv, zpdlr, zri, zsqen ! - - |
|---|
| 167 | REAL(wp) :: zemlm, zemlp ! - - |
|---|
| 168 | !REAL(wp), POINTER, DIMENSION(:,:,:) :: zmpdl, zmxlm, zmxld |
|---|
| 169 | REAL(wp), POINTER, DIMENSION(:,:,:) :: zmxlm_crs,zmxld_crs |
|---|
| 170 | REAL(wp), POINTER, DIMENSION(:,:,:) :: avm_crs, avmu_crs,avmv_crs |
|---|
| 171 | !!-------------------------------------------------------------------- |
|---|
| 172 | ! |
|---|
| 173 | IF( nn_timing == 1 ) CALL timing_start('tke_avn') |
|---|
| 174 | |
|---|
| 175 | !CALL wrk_alloc( jpi,jpj,jpk, zmpdl, zmxlm, zmxld ) |
|---|
| 176 | CALL wrk_alloc( jpi_crs,jpj_crs,jpk, zmxlm_crs, zmxld_crs, avm_crs, avmu_crs, avmv_crs ) |
|---|
| 177 | |
|---|
| 178 | ! !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< |
|---|
| 179 | ! ! Mixing length |
|---|
| 180 | ! !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< |
|---|
| 181 | ! |
|---|
| 182 | ! !* Buoyancy length scale: l=sqrt(2*e/n**2) |
|---|
| 183 | ! |
|---|
| 184 | IF( ln_mxl0 ) THEN ! surface mixing length = F(stress) : l=vkarmn*2.e5*taum/(rau0*g) |
|---|
| 185 | DO jj = 2, jpj_crs-1 |
|---|
| 186 | DO ji = 2, jpi_crs-1 |
|---|
| 187 | !IF (mikt(ji,jj) .GT. 1) THEN |
|---|
| 188 | ! zmxlm(ji,jj,mikt(ji,jj)) = rmxl_min |
|---|
| 189 | !ELSE |
|---|
| 190 | zraug = vkarmn * 2.e5_wp / ( rau0 * grav ) |
|---|
| 191 | zmxlm_crs(ji,jj,1) = MAX( rn_mxl0, zraug * taum_crs(ji,jj) ) |
|---|
| 192 | !END IF |
|---|
| 193 | END DO |
|---|
| 194 | END DO |
|---|
| 195 | !ELSE |
|---|
| 196 | ! DO jj = 2, jpjm1 |
|---|
| 197 | ! DO ji = fs_2, fs_jpim1 ! surface set to the minimum value |
|---|
| 198 | ! zmxlm(ji,jj,mikt(ji,jj)) = MAX( tmask(ji,jj,1) * rn_mxl0, rmxl_min) |
|---|
| 199 | ! END DO |
|---|
| 200 | ! END DO |
|---|
| 201 | ENDIF |
|---|
| 202 | zmxlm_crs(:,:,jpk) = rmxl_min ! last level set to the interior minium value |
|---|
| 203 | ! |
|---|
| 204 | DO jj = 2, jpj_crs-1 |
|---|
| 205 | DO ji = 2, jpi_crs-1 ! vector opt. |
|---|
| 206 | DO jk = 2, jpk-1 ! interior value : l=sqrt(2*e/n^2) |
|---|
| 207 | zrn2 = MAX( rn2_crs(ji,jj,jk), rsmall ) |
|---|
| 208 | zmxlm_crs(ji,jj,jk) = MAX( rmxl_min, SQRT( 2._wp * en_crs(ji,jj,jk) / zrn2 ) ) |
|---|
| 209 | END DO |
|---|
| 210 | zmxld_crs(ji,jj,1) = zmxlm_crs(ji,jj,1) ! surface set to the minimum value |
|---|
| 211 | END DO |
|---|
| 212 | END DO |
|---|
| 213 | ! |
|---|
| 214 | ! !* Physical limits for the mixing length |
|---|
| 215 | ! |
|---|
| 216 | zmxld_crs(:,:, 1 ) = zmxlm_crs(:,:,1) ! surface set to the zmxlm value |
|---|
| 217 | zmxld_crs(:,:,jpk) = rmxl_min ! last level set to the minimum value |
|---|
| 218 | ! |
|---|
| 219 | SELECT CASE ( nn_mxl ) |
|---|
| 220 | ! |
|---|
| 221 | !CASE ( 0 ) ! bounded by the distance to surface and bottom |
|---|
| 222 | ! DO jj = 2, jpjm1 |
|---|
| 223 | ! DO ji = fs_2, fs_jpim1 ! vector opt. |
|---|
| 224 | ! DO jk = mikt(ji,jj)+1, jpkm1 |
|---|
| 225 | ! zemxl = MIN( fsdepw(ji,jj,jk) - fsdepw(ji,jj,mikt(ji,jj)), zmxlm(ji,jj,jk), & |
|---|
| 226 | ! & fsdepw(ji,jj,mbkt(ji,jj)+1) - fsdepw(ji,jj,jk) ) |
|---|
| 227 | ! zmxlm(ji,jj,jk) = zemxl |
|---|
| 228 | ! zmxld(ji,jj,jk) = zemxl |
|---|
| 229 | ! END DO |
|---|
| 230 | ! END DO |
|---|
| 231 | ! END DO |
|---|
| 232 | ! ! |
|---|
| 233 | !CASE ( 1 ) ! bounded by the vertical scale factor |
|---|
| 234 | ! DO jj = 2, jpjm1 |
|---|
| 235 | ! DO ji = fs_2, fs_jpim1 ! vector opt. |
|---|
| 236 | ! DO jk = mikt(ji,jj)+1, jpkm1 |
|---|
| 237 | ! zemxl = MIN( fse3w(ji,jj,jk), zmxlm(ji,jj,jk) ) |
|---|
| 238 | ! zmxlm(ji,jj,jk) = zemxl |
|---|
| 239 | ! zmxld(ji,jj,jk) = zemxl |
|---|
| 240 | ! END DO |
|---|
| 241 | ! END DO |
|---|
| 242 | ! END DO |
|---|
| 243 | ! ! |
|---|
| 244 | !CASE ( 2 ) ! |dk[xml]| bounded by e3t : |
|---|
| 245 | ! DO jj = 2, jpjm1 |
|---|
| 246 | ! DO ji = fs_2, fs_jpim1 ! vector opt. |
|---|
| 247 | ! DO jk = mikt(ji,jj)+1, jpkm1 ! from the surface to the bottom : |
|---|
| 248 | ! zmxlm(ji,jj,jk) = MIN( zmxlm(ji,jj,jk-1) + fse3t(ji,jj,jk-1), zmxlm(ji,jj,jk) ) |
|---|
| 249 | ! END DO |
|---|
| 250 | ! DO jk = jpkm1, mikt(ji,jj)+1, -1 ! from the bottom to the surface : |
|---|
| 251 | ! zemxl = MIN( zmxlm(ji,jj,jk+1) + fse3t(ji,jj,jk+1), zmxlm(ji,jj,jk) ) |
|---|
| 252 | ! zmxlm(ji,jj,jk) = zemxl |
|---|
| 253 | ! zmxld(ji,jj,jk) = zemxl |
|---|
| 254 | ! END DO |
|---|
| 255 | ! END DO |
|---|
| 256 | ! END DO |
|---|
| 257 | ! ! |
|---|
| 258 | CASE ( 3 ) ! lup and ldown, |dk[xml]| bounded by e3t : |
|---|
| 259 | DO jj = 2, jpj_crs-1 |
|---|
| 260 | DO ji = 2, jpi_crs-1 ! vector opt. |
|---|
| 261 | DO jk = 2, jpkm1 ! from the surface to the bottom : lup |
|---|
| 262 | zmxld_crs(ji,jj,jk) = MIN( zmxld_crs(ji,jj,jk-1) + e3t_crs(ji,jj,jk-1), zmxlm_crs(ji,jj,jk) ) |
|---|
| 263 | END DO |
|---|
| 264 | DO jk = jpkm1, 2 , -1 ! from the bottom to the surface : ldown |
|---|
| 265 | zmxlm_crs(ji,jj,jk) = MIN( zmxlm_crs(ji,jj,jk+1) + e3t_crs(ji,jj,jk+1), zmxlm_crs(ji,jj,jk) ) |
|---|
| 266 | END DO |
|---|
| 267 | END DO |
|---|
| 268 | END DO |
|---|
| 269 | DO jk = 2, jpkm1 |
|---|
| 270 | DO jj = 2, jpj_crs-1 |
|---|
| 271 | DO ji = 2, jpi_crs-1 ! vector opt. |
|---|
| 272 | zemlm = MIN ( zmxld_crs(ji,jj,jk), zmxlm_crs(ji,jj,jk) ) |
|---|
| 273 | zemlp = SQRT( zmxld_crs(ji,jj,jk) * zmxlm_crs(ji,jj,jk) ) |
|---|
| 274 | zmxlm_crs(ji,jj,jk) = zemlm |
|---|
| 275 | zmxld_crs(ji,jj,jk) = zemlp |
|---|
| 276 | END DO |
|---|
| 277 | END DO |
|---|
| 278 | END DO |
|---|
| 279 | ! |
|---|
| 280 | END SELECT |
|---|
| 281 | ! |
|---|
| 282 | ! !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< |
|---|
| 283 | ! ! Vertical eddy viscosity and diffusivity (avmu, avmv, avt) |
|---|
| 284 | ! !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< |
|---|
| 285 | DO jk = 1, jpkm1 !* vertical eddy viscosity & diffivity at w-points |
|---|
| 286 | DO jj = 2, jpj_crs-1 |
|---|
| 287 | DO ji = 2, jpi_crs-1 ! vector opt. |
|---|
| 288 | zsqen = SQRT( en_crs(ji,jj,jk) ) |
|---|
| 289 | zav = rn_ediff * zmxlm_crs(ji,jj,jk) * zsqen |
|---|
| 290 | avm_crs(ji,jj,jk) = MAX( zav, avmb(jk) ) * tmask_crs(ji,jj,jk) |
|---|
| 291 | avt_crs(ji,jj,jk) = MAX( zav, avtb_2d(ji,jj) * avtb(jk) ) * tmask_crs(ji,jj,jk) |
|---|
| 292 | !dissl(ji,jj,jk) = zsqen / zmxld(ji,jj,jk) |
|---|
| 293 | END DO |
|---|
| 294 | END DO |
|---|
| 295 | END DO |
|---|
| 296 | CALL crs_lbc_lnk( avm_crs, 'W', 1. ) ! Lateral boundary conditions (sign unchanged) |
|---|
| 297 | ! |
|---|
| 298 | DO jj = 2, jpj_crs-1 |
|---|
| 299 | DO ji = 2, jpi_crs-1 ! vector opt. |
|---|
| 300 | DO jk = 2, jpkm1 !* vertical eddy viscosity at u- and v-points |
|---|
| 301 | avmu_crs(ji,jj,jk) = 0.5 * ( avm_crs(ji,jj,jk) + avm_crs(ji+1,jj ,jk) ) * umask_crs(ji,jj,jk) |
|---|
| 302 | END DO |
|---|
| 303 | DO jk = 2, jpkm1 |
|---|
| 304 | avmv_crs(ji,jj,jk) = 0.5 * ( avm_crs(ji,jj,jk) + avm_crs(ji ,jj+1,jk) ) * vmask_crs(ji,jj,jk) |
|---|
| 305 | END DO |
|---|
| 306 | END DO |
|---|
| 307 | END DO |
|---|
| 308 | CALL crs_lbc_lnk( avmu_crs, 'U', 1. ) ; CALL crs_lbc_lnk( avmv_crs, 'V', 1. ) ! Lateral boundary conditions |
|---|
| 309 | ! |
|---|
| 310 | IF( nn_pdl == 1 ) THEN !* Prandtl number case: update avt |
|---|
| 311 | DO jj = 2, jpj_crs-1 |
|---|
| 312 | DO ji = 2, jpi_crs-1 ! vector opt. |
|---|
| 313 | DO jk = 2, jpkm1 |
|---|
| 314 | zcoef = avm_crs(ji,jj,jk) * 2._wp * e3w_crs(ji,jj,jk) * e3w_crs(ji,jj,jk) |
|---|
| 315 | ! ! shear |
|---|
| 316 | zdku = avmu_crs(ji-1,jj,jk) * ( un_crs(ji-1,jj,jk-1) - un_crs(ji-1,jj,jk) ) * ( ub_crs(ji-1,jj,jk-1) - ub_crs(ji-1,jj,jk) ) & |
|---|
| 317 | & + avmu_crs(ji ,jj,jk) * ( un_crs(ji ,jj,jk-1) - un_crs(ji ,jj,jk) ) * ( ub_crs(ji ,jj,jk-1) - ub_crs(ji ,jj,jk) ) |
|---|
| 318 | zdkv = avmv_crs(ji,jj-1,jk) * ( vn_crs(ji,jj-1,jk-1) - vn_crs(ji,jj-1,jk) ) * ( vb_crs(ji,jj-1,jk-1) - vb_crs(ji,jj-1,jk) ) & |
|---|
| 319 | & + avmv_crs(ji,jj ,jk) * ( vn_crs(ji,jj ,jk-1) - vn_crs(ji,jj ,jk) ) * ( vb_crs(ji,jj ,jk-1) - vb_crs(ji,jj ,jk) ) |
|---|
| 320 | ! ! local Richardson number |
|---|
| 321 | zri = MAX( rb2_crs(ji,jj,jk), 0._wp ) * zcoef / (zdku + zdkv + rn_bshear ) |
|---|
| 322 | zpdlr = MAX( 0.1_wp, 0.2 / MAX( 0.2 , zri ) ) |
|---|
| 323 | avt_crs(ji,jj,jk) = MAX( zpdlr * avt_crs(ji,jj,jk), avtb_2d(ji,jj) * avtb(jk) ) * tmask_crs(ji,jj,jk) |
|---|
| 324 | END DO |
|---|
| 325 | END DO |
|---|
| 326 | END DO |
|---|
| 327 | ENDIF |
|---|
| 328 | CALL crs_lbc_lnk( avt_crs, 'W', 1. ) ! Lateral boundary conditions on avt (sign unchanged) |
|---|
| 329 | |
|---|
| 330 | !IF(ln_ctl) THEN |
|---|
| 331 | ! CALL prt_ctl( tab3d_1=en , clinfo1=' tke - e: ', tab3d_2=avt, clinfo2=' t: ', ovlap=1, kdim=jpk) |
|---|
| 332 | ! CALL prt_ctl( tab3d_1=avmu, clinfo1=' tke - u: ', mask1=umask, & |
|---|
| 333 | ! & tab3d_2=avmv, clinfo2= ' v: ', mask2=vmask, ovlap=1, kdim=jpk ) |
|---|
| 334 | !ENDIF |
|---|
| 335 | ! |
|---|
| 336 | CALL wrk_dealloc( jpi_crs,jpj_crs,jpk, zmxlm_crs, zmxld_crs, avm_crs, avmu_crs, avmv_crs ) |
|---|
| 337 | ! |
|---|
| 338 | IF( nn_timing == 1 ) CALL timing_stop('tke_avn') |
|---|
| 339 | ! |
|---|
| 340 | END SUBROUTINE tke_avn_crs |
|---|
| 341 | |
|---|
| 342 | SUBROUTINE tke_avn_ini_crs |
|---|
| 343 | !!---------------------------------------------------------------------- |
|---|
| 344 | !! |
|---|
| 345 | !! |
|---|
| 346 | !! |
|---|
| 347 | !! |
|---|
| 348 | !!---------------------------------------------------------------------- |
|---|
| 349 | !IF( nn_avb == 0 ) THEN ! Define avmb, avtb from namelist parameter |
|---|
| 350 | ! avtb(:) = rn_avt0 |
|---|
| 351 | !ELSE ! Background profile of avt (fit a theoretical/observational profile (Krauss 1990) |
|---|
| 352 | ! avmb(:) = rn_avm0 |
|---|
| 353 | ! avtb(:) = rn_avt0 + ( 3.e-4_wp - 2._wp * rn_avt0 ) * 1.e-4_wp * gdepw_1d(:) ! m2/s |
|---|
| 354 | ! IF(ln_sco .AND. lwp) CALL ctl_warn( 'avtb profile not valid in sco' ) |
|---|
| 355 | !ENDIF |
|---|
| 356 | |
|---|
| 357 | |
|---|
| 358 | |
|---|
| 359 | avtb_2d_crs(:,:) = 1.e0 |
|---|
| 360 | !IF( nn_havtb == 1 ) THEN ! decrease avtb in the equatorial band |
|---|
| 361 | ! ! -15S -5S : linear decrease from avt0 to avt0/10. |
|---|
| 362 | ! ! -5S +5N : cst value avt0/10. |
|---|
| 363 | ! ! 5N 15N : linear increase from avt0/10, to avt0 |
|---|
| 364 | ! WHERE(-15. <= gphit .AND. gphit < -5 ) avtb_2d = (1. - 0.09 * (gphit + 15.)) |
|---|
| 365 | ! WHERE( -5. <= gphit .AND. gphit < 5 ) avtb_2d = 0.1 |
|---|
| 366 | ! WHERE( 5. <= gphit .AND. gphit < 15 ) avtb_2d = (0.1 + 0.09 * (gphit - 5.)) |
|---|
| 367 | !ENDIF |
|---|
| 368 | |
|---|
| 369 | |
|---|
| 370 | |
|---|
| 371 | END SUBROUTINE tke_avn_ini_crs |
|---|
| 372 | #else |
|---|
| 373 | !!---------------------------------------------------------------------- |
|---|
| 374 | !! Dummy module : NO TKE scheme |
|---|
| 375 | !!---------------------------------------------------------------------- |
|---|
| 376 | LOGICAL, PUBLIC, PARAMETER :: lk_zdftke = .FALSE. !: TKE flag |
|---|
| 377 | #endif |
|---|
| 378 | |
|---|
| 379 | !!====================================================================== |
|---|
| 380 | END MODULE zdftke_crs |
|---|
