- Timestamp:
- 2017-10-18T19:14:32+02:00 (7 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfiwm.F90
r8568 r8637 33 33 PRIVATE 34 34 35 PUBLIC zdf_iwm ! called in step module 36 PUBLIC zdf_iwm_init ! called in nemogcm module 37 38 ! !!* Namelist namzdf_iwm : internal wave-driven mixing * 39 INTEGER :: nn_zpyc ! pycnocline-intensified mixing energy proportional to N (=1) or N^2 (=2) 40 LOGICAL :: ln_mevar ! variable (=T) or constant (=F) mixing efficiency 41 LOGICAL :: ln_tsdiff ! account for differential T/S wave-driven mixing (=T) or not (=F) 42 43 REAL(wp) :: r1_6 = 1._wp / 6._wp 44 45 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: ebot_iwm ! power available from high-mode wave breaking (W/m2) 46 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: epyc_iwm ! power available from low-mode, pycnocline-intensified wave breaking (W/m2) 47 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: ecri_iwm ! power available from low-mode, critical slope wave breaking (W/m2) 48 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: hbot_iwm ! WKB decay scale for high-mode energy dissipation (m) 49 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: hcri_iwm ! decay scale for low-mode critical slope dissipation (m) 50 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: emix_iwm ! local energy density available for mixing (W/kg) 51 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: bflx_iwm ! buoyancy flux Kz * N^2 (W/kg) 52 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: pcmap_iwm ! vertically integrated buoyancy flux (W/m2) 53 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: zav_ratio ! S/T diffusivity ratio (only for ln_tsdiff=T) 54 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: zav_wave ! Internal wave-induced diffusivity 35 PUBLIC zdf_iwm ! called in step module 36 PUBLIC zdf_iwm_init ! called in nemogcm module 37 38 ! !!* Namelist namzdf_iwm : internal wave-driven mixing * 39 INTEGER :: nn_zpyc ! pycnocline-intensified mixing energy proportional to N (=1) or N^2 (=2) 40 LOGICAL :: ln_mevar ! variable (=T) or constant (=F) mixing efficiency 41 LOGICAL :: ln_tsdiff ! account for differential T/S wave-driven mixing (=T) or not (=F) 42 43 REAL(wp):: r1_6 = 1._wp / 6._wp 44 45 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: ebot_iwm ! power available from high-mode wave breaking (W/m2) 46 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: epyc_iwm ! power available from low-mode, pycnocline-intensified wave breaking (W/m2) 47 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: ecri_iwm ! power available from low-mode, critical slope wave breaking (W/m2) 48 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: hbot_iwm ! WKB decay scale for high-mode energy dissipation (m) 49 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: hcri_iwm ! decay scale for low-mode critical slope dissipation (m) 55 50 56 51 !! * Substitutions … … 67 62 !! *** FUNCTION zdf_iwm_alloc *** 68 63 !!---------------------------------------------------------------------- 69 ALLOCATE( ebot_iwm(jpi,jpj), epyc_iwm(jpi,jpj), ecri_iwm(jpi,jpj) , & 70 & hbot_iwm(jpi,jpj), hcri_iwm(jpi,jpj), emix_iwm(jpi,jpj,jpk), & 71 & bflx_iwm(jpi,jpj,jpk), pcmap_iwm(jpi,jpj), zav_ratio(jpi,jpj,jpk), & 72 & zav_wave(jpi,jpj,jpk), STAT=zdf_iwm_alloc ) 64 ALLOCATE( ebot_iwm(jpi,jpj), epyc_iwm(jpi,jpj), ecri_iwm(jpi,jpj) , & 65 & hbot_iwm(jpi,jpj), hcri_iwm(jpi,jpj) , STAT=zdf_iwm_alloc ) 73 66 ! 74 67 IF( lk_mpp ) CALL mpp_sum ( zdf_iwm_alloc ) … … 85 78 !! 86 79 !! ** Method : - internal wave-driven vertical mixing is given by: 87 !! Kz_wave = min( 100 cm2/s, f( Reb = emix_iwm /( Nu * N^2 ) )88 !! where emix_iwm is the 3D space distribution of the wave-breaking80 !! Kz_wave = min( 100 cm2/s, f( Reb = zemx_iwm /( Nu * N^2 ) ) 81 !! where zemx_iwm is the 3D space distribution of the wave-breaking 89 82 !! energy and Nu the molecular kinematic viscosity. 90 83 !! The function f(Reb) is linear (constant mixing efficiency) 91 84 !! if the namelist parameter ln_mevar = F and nonlinear if ln_mevar = T. 92 85 !! 93 !! - Compute emix_iwm, the 3D power density that allows to compute86 !! - Compute zemx_iwm, the 3D power density that allows to compute 94 87 !! Reb and therefrom the wave-induced vertical diffusivity. 95 88 !! This is divided into three components: 96 89 !! 1. Bottom-intensified low-mode dissipation at critical slopes 97 !! emix_iwm(z) = ( ecri_iwm / rau0 ) * EXP( -(H-z)/hcri_iwm )90 !! zemx_iwm(z) = ( ecri_iwm / rau0 ) * EXP( -(H-z)/hcri_iwm ) 98 91 !! / ( 1. - EXP( - H/hcri_iwm ) ) * hcri_iwm 99 92 !! where hcri_iwm is the characteristic length scale of the bottom 100 93 !! intensification, ecri_iwm a map of available power, and H the ocean depth. 101 94 !! 2. Pycnocline-intensified low-mode dissipation 102 !! emix_iwm(z) = ( epyc_iwm / rau0 ) * ( sqrt(rn2(z))^nn_zpyc )95 !! zemx_iwm(z) = ( epyc_iwm / rau0 ) * ( sqrt(rn2(z))^nn_zpyc ) 103 96 !! / SUM( sqrt(rn2(z))^nn_zpyc * e3w(z) ) 104 97 !! where epyc_iwm is a map of available power, and nn_zpyc … … 106 99 !! energy dissipation. 107 100 !! 3. WKB-height dependent high mode dissipation 108 !! emix_iwm(z) = ( ebot_iwm / rau0 ) * rn2(z) * EXP(-z_wkb(z)/hbot_iwm)101 !! zemx_iwm(z) = ( ebot_iwm / rau0 ) * rn2(z) * EXP(-z_wkb(z)/hbot_iwm) 109 102 !! / SUM( rn2(z) * EXP(-z_wkb(z)/hbot_iwm) * e3w(z) ) 110 103 !! where hbot_iwm is the characteristic length scale of the WKB bottom … … 121 114 !! avs = avt + av_wave * diffusivity_ratio(Reb) 122 115 !! 123 !! ** Action : - Define emix_iwm used to compute internal wave-induced mixing 124 !! - avt, avs, avm, increased by internal wave-driven mixing 116 !! ** Action : - avt, avs, avm, increased by tide internal wave-driven mixing 125 117 !! 126 118 !! References : de Lavergne et al. 2015, JPO; 2016, in prep. … … 132 124 INTEGER :: ji, jj, jk ! dummy loop indices 133 125 REAL(wp) :: zztmp ! scalar workspace 134 REAL(wp), DIMENSION(jpi,jpj) :: zfact ! Used for vertical structure 135 REAL(wp), DIMENSION(jpi,jpj) :: zhdep ! Ocean depth 136 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zwkb ! WKB-stretched height above bottom 137 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zweight ! Weight for high mode vertical distribution 138 REAL(wp), DIMENSION(jpi,jpj,jpk) :: znu_t ! Molecular kinematic viscosity (T grid) 139 REAL(wp), DIMENSION(jpi,jpj,jpk) :: znu_w ! Molecular kinematic viscosity (W grid) 140 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zReb ! Turbulence intensity parameter 126 REAL(wp), DIMENSION(jpi,jpj) :: zfact ! Used for vertical structure 127 REAL(wp), DIMENSION(jpi,jpj) :: zhdep ! Ocean depth 128 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zwkb ! WKB-stretched height above bottom 129 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zweight ! Weight for high mode vertical distribution 130 REAL(wp), DIMENSION(jpi,jpj,jpk) :: znu_t ! Molecular kinematic viscosity (T grid) 131 REAL(wp), DIMENSION(jpi,jpj,jpk) :: znu_w ! Molecular kinematic viscosity (W grid) 132 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zReb ! Turbulence intensity parameter 133 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zemx_iwm ! local energy density available for mixing (W/kg) 134 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zav_ratio ! S/T diffusivity ratio (only for ln_tsdiff=T) 135 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zav_wave ! Internal wave-induced diffusivity 136 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: z3d ! 3D workspace used for iom_put 137 REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: z2d ! 2D - - - - 141 138 !!---------------------------------------------------------------------- 142 139 ! 143 140 IF( ln_timing ) CALL timing_start('zdf_iwm') 144 141 ! 145 ! ! ----------------------------- ! 146 ! ! Internal wave-driven mixing ! (compute zav_wave) 147 ! ! ----------------------------- ! 142 ! !* Set to zero the 1st and last vertical levels of appropriate variables 143 zemx_iwm (:,:,1) = 0._wp ; zemx_iwm (:,:,jpk) = 0._wp 144 zav_ratio(:,:,1) = 0._wp ; zav_ratio(:,:,jpk) = 0._wp 145 zav_wave (:,:,1) = 0._wp ; zav_wave (:,:,jpk) = 0._wp 146 ! 147 ! ! ----------------------------- ! 148 ! ! Internal wave-driven mixing ! (compute zav_wave) 149 ! ! ----------------------------- ! 148 150 ! 149 ! 151 ! !* Critical slope mixing: distribute energy over the time-varying ocean depth, 150 152 ! using an exponential decay from the seafloor. 151 153 DO jj = 1, jpj ! part independent of the level … … 156 158 END DO 157 159 END DO 158 160 !!gm gde3w ==>>> check for ssh taken into account.... seem OK gde3w_n=gdept_n - sshn 159 161 DO jk = 2, jpkm1 ! complete with the level-dependent part 160 emix_iwm(:,:,jk) = zfact(:,:) * ( EXP( ( gde3w_n(:,:,jk ) - zhdep(:,:) ) / hcri_iwm(:,:) ) &162 zemx_iwm(:,:,jk) = zfact(:,:) * ( EXP( ( gde3w_n(:,:,jk ) - zhdep(:,:) ) / hcri_iwm(:,:) ) & 161 163 & - EXP( ( gde3w_n(:,:,jk-1) - zhdep(:,:) ) / hcri_iwm(:,:) ) ) * wmask(:,:,jk) & 162 164 & / ( gde3w_n(:,:,jk) - gde3w_n(:,:,jk-1) ) … … 186 188 ! 187 189 DO jk = 2, jpkm1 ! complete with the level-dependent part 188 emix_iwm(:,:,jk) = emix_iwm(:,:,jk) + zfact(:,:) * SQRT( MAX( 0._wp, rn2(:,:,jk) ) ) * wmask(:,:,jk)190 zemx_iwm(:,:,jk) = zemx_iwm(:,:,jk) + zfact(:,:) * SQRT( MAX( 0._wp, rn2(:,:,jk) ) ) * wmask(:,:,jk) 189 191 END DO 190 192 ! … … 203 205 ! 204 206 DO jk = 2, jpkm1 ! complete with the level-dependent part 205 emix_iwm(:,:,jk) = emix_iwm(:,:,jk) + zfact(:,:) * MAX( 0._wp, rn2(:,:,jk) ) * wmask(:,:,jk)207 zemx_iwm(:,:,jk) = zemx_iwm(:,:,jk) + zfact(:,:) * MAX( 0._wp, rn2(:,:,jk) ) * wmask(:,:,jk) 206 208 END DO 207 209 ! … … 253 255 ! 254 256 DO jk = 2, jpkm1 ! complete with the level-dependent part 255 emix_iwm(:,:,jk) = emix_iwm(:,:,jk) + zweight(:,:,jk) * zfact(:,:) * wmask(:,:,jk) &257 zemx_iwm(:,:,jk) = zemx_iwm(:,:,jk) + zweight(:,:,jk) * zfact(:,:) * wmask(:,:,jk) & 256 258 & / ( gde3w_n(:,:,jk) - gde3w_n(:,:,jk-1) ) 257 259 !!gm use of e3t_n just above? … … 269 271 ! Calculate turbulence intensity parameter Reb 270 272 DO jk = 2, jpkm1 271 zReb(:,:,jk) = emix_iwm(:,:,jk) / MAX( 1.e-20_wp, znu_w(:,:,jk) * rn2(:,:,jk) )273 zReb(:,:,jk) = zemx_iwm(:,:,jk) / MAX( 1.e-20_wp, znu_w(:,:,jk) * rn2(:,:,jk) ) 272 274 END DO 273 275 ! … … 349 351 ! !* output internal wave-driven mixing coefficient 350 352 CALL iom_put( "av_wave", zav_wave ) 351 !* output useful diagnostics: N^2, Kz * N^2 (bflx_iwm), 352 ! vertical integral of rau0 * Kz * N^2 (pcmap_iwm), energy density (emix_iwm) 353 !* output useful diagnostics: Kz*N^2 , 354 !!gm Kz*N2 should take into account the ratio avs/avt if it is used.... (see diaar5) 355 ! vertical integral of rau0 * Kz * N^2 , energy density (zemx_iwm) 353 356 IF( iom_use("bflx_iwm") .OR. iom_use("pcmap_iwm") ) THEN 354 bflx_iwm(:,:,:) = MAX( 0._wp, rn2(:,:,:) ) * zav_wave(:,:,:) 355 pcmap_iwm(:,:) = 0._wp 357 ALLOCATE( z2d(jpi,jpj) , z3d(jpi,jpj,jpk) ) 358 z3d(:,:,:) = MAX( 0._wp, rn2(:,:,:) ) * zav_wave(:,:,:) 359 z2d(:,:) = 0._wp 356 360 DO jk = 2, jpkm1 357 pcmap_iwm(:,:) = pcmap_iwm(:,:) + e3w_n(:,:,jk) * bflx_iwm(:,:,jk) * wmask(:,:,jk)358 END DO 359 pcmap_iwm(:,:) = rau0 * pcmap_iwm(:,:)360 CALL iom_put( "bflx_iwm", bflx_iwm)361 CALL iom_put( "pcmap_iwm", pcmap_iwm)362 ENDIF363 CALL iom_put( "bn2", rn2 )364 CALL iom_put( "emix_iwm", emix_iwm )361 z2d(:,:) = z2d(:,:) + e3w_n(:,:,jk) * z3d(:,:,jk) * wmask(:,:,jk) 362 END DO 363 z2d(:,:) = rau0 * z2d(:,:) 364 CALL iom_put( "bflx_iwm", z3d ) 365 CALL iom_put( "pcmap_iwm", z2d ) 366 DEALLOCATE( z2d , z3d ) 367 ENDIF 368 CALL iom_put( "emix_iwm", zemx_iwm ) 365 369 366 370 IF(ln_ctl) CALL prt_ctl(tab3d_1=zav_wave , clinfo1=' iwm - av_wave: ', tab3d_2=avt, clinfo2=' avt: ', ovlap=1, kdim=jpk) … … 466 470 ecri_iwm(:,:) = ecri_iwm(:,:) * ssmask(:,:) 467 471 468 ! Set once for all to zero the first and last vertical levels of appropriate variables469 emix_iwm (:,:, 1 ) = 0._wp470 emix_iwm (:,:,jpk) = 0._wp471 zav_ratio(:,:, 1 ) = 0._wp472 zav_ratio(:,:,jpk) = 0._wp473 zav_wave (:,:, 1 ) = 0._wp474 zav_wave (:,:,jpk) = 0._wp475 476 472 zbot = glob_sum( e1e2t(:,:) * ebot_iwm(:,:) ) 477 473 zpyc = glob_sum( e1e2t(:,:) * epyc_iwm(:,:) )
Note: See TracChangeset
for help on using the changeset viewer.
