Changeset 12488
 Timestamp:
 20200228T16:52:11+01:00 (3 years ago)
 Location:
 NEMO/branches/2020/KERNEL03_Storkey_Coward_RK3_stage2/tests/BENCH/MY_SRC
 Files:

 1 deleted
 1 edited
Legend:
 Unmodified
 Added
 Removed

NEMO/branches/2020/KERNEL03_Storkey_Coward_RK3_stage2/tests/BENCH/MY_SRC/zdfiwm.F90
r12443 r12488 48 48 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: hcri_iwm ! decay scale for lowmode critical slope dissipation (m) 49 49 50 !! * Substitutions 51 # include "do_loop_substitute.h90" 50 52 !! 51 53 !! NEMO/OCE 4.0 , NEMO Consortium (2018) 52 !! $Id: zdfiwm.F90 8093 20170530 08:13:14Z gm$53 !! Software governed by the CeCILL licen ce (./LICENSE)54 !! $Id: zdfiwm.F90 12377 20200212 14:39:06Z acc $ 55 !! Software governed by the CeCILL license (see ./LICENSE) 54 56 !! 55 57 CONTAINS … … 121 123 ! 122 124 INTEGER :: ji, jj, jk ! dummy loop indices 123 REAL(wp) :: zztmp ! scalar workspace125 REAL(wp) :: zztmp, ztmp1, ztmp2 ! scalar workspace 124 126 REAL(wp), DIMENSION(jpi,jpj) :: zfact ! Used for vertical structure 125 127 REAL(wp), DIMENSION(jpi,jpj) :: zhdep ! Ocean depth … … 147 149 ! !* Critical slope mixing: distribute energy over the timevarying ocean depth, 148 150 ! using an exponential decay from the seafloor. 149 DO jj = 1, jpj ! part independent of the level 150 DO ji = 1, jpi 151 zhdep(ji,jj) = gdepw_0(ji,jj,mbkt(ji,jj)+1) ! depth of the ocean 152 zfact(ji,jj) = rho0 * ( 1._wp  EXP( zhdep(ji,jj) / hcri_iwm(ji,jj) ) ) 153 IF( zfact(ji,jj) /= 0._wp ) zfact(ji,jj) = ecri_iwm(ji,jj) / zfact(ji,jj) 154 END DO 155 END DO 156 !!gm gde3w ==>>> check for ssh taken into account.... seem OK gde3w_n=gdept(Kmm)  ssh(Kmm) 157 DO jk = 2, jpkm1 ! complete with the leveldependent part 158 zemx_iwm(:,:,jk) = zfact(:,:) * ( EXP( ( gde3w(:,:,jk )  zhdep(:,:) ) / hcri_iwm(:,:) ) & 159 &  EXP( ( gde3w(:,:,jk1)  zhdep(:,:) ) / hcri_iwm(:,:) ) ) * wmask(:,:,jk) & 160 & / ( gde3w(:,:,jk)  gde3w(:,:,jk1) ) 161 162 !!gm delta(gde3w) = e3t(Kmm) !! Please verify the gridpoint position w versus tpoint 151 DO_2D_11_11 152 zhdep(ji,jj) = gdepw_0(ji,jj,mbkt(ji,jj)+1) ! depth of the ocean 153 zfact(ji,jj) = rho0 * ( 1._wp  EXP( zhdep(ji,jj) / hcri_iwm(ji,jj) ) ) 154 IF( zfact(ji,jj) /= 0._wp ) zfact(ji,jj) = ecri_iwm(ji,jj) / zfact(ji,jj) 155 END_2D 156 !!gm gde3w ==>>> check for ssh taken into account.... seem OK gde3w_n=gdept(:,:,:,Kmm)  ssh(:,:,Kmm) 157 DO_3D_11_11( 2, jpkm1 ) 158 IF ( zfact(ji,jj) == 0._wp .OR. wmask(ji,jj,jk) == 0._wp ) THEN ! optimization 159 zemx_iwm(ji,jj,jk) = 0._wp 160 ELSE 161 zemx_iwm(ji,jj,jk) = zfact(ji,jj) * ( EXP( ( gde3w(ji,jj,jk )  zhdep(ji,jj) ) / hcri_iwm(ji,jj) ) & 162 &  EXP( ( gde3w(ji,jj,jk1)  zhdep(ji,jj) ) / hcri_iwm(ji,jj) ) ) & 163 & / ( gde3w(ji,jj,jk)  gde3w(ji,jj,jk1) ) 164 ENDIF 165 END_3D 166 !!gm delta(gde3w) = e3t(:,:,:,Kmm) !! Please verify the gridpoint position w versus tpoint 163 167 !!gm it seems to me that only 1/hcri_iwm is used ==> compute it one for all 164 168 165 END DO166 169 167 170 ! !* Pycnoclineintensified mixing: distribute energy over the timevarying … … 177 180 END DO 178 181 ! 179 DO jj = 1, jpj 180 DO ji = 1, jpi 181 IF( zfact(ji,jj) /= 0 ) zfact(ji,jj) = epyc_iwm(ji,jj) / ( rho0 * zfact(ji,jj) ) 182 END DO 183 END DO 182 DO_2D_11_11 183 IF( zfact(ji,jj) /= 0 ) zfact(ji,jj) = epyc_iwm(ji,jj) / ( rho0 * zfact(ji,jj) ) 184 END_2D 184 185 ! 185 186 DO jk = 2, jpkm1 ! complete with the leveldependent part … … 194 195 END DO 195 196 ! 196 DO jj= 1, jpj 197 DO ji = 1, jpi 198 IF( zfact(ji,jj) /= 0 ) zfact(ji,jj) = epyc_iwm(ji,jj) / ( rho0 * zfact(ji,jj) ) 199 END DO 200 END DO 197 DO_2D_11_11 198 IF( zfact(ji,jj) /= 0 ) zfact(ji,jj) = epyc_iwm(ji,jj) / ( rho0 * zfact(ji,jj) ) 199 END_2D 201 200 ! 202 201 DO jk = 2, jpkm1 ! complete with the leveldependent part … … 223 222 !!gm 224 223 ! 225 DO jk = 2, jpkm1 226 DO jj = 1, jpj 227 DO ji = 1, jpi 228 IF( zfact(ji,jj) /= 0 ) zwkb(ji,jj,jk) = zhdep(ji,jj) * ( zfact(ji,jj)  zwkb(ji,jj,jk) ) & 229 & * wmask(ji,jj,jk) / zfact(ji,jj) 230 END DO 231 END DO 232 END DO 224 DO_3D_11_11( 2, jpkm1 ) 225 IF( zfact(ji,jj) /= 0 ) zwkb(ji,jj,jk) = zhdep(ji,jj) * ( zfact(ji,jj)  zwkb(ji,jj,jk) ) & 226 & * wmask(ji,jj,jk) / zfact(ji,jj) 227 END_3D 233 228 zwkb(:,:,1) = zhdep(:,:) * wmask(:,:,1) 234 229 ! 235 zweight(:,:,:) = 0._wp 236 DO jk = 2, jpkm1 237 zweight(:,:,jk) = MAX( 0._wp, rn2(:,:,jk) ) * hbot_iwm(:,:) * wmask(:,:,jk) & 238 & * ( EXP( zwkb(:,:,jk) / hbot_iwm(:,:) )  EXP( zwkb(:,:,jk1) / hbot_iwm(:,:) ) ) 239 END DO 230 DO_3D_11_11( 2, jpkm1 ) 231 IF ( rn2(ji,jj,jk) <= 0._wp .OR. wmask(ji,jj,jk) == 0._wp ) THEN ! optimization 232 zweight(ji,jj,jk) = 0._wp 233 ELSE 234 zweight(ji,jj,jk) = rn2(ji,jj,jk) * hbot_iwm(ji,jj) & 235 & * ( EXP( zwkb(ji,jj,jk) / hbot_iwm(ji,jj) )  EXP( zwkb(ji,jj,jk1) / hbot_iwm(ji,jj) ) ) 236 ENDIF 237 END_3D 240 238 ! 241 239 zfact(:,:) = 0._wp … … 244 242 END DO 245 243 ! 246 DO jj = 1, jpj 247 DO ji = 1, jpi 248 IF( zfact(ji,jj) /= 0 ) zfact(ji,jj) = ebot_iwm(ji,jj) / ( rho0 * zfact(ji,jj) ) 249 END DO 250 END DO 244 DO_2D_11_11 245 IF( zfact(ji,jj) /= 0 ) zfact(ji,jj) = ebot_iwm(ji,jj) / ( rho0 * zfact(ji,jj) ) 246 END_2D 251 247 ! 252 248 DO jk = 2, jpkm1 ! complete with the leveldependent part … … 276 272 ! 277 273 IF( ln_mevar ) THEN ! Variable mixing efficiency case : modify zav_wave in the 278 DO jk = 2, jpkm1 ! energetic (Reb > 480) and buoyancycontrolled (Reb <10.224 ) regimes 279 DO jj = 1, jpj 280 DO ji = 1, jpi 281 IF( zReb(ji,jj,jk) > 480.00_wp ) THEN 282 zav_wave(ji,jj,jk) = 3.6515_wp * znu_w(ji,jj,jk) * SQRT( zReb(ji,jj,jk) ) 283 ELSEIF( zReb(ji,jj,jk) < 10.224_wp ) THEN 284 zav_wave(ji,jj,jk) = 0.052125_wp * znu_w(ji,jj,jk) * zReb(ji,jj,jk) * SQRT( zReb(ji,jj,jk) ) 285 ENDIF 286 END DO 287 END DO 288 END DO 274 DO_3D_11_11( 2, jpkm1 ) 275 IF( zReb(ji,jj,jk) > 480.00_wp ) THEN 276 zav_wave(ji,jj,jk) = 3.6515_wp * znu_w(ji,jj,jk) * SQRT( zReb(ji,jj,jk) ) 277 ELSEIF( zReb(ji,jj,jk) < 10.224_wp ) THEN 278 zav_wave(ji,jj,jk) = 0.052125_wp * znu_w(ji,jj,jk) * zReb(ji,jj,jk) * SQRT( zReb(ji,jj,jk) ) 279 ENDIF 280 END_3D 289 281 ENDIF 290 282 ! … … 296 288 zztmp = 0._wp 297 289 !!gm used of glosum 3D.... 298 DO jk = 2, jpkm1 299 DO jj = 1, jpj 300 DO ji = 1, jpi 301 zztmp = zztmp + e3w(ji,jj,jk,Kmm) * e1e2t(ji,jj) & 302 & * MAX( 0._wp, rn2(ji,jj,jk) ) * zav_wave(ji,jj,jk) * wmask(ji,jj,jk) * tmask_i(ji,jj) 303 END DO 304 END DO 305 END DO 290 DO_3D_11_11( 2, jpkm1 ) 291 zztmp = zztmp + e3w(ji,jj,jk,Kmm) * e1e2t(ji,jj) & 292 & * MAX( 0._wp, rn2(ji,jj,jk) ) * zav_wave(ji,jj,jk) * wmask(ji,jj,jk) * tmask_i(ji,jj) 293 END_3D 306 294 CALL mpp_sum( 'zdfiwm', zztmp ) 307 295 zztmp = rho0 * zztmp ! Global integral of rauo * Kz * N^2 = power contributing to mixing … … 321 309 ! 322 310 IF( ln_tsdiff ) THEN !* Option for differential mixing of salinity and temperature 323 DO jk = 2, jpkm1 ! Calculate S/T diffusivity ratio as a function of Reb324 DO jj = 1, jpj325 DO ji = 1, jpi326 zav_ratio(ji,jj,jk) = ( 0.505_wp + 0.495_wp * &327 & TANH( 0.92_wp * ( LOG10( MAX( 1.e20_wp, zReb(ji,jj,jk) * 5._wp * r1_6 ) )  0.60_wp ) ) &328 & ) * wmask(ji,jj,jk)329 END DO330 END DO331 END DO311 ztmp1 = 0.505_wp + 0.495_wp * TANH( 0.92_wp * ( LOG10( 1.e20_wp )  0.60_wp ) ) 312 DO_3D_11_11( 2, jpkm1 ) 313 ztmp2 = zReb(ji,jj,jk) * 5._wp * r1_6 314 IF ( ztmp2 > 1.e20_wp .AND. wmask(ji,jj,jk) == 1._wp ) THEN 315 zav_ratio(ji,jj,jk) = 0.505_wp + 0.495_wp * TANH( 0.92_wp * ( LOG10(ztmp2)  0.60_wp ) ) 316 ELSE 317 zav_ratio(ji,jj,jk) = ztmp1 * wmask(ji,jj,jk) 318 ENDIF 319 END_3D 332 320 CALL iom_put( "av_ratio", zav_ratio ) 333 321 DO jk = 2, jpkm1 !* update momentum & tracer diffusivity with wavedriven mixing … … 364 352 CALL iom_put( "emix_iwm", zemx_iwm ) 365 353 366 IF( ln_ctl) CALL prt_ctl(tab3d_1=zav_wave , clinfo1=' iwm  av_wave: ', tab3d_2=avt, clinfo2=' avt: ', kdim=jpk)354 IF(sn_cfctl%l_prtctl) CALL prt_ctl(tab3d_1=zav_wave , clinfo1=' iwm  av_wave: ', tab3d_2=avt, clinfo2=' avt: ', kdim=jpk) 367 355 ! 368 356 END SUBROUTINE zdf_iwm … … 395 383 !! de Lavergne et al. in prep., 2017 396 384 !! 397 INTEGER :: ji, jj, jk ! dummy loop indices398 385 INTEGER :: inum ! local integer 399 386 INTEGER :: ios
Note: See TracChangeset
for help on using the changeset viewer.