- Timestamp:
- 2017-05-30T10:13:14+02:00 (7 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfiwm.F90
r7990 r8093 35 35 PUBLIC zdf_iwm ! called in step module 36 36 PUBLIC zdf_iwm_init ! called in nemogcm module 37 PUBLIC zdf_iwm_alloc ! called in nemogcm module38 37 39 38 ! !!* Namelist namzdf_iwm : internal wave-driven mixing * … … 78 77 79 78 80 SUBROUTINE zdf_iwm( kt )79 SUBROUTINE zdf_iwm( kt, p_avm, p_avt, p_avs ) 81 80 !!---------------------------------------------------------------------- 82 81 !! *** ROUTINE zdf_iwm *** … … 127 126 !! References : de Lavergne et al. 2015, JPO; 2016, in prep. 128 127 !!---------------------------------------------------------------------- 129 INTEGER, INTENT(in) :: kt ! ocean time-step 128 INTEGER , INTENT(in ) :: kt ! ocean time step 129 REAL(wp), DIMENSION(:,:,:) , INTENT(inout) :: p_avm ! momentum Kz (w-points) 130 REAL(wp), DIMENSION(:,:,:) , INTENT(inout) :: p_avt, p_avs ! tracer Kz (w-points) 130 131 ! 131 132 INTEGER :: ji, jj, jk ! dummy loop indices 132 REAL(wp) :: z tpc! scalar workspace133 REAL(wp), DIMENSION( :,:) , POINTER:: zfact ! Used for vertical structure134 REAL(wp), DIMENSION( :,:) , POINTER:: zhdep ! Ocean depth135 REAL(wp), DIMENSION( :,:,:), POINTER:: zwkb ! WKB-stretched height above bottom136 REAL(wp), DIMENSION( :,:,:), POINTER:: zweight ! Weight for high mode vertical distribution137 REAL(wp), DIMENSION( :,:,:), POINTER:: znu_t ! Molecular kinematic viscosity (T grid)138 REAL(wp), DIMENSION( :,:,:), POINTER:: znu_w ! Molecular kinematic viscosity (W grid)139 REAL(wp), DIMENSION( :,:,:), POINTER:: zReb ! Turbulence intensity parameter133 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 140 141 !!---------------------------------------------------------------------- 141 142 ! 142 143 IF( nn_timing == 1 ) CALL timing_start('zdf_iwm') 143 144 ! 144 CALL wrk_alloc( jpi,jpj, zfact, zhdep ) 145 CALL wrk_alloc( jpi,jpj,jpk, zwkb, zweight, znu_t, znu_w, zReb ) 146 147 ! ! ----------------------------- ! 148 ! ! Internal wave-driven mixing ! (compute zav_wave) 149 ! ! ----------------------------- ! 145 ! ! ----------------------------- ! 146 ! ! Internal wave-driven mixing ! (compute zav_wave) 147 ! ! ----------------------------- ! 150 148 ! 151 149 ! !* Critical slope mixing: distribute energy over the time-varying ocean depth, … … 162 160 emix_iwm(:,:,jk) = zfact(:,:) * ( EXP( ( gde3w_n(:,:,jk ) - zhdep(:,:) ) / hcri_iwm(:,:) ) & 163 161 & - EXP( ( gde3w_n(:,:,jk-1) - zhdep(:,:) ) / hcri_iwm(:,:) ) ) * wmask(:,:,jk) & 162 & / ( gde3w_n(:,:,jk) - gde3w_n(:,:,jk-1) ) 163 164 164 !!gm delta(gde3w_n) = e3t_n !! Please verify the grid-point position w versus t-point 165 & / ( gde3w_n(:,:,jk) - gde3w_n(:,:,jk-1) ) 165 !!gm it seems to me that only 1/hcri_iwm is used ==> compute it one for all 166 166 167 END DO 167 168 168 169 ! !* Pycnocline-intensified mixing: distribute energy over the time-varying 169 170 ! !* ocean depth as proportional to sqrt(rn2)^nn_zpyc 170 ! 171 ! ! (NB: N2 is masked, so no use of wmask here) 171 172 SELECT CASE ( nn_zpyc ) 172 173 ! … … 210 211 ! !* ocean depth as proportional to rn2 * exp(-z_wkb/rn_hbot) 211 212 ! 212 zwkb (:,:,:) = 0._wp213 zfact(:,:) = 0._wp213 zwkb (:,:,:) = 0._wp 214 zfact(:,:) = 0._wp 214 215 DO jk = 2, jpkm1 215 216 zfact(:,:) = zfact(:,:) + e3w_n(:,:,jk) * SQRT( MAX( 0._wp, rn2(:,:,jk) ) ) * wmask(:,:,jk) 216 217 zwkb(:,:,jk) = zfact(:,:) 217 218 END DO 219 !!gm even better: 220 ! DO jk = 2, jpkm1 221 ! zwkb(:,:) = zwkb(:,:) + e3w_n(:,:,jk) * SQRT( MAX( 0._wp, rn2(:,:,jk) ) ) 222 ! END DO 223 ! zfact(:,:) = zwkb(:,:,jpkm1) 224 !!gm or just use zwkb(k=jpk-1) instead of zfact... 225 !!gm 218 226 ! 219 227 DO jk = 2, jpkm1 … … 221 229 DO ji = 1, jpi 222 230 IF( zfact(ji,jj) /= 0 ) zwkb(ji,jj,jk) = zhdep(ji,jj) * ( zfact(ji,jj) - zwkb(ji,jj,jk) ) & 223 & * tmask(ji,jj,jk) / zfact(ji,jj)224 END DO 225 END DO 226 END DO 227 zwkb(:,:,1) = zhdep(:,:) * tmask(:,:,1)231 & * wmask(ji,jj,jk) / zfact(ji,jj) 232 END DO 233 END DO 234 END DO 235 zwkb(:,:,1) = zhdep(:,:) * wmask(:,:,1) 228 236 ! 229 237 zweight(:,:,:) = 0._wp … … 247 255 emix_iwm(:,:,jk) = emix_iwm(:,:,jk) + zweight(:,:,jk) * zfact(:,:) * wmask(:,:,jk) & 248 256 & / ( gde3w_n(:,:,jk) - gde3w_n(:,:,jk-1) ) 249 END DO 250 ! 257 !!gm use of e3t_n just above? 258 END DO 259 ! 260 !!gm this is to be replaced by just a constant value znu=1.e-6 m2/s 251 261 ! Calculate molecular kinematic viscosity 252 262 znu_t(:,:,:) = 1.e-4_wp * ( 17.91_wp - 0.53810_wp * tsn(:,:,:,jp_tem) + 0.00694_wp * tsn(:,:,:,jp_tem) * tsn(:,:,:,jp_tem) & … … 255 265 znu_w(:,:,jk) = 0.5_wp * ( znu_t(:,:,jk-1) + znu_t(:,:,jk) ) * wmask(:,:,jk) 256 266 END DO 267 !!gm end 257 268 ! 258 269 ! Calculate turbulence intensity parameter Reb … … 285 296 ! 286 297 IF( kt == nit000 ) THEN !* Control print at first time-step: diagnose the energy consumed by zav_wave 287 z tpc= 0._wp298 zztmp = 0._wp 288 299 !!gm used of glosum 3D.... 289 300 DO jk = 2, jpkm1 290 301 DO jj = 1, jpj 291 302 DO ji = 1, jpi 292 z tpc = ztpc+ e3w_n(ji,jj,jk) * e1e2t(ji,jj) &293 & * MAX( 0._wp, rn2(ji,jj,jk) ) * zav_wave(ji,jj,jk) * wmask(ji,jj,jk) * tmask_i(ji,jj)303 zztmp = zztmp + e3w_n(ji,jj,jk) * e1e2t(ji,jj) & 304 & * MAX( 0._wp, rn2(ji,jj,jk) ) * zav_wave(ji,jj,jk) * wmask(ji,jj,jk) * tmask_i(ji,jj) 294 305 END DO 295 306 END DO 296 307 END DO 297 IF( lk_mpp ) CALL mpp_sum( z tpc)298 z tpc = rau0 * ztpc! Global integral of rauo * Kz * N^2 = power contributing to mixing308 IF( lk_mpp ) CALL mpp_sum( zztmp ) 309 zztmp = rau0 * zztmp ! Global integral of rauo * Kz * N^2 = power contributing to mixing 299 310 ! 300 311 IF(lwp) THEN … … 303 314 WRITE(numout,*) '~~~~~~~ ' 304 315 WRITE(numout,*) 305 WRITE(numout,*) ' Total power consumption by av_wave : ztpc = ', ztpc* 1.e-12_wp, 'TW'316 WRITE(numout,*) ' Total power consumption by av_wave = ', zztmp * 1.e-12_wp, 'TW' 306 317 ENDIF 307 318 ENDIF … … 323 334 CALL iom_put( "av_ratio", zav_ratio ) 324 335 DO jk = 2, jpkm1 !* update momentum & tracer diffusivity with wave-driven mixing 325 avs(:,:,jk) =avs(:,:,jk) + zav_wave(:,:,jk) * zav_ratio(:,:,jk)326 avt(:,:,jk) =avt(:,:,jk) + zav_wave(:,:,jk)327 avm(:,:,jk) =avm(:,:,jk) + zav_wave(:,:,jk)336 p_avs(:,:,jk) = p_avs(:,:,jk) + zav_wave(:,:,jk) * zav_ratio(:,:,jk) 337 p_avt(:,:,jk) = p_avt(:,:,jk) + zav_wave(:,:,jk) 338 p_avm(:,:,jk) = p_avm(:,:,jk) + zav_wave(:,:,jk) 328 339 END DO 329 340 ! 330 341 ELSE !* update momentum & tracer diffusivity with wave-driven mixing 331 342 DO jk = 2, jpkm1 332 avs(:,:,jk) =avs(:,:,jk) + zav_wave(:,:,jk)333 avt(:,:,jk) =avt(:,:,jk) + zav_wave(:,:,jk)334 avm(:,:,jk) =avm(:,:,jk) + zav_wave(:,:,jk)343 p_avs(:,:,jk) = p_avs(:,:,jk) + zav_wave(:,:,jk) 344 p_avt(:,:,jk) = p_avt(:,:,jk) + zav_wave(:,:,jk) 345 p_avm(:,:,jk) = p_avm(:,:,jk) + zav_wave(:,:,jk) 335 346 END DO 336 347 ENDIF … … 353 364 CALL iom_put( "emix_iwm", emix_iwm ) 354 365 355 CALL wrk_dealloc( jpi,jpj, zfact, zhdep )356 CALL wrk_dealloc( jpi,jpj,jpk, zwkb, zweight, znu_t, znu_w, zReb )357 358 366 IF(ln_ctl) CALL prt_ctl(tab3d_1=zav_wave , clinfo1=' iwm - av_wave: ', tab3d_2=avt, clinfo2=' avt: ', ovlap=1, kdim=jpk) 359 367 !
Note: See TracChangeset
for help on using the changeset viewer.
