- Timestamp:
- 2015-12-08T12:39:53+01:00 (8 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/UKMO/icebergs_restart_single_file/NEMOGCM/NEMO/OPA_SRC/ZDF/zdftmx.F90
r6019 r6020 8 8 !! 3.3 ! 2010-10 (C. Ethe, G. Madec) reorganisation of initialisation phase 9 9 !!---------------------------------------------------------------------- 10 #if defined key_zdftmx || defined key_esopa10 #if defined key_zdftmx 11 11 !!---------------------------------------------------------------------- 12 12 !! 'key_zdftmx' Tidal vertical mixing … … 54 54 # include "vectopt_loop_substitute.h90" 55 55 !!---------------------------------------------------------------------- 56 !! NEMO/OPA 4.0 , NEMO Consortium (2011)56 !! NEMO/OPA 3.7 , NEMO Consortium (2014) 57 57 !! $Id$ 58 58 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) … … 105 105 !! Koch-Larrouy et al. 2007, GRL. 106 106 !!---------------------------------------------------------------------- 107 USE oce, zav_tide => ua ! use ua as workspace108 !!109 107 INTEGER, INTENT(in) :: kt ! ocean time-step 110 ! !108 ! 111 109 INTEGER :: ji, jj, jk ! dummy loop indices 112 110 REAL(wp) :: ztpc ! scalar workspace 113 REAL(wp), POINTER, DIMENSION(:,:) :: zkz 111 REAL(wp), POINTER, DIMENSION(:,:) :: zkz 112 REAL(wp), POINTER, DIMENSION(:,:,:) :: zav_tide 114 113 !!---------------------------------------------------------------------- 115 114 ! 116 115 IF( nn_timing == 1 ) CALL timing_start('zdf_tmx') 117 116 ! 118 CALL wrk_alloc( jpi,jpj, zkz ) 119 117 CALL wrk_alloc( jpi,jpj, zkz ) 118 CALL wrk_alloc( jpi,jpj,jpk, zav_tide ) 119 ! 120 120 ! ! ----------------------- ! 121 121 ! ! Standard tidal mixing ! (compute zav_tide) … … 136 136 137 137 DO jk = 2, jpkm1 !* Mutiply by zkz to recover en_tmx, BUT bound by 30/6 ==> zav_tide bound by 300 cm2/s 138 DO jj = 1, jpj !* Here zkz should be equal to en_tmx ==> multiply by en_tmx/zkz to recover en_tmx 139 DO ji = 1, jpi 140 zav_tide(ji,jj,jk) = zav_tide(ji,jj,jk) * MIN( zkz(ji,jj), 30./6. ) * wmask(ji,jj,jk) !kz max = 300 cm2/s 141 END DO 142 END DO 138 zav_tide(:,:,jk) = zav_tide(:,:,jk) * MIN( zkz(:,:), 30./6. ) * wmask(:,:,jk) !kz max = 300 cm2/s 143 139 END DO 144 140 145 141 IF( kt == nit000 ) THEN !* check at first time-step: diagnose the energy consumed by zav_tide 146 ztpc = 0. e0142 ztpc = 0._wp 147 143 DO jk= 1, jpk 148 144 DO jj= 1, jpj 149 145 DO ji= 1, jpi 150 ztpc = ztpc + fse3w(ji,jj,jk) * e1 t(ji,jj) * e2t(ji,jj)&151 & 146 ztpc = ztpc + fse3w(ji,jj,jk) * e1e2t(ji,jj) & 147 & * MAX( 0.e0, rn2(ji,jj,jk) ) * zav_tide(ji,jj,jk) * tmask(ji,jj,jk) * tmask_i(ji,jj) 152 148 END DO 153 149 END DO 154 150 END DO 155 151 ztpc= rau0 / ( rn_tfe * rn_me ) * ztpc 152 IF( lk_mpp ) CALL mpp_sum( ztpc ) 156 153 IF(lwp) WRITE(numout,*) 157 154 IF(lwp) WRITE(numout,*) ' N Total power consumption by av_tide : ztpc = ', ztpc * 1.e-12 ,'TW' … … 167 164 ! ! ----------------------- ! 168 165 DO jk = 2, jpkm1 !* update momentum & tracer diffusivity with tidal mixing 169 DO jj = 1, jpj !* Here zkz should be equal to en_tmx ==> multiply by en_tmx/zkz to recover en_tmx 170 DO ji = 1, jpi 171 avt(ji,jj,jk) = avt(ji,jj,jk) + zav_tide(ji,jj,jk) * wmask(ji,jj,jk) 172 avm(ji,jj,jk) = avm(ji,jj,jk) + zav_tide(ji,jj,jk) * wmask(ji,jj,jk) 173 END DO 174 END DO 175 END DO 176 177 DO jk = 2, jpkm1 !* update momentum & tracer diffusivity with tidal mixing 166 avt(:,:,jk) = avt(:,:,jk) + zav_tide(:,:,jk) * wmask(:,:,jk) 167 avm(:,:,jk) = avm(:,:,jk) + zav_tide(:,:,jk) * wmask(:,:,jk) 178 168 DO jj = 2, jpjm1 179 169 DO ji = fs_2, fs_jpim1 ! vector opt. … … 190 180 IF(ln_ctl) CALL prt_ctl(tab3d_1=zav_tide , clinfo1=' tmx - av_tide: ', tab3d_2=avt, clinfo2=' avt: ', ovlap=1, kdim=jpk) 191 181 ! 192 CALL wrk_dealloc( jpi,jpj, zkz ) 182 CALL wrk_dealloc( jpi,jpj, zkz ) 183 CALL wrk_dealloc( jpi,jpj,jpk, zav_tide ) 193 184 ! 194 185 IF( nn_timing == 1 ) CALL timing_stop('zdf_tmx') … … 239 230 DO jk = 1, jpkm1 240 231 zdn2dz (:,:,jk) = rn2(:,:,jk) - rn2(:,:,jk+1) ! Vertical profile of dN2/dz 241 !CDIR NOVERRCHK242 232 zempba_3d_1(:,:,jk) = SQRT( MAX( 0.e0, rn2(:,:,jk) ) ) ! - - of N 243 233 zempba_3d_2(:,:,jk) = MAX( 0.e0, rn2(:,:,jk) ) ! - - of N^2 … … 248 238 zsum2(:,:) = 0.e0 249 239 DO jk= 2, jpk 250 zsum1(:,:) = zsum1(:,:) + zempba_3d_1(:,:,jk) * fse3w(:,:,jk) * tmask(:,:,jk) * tmask(:,:,jk-1)251 zsum2(:,:) = zsum2(:,:) + zempba_3d_2(:,:,jk) * fse3w(:,:,jk) * tmask(:,:,jk) * tmask(:,:,jk-1)240 zsum1(:,:) = zsum1(:,:) + zempba_3d_1(:,:,jk) * fse3w(:,:,jk) * wmask(:,:,jk) 241 zsum2(:,:) = zsum2(:,:) + zempba_3d_2(:,:,jk) * fse3w(:,:,jk) * wmask(:,:,jk) 252 242 END DO 253 243 DO jj = 1, jpj … … 285 275 zkz(:,:) = 0.e0 ! Associated potential energy consummed over the whole water column 286 276 DO jk = 2, jpkm1 287 zkz(:,:) = zkz(:,:) + fse3w(:,:,jk) * MAX( 0.e0, rn2(:,:,jk) ) * rau0 * zavt_itf(:,:,jk) * tmask(:,:,jk) * tmask(:,:,jk-1)277 zkz(:,:) = zkz(:,:) + fse3w(:,:,jk) * MAX( 0.e0, rn2(:,:,jk) ) * rau0 * zavt_itf(:,:,jk) * wmask(:,:,jk) 288 278 END DO 289 279 … … 295 285 296 286 DO jk = 2, jpkm1 ! Mutiply by zkz to recover en_tmx, BUT bound by 30/6 ==> zavt_itf bound by 300 cm2/s 297 zavt_itf(:,:,jk) = zavt_itf(:,:,jk) * MIN( zkz(:,:), 120./10. ) * tmask(:,:,jk) * tmask(:,:,jk-1) ! kz max = 120 cm2/s287 zavt_itf(:,:,jk) = zavt_itf(:,:,jk) * MIN( zkz(:,:), 120./10. ) * wmask(:,:,jk) ! kz max = 120 cm2/s 298 288 END DO 299 289 … … 303 293 DO jj= 1, jpj 304 294 DO ji= 1, jpi 305 ztpc = ztpc + e1 t(ji,jj) *e2t(ji,jj) * fse3w(ji,jj,jk) * MAX( 0.e0, rn2(ji,jj,jk) ) &306 & * zavt_itf(ji,jj,jk) * tmask(ji,jj,jk) * tmask_i(ji,jj)295 ztpc = ztpc + e1e2t(ji,jj) * fse3w(ji,jj,jk) * MAX( 0.e0, rn2(ji,jj,jk) ) & 296 & * zavt_itf(ji,jj,jk) * tmask(ji,jj,jk) * tmask_i(ji,jj) 307 297 END DO 308 298 END DO 309 299 END DO 300 IF( lk_mpp ) CALL mpp_sum( ztpc ) 310 301 ztpc= rau0 * ztpc / ( rn_me * rn_tfe_itf ) 311 302 IF(lwp) WRITE(numout,*) ' N Total power consumption by zavt_itf: ztpc = ', ztpc * 1.e-12 ,'TW' … … 361 352 !! Koch-Larrouy et al. 2007, GRL. 362 353 !!---------------------------------------------------------------------- 363 USE oce , zav_tide => ua ! ua used as workspace364 !!365 354 INTEGER :: ji, jj, jk ! dummy loop indices 366 355 INTEGER :: inum ! local integer 367 356 INTEGER :: ios 368 357 REAL(wp) :: ztpc, ze_z ! local scalars 369 REAL(wp), DIMENSION(:,:) , POINTER :: zem2, zek1 ! read M2 and K1 tidal energy370 REAL(wp), DIMENSION(:,:) , POINTER :: zkz ! total M2, K1 and S2 tidal energy371 REAL(wp), DIMENSION(:,:) , POINTER :: zfact ! used for vertical structure function372 REAL(wp), DIMENSION(:,:) , POINTER :: zhdep ! Ocean depth373 REAL(wp), DIMENSION(:,:,:), POINTER :: zpc 358 REAL(wp), DIMENSION(:,:) , POINTER :: zem2, zek1 ! read M2 and K1 tidal energy 359 REAL(wp), DIMENSION(:,:) , POINTER :: zkz ! total M2, K1 and S2 tidal energy 360 REAL(wp), DIMENSION(:,:) , POINTER :: zfact ! used for vertical structure function 361 REAL(wp), DIMENSION(:,:) , POINTER :: zhdep ! Ocean depth 362 REAL(wp), DIMENSION(:,:,:), POINTER :: zpc, zav_tide ! power consumption 374 363 !! 375 364 NAMELIST/namzdf_tmx/ rn_htmx, rn_n2min, rn_tfe, rn_me, ln_tmx_itf, rn_tfe_itf … … 378 367 IF( nn_timing == 1 ) CALL timing_start('zdf_tmx_init') 379 368 ! 380 CALL wrk_alloc( jpi,jpj, zem2, zek1, zkz, zfact, zhdep )381 CALL wrk_alloc( jpi,jpj,jpk, zpc)382 383 REWIND( numnam_ref ) 369 CALL wrk_alloc( jpi,jpj, zem2, zek1, zkz, zfact, zhdep ) 370 CALL wrk_alloc( jpi,jpj,jpk, zpc, zav_tide ) 371 ! 372 REWIND( numnam_ref ) ! Namelist namzdf_tmx in reference namelist : Tidal Mixing 384 373 READ ( numnam_ref, namzdf_tmx, IOSTAT = ios, ERR = 901) 385 374 901 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namzdf_tmx in reference namelist', lwp ) 386 387 REWIND( numnam_cfg ) 375 ! 376 REWIND( numnam_cfg ) ! Namelist namzdf_tmx in configuration namelist : Tidal Mixing 388 377 READ ( numnam_cfg, namzdf_tmx, IOSTAT = ios, ERR = 902 ) 389 378 902 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namzdf_tmx in configuration namelist', lwp ) 390 379 IF(lwm) WRITE ( numond, namzdf_tmx ) 391 392 IF(lwp) THEN ! Control print380 ! 381 IF(lwp) THEN ! Control print 393 382 WRITE(numout,*) 394 383 WRITE(numout,*) 'zdf_tmx_init : tidal mixing' … … 402 391 WRITE(numout,*) ' ITF tidal dissipation efficiency = ', rn_tfe_itf 403 392 ENDIF 404 405 ! ! allocate tmx arrays 393 ! ! allocate tmx arrays 406 394 IF( zdf_tmx_alloc() /= 0 ) CALL ctl_stop( 'STOP', 'zdf_tmx_init : unable to allocate tmx arrays' ) 407 395 408 IF( ln_tmx_itf ) THEN ! read the Indonesian Through Flow mask396 IF( ln_tmx_itf ) THEN ! read the Indonesian Through Flow mask 409 397 CALL iom_open('mask_itf',inum) 410 398 CALL iom_get (inum, jpdom_data, 'tmaskitf',mask_itf,1) ! 411 399 CALL iom_close(inum) 412 400 ENDIF 413 414 ! read M2 tidal energy flux : W/m2 ( zem2 < 0 ) 401 ! ! read M2 tidal energy flux : W/m2 ( zem2 < 0 ) 415 402 CALL iom_open('M2rowdrg',inum) 416 403 CALL iom_get (inum, jpdom_data, 'field',zem2,1) ! 417 404 CALL iom_close(inum) 418 419 ! read K1 tidal energy flux : W/m2 ( zek1 < 0 ) 405 ! ! read K1 tidal energy flux : W/m2 ( zek1 < 0 ) 420 406 CALL iom_open('K1rowdrg',inum) 421 407 CALL iom_get (inum, jpdom_data, 'field',zek1,1) ! 422 408 CALL iom_close(inum) 423 424 ! Total tidal energy ( M2, S2 and K1 with S2=(1/2)^2 * M2 ) 425 ! only the energy available for mixing is taken into account, 426 ! (mixing efficiency tidal dissipation efficiency) 409 ! ! Total tidal energy ( M2, S2 and K1 with S2=(1/2)^2 * M2 ) 410 ! ! only the energy available for mixing is taken into account, 411 ! ! (mixing efficiency tidal dissipation efficiency) 427 412 en_tmx(:,:) = - rn_tfe * rn_me * ( zem2(:,:) * 1.25 + zek1(:,:) ) * ssmask(:,:) 428 413 429 414 !============ 430 415 !TG: Bug for VVL? Should this section be moved out of _init and be updated at every timestep? 431 ! Vertical structure (az_tmx) 432 DO jj = 1, jpj ! part independent of the level 416 !!gm : you are right, but tidal mixing acts in deep ocean (H>500m) where e3 is O(100m) 417 !! the error is thus ~1% which I feel comfortable with, compared to uncertainties in tidal energy dissipation. 418 ! ! Vertical structure (az_tmx) 419 DO jj = 1, jpj ! part independent of the level 433 420 DO ji = 1, jpi 434 421 zhdep(ji,jj) = gdepw_0(ji,jj,mbkt(ji,jj)+1) ! depth of the ocean … … 437 424 END DO 438 425 END DO 439 DO jk= 1, jpk ! complete with the level-dependent part426 DO jk= 1, jpk ! complete with the level-dependent part 440 427 DO jj = 1, jpj 441 428 DO ji = 1, jpi … … 445 432 END DO 446 433 !=========== 447 434 ! 448 435 IF( nprint == 1 .AND. lwp ) THEN 449 436 ! Control print … … 454 441 zav_tide(:,:,jk) = az_tmx(:,:,jk) / MAX( rn_n2min, rn2(:,:,jk) ) 455 442 END DO 456 457 ztpc = 0. e0443 ! 444 ztpc = 0._wp 458 445 zpc(:,:,:) = MAX(rn_n2min,rn2(:,:,:)) * zav_tide(:,:,:) 459 446 DO jk= 2, jpkm1 460 447 DO jj = 1, jpj 461 448 DO ji = 1, jpi 462 ztpc = ztpc + fse3w(ji,jj,jk) * e1 t(ji,jj) *e2t(ji,jj) * zpc(ji,jj,jk) * wmask(ji,jj,jk) * tmask_i(ji,jj)449 ztpc = ztpc + fse3w(ji,jj,jk) * e1e2t(ji,jj) * zpc(ji,jj,jk) * wmask(ji,jj,jk) * tmask_i(ji,jj) 463 450 END DO 464 451 END DO 465 452 END DO 453 IF( lk_mpp ) CALL mpp_sum( ztpc ) 466 454 ztpc= rau0 * 1/(rn_tfe * rn_me) * ztpc 467 455 ! 468 456 WRITE(numout,*) 469 457 WRITE(numout,*) ' Total power consumption of the tidally driven part of Kz : ztpc = ', ztpc * 1.e-12 ,'TW' 470 471 458 ! 472 459 ! control print 2 473 460 zav_tide(:,:,:) = MIN( zav_tide(:,:,:), 60.e-4 ) 474 zkz(:,:) = 0. e0461 zkz(:,:) = 0._wp 475 462 DO jk = 2, jpkm1 476 DO jj = 1, jpj 477 DO ji = 1, jpi 478 zkz(ji,jj) = zkz(ji,jj) + fse3w(ji,jj,jk) * MAX(0.e0, rn2(ji,jj,jk)) * rau0 * zav_tide(ji,jj,jk) * wmask(ji,jj,jk) 479 END DO 480 END DO 463 zkz(:,:) = zkz(:,:) + fse3w(:,:,jk) * MAX(0.e0, rn2(:,:,jk)) * rau0 * zav_tide(:,:,jk) * wmask(:,:,jk) 481 464 END DO 482 465 ! Here zkz should be equal to en_tmx ==> multiply by en_tmx/zkz … … 497 480 END DO 498 481 WRITE(numout,*) ' Min de zkz ', ztpc, ' Max = ', maxval(zkz(:,:) ) 499 482 ! 500 483 DO jk = 2, jpkm1 501 DO jj = 1, jpj 502 DO ji = 1, jpi 503 zav_tide(ji,jj,jk) = zav_tide(ji,jj,jk) * MIN( zkz(ji,jj), 30./6. ) * wmask(ji,jj,jk) !kz max = 300 cm2/s 504 END DO 505 END DO 506 END DO 507 ztpc = 0.e0 484 zav_tide(:,:,jk) = zav_tide(:,:,jk) * MIN( zkz(:,:), 30./6. ) * wmask(:,:,jk) !kz max = 300 cm2/s 485 END DO 486 ztpc = 0._wp 508 487 zpc(:,:,:) = Max(0.e0,rn2(:,:,:)) * zav_tide(:,:,:) 509 488 DO jk= 1, jpk 510 489 DO jj = 1, jpj 511 490 DO ji = 1, jpi 512 ztpc = ztpc + fse3w(ji,jj,jk) * e1 t(ji,jj) *e2t(ji,jj) * zpc(ji,jj,jk) * wmask(ji,jj,jk) * tmask_i(ji,jj)491 ztpc = ztpc + fse3w(ji,jj,jk) * e1e2t(ji,jj) * zpc(ji,jj,jk) * wmask(ji,jj,jk) * tmask_i(ji,jj) 513 492 END DO 514 493 END DO 515 494 END DO 495 IF( lk_mpp ) CALL mpp_sum( ztpc ) 516 496 ztpc= rau0 * 1/(rn_tfe * rn_me) * ztpc 517 497 WRITE(numout,*) ' 2 Total power consumption of the tidally driven part of Kz : ztpc = ', ztpc * 1.e-12 ,'TW' 518 498 !!gm bug mpp in these diagnostics 519 499 DO jk = 1, jpk 520 ze_z = SUM( e1 t(:,:) * e2t(:,:) * zav_tide(:,:,jk)* tmask_i(:,:) ) &521 & / MAX( 1.e-20, SUM( e1 t(:,:) * e2t(:,:) * wmask(:,:,jk) * tmask_i(:,:) ) )522 ztpc = 1. E50500 ze_z = SUM( e1e2t(:,:) * zav_tide(:,:,jk) * tmask_i(:,:) ) & 501 & / MAX( 1.e-20, SUM( e1e2t(:,:) * wmask (:,:,jk) * tmask_i(:,:) ) ) 502 ztpc = 1.e50 523 503 DO jj = 1, jpj 524 504 DO ji = 1, jpi 525 IF( zav_tide(ji,jj,jk) /= 0.e0 ) ztpc = Min( ztpc, zav_tide(ji,jj,jk) )505 IF( zav_tide(ji,jj,jk) /= 0.e0 ) ztpc = MIN( ztpc, zav_tide(ji,jj,jk) ) 526 506 END DO 527 507 END DO … … 530 510 END DO 531 511 532 WRITE(numout,*) ' e_tide : ', SUM( e1 t*e2t*en_tmx ) / ( rn_tfe * rn_me ) * 1.e-12, 'TW'512 WRITE(numout,*) ' e_tide : ', SUM( e1e2t*en_tmx ) / ( rn_tfe * rn_me ) * 1.e-12, 'TW' 533 513 WRITE(numout,*) 534 514 WRITE(numout,*) ' Initial profile of tidal vertical mixing' … … 539 519 END DO 540 520 END DO 541 ze_z = SUM( e1 t(:,:) * e2t(:,:) * zkz(:,:)* tmask_i(:,:) ) &542 & / MAX( 1.e-20, SUM( e1 t(:,:) * e2t(:,:) * wmask(:,:,jk) * tmask_i(:,:) ) )521 ze_z = SUM( e1e2t(:,:) * zkz (:,:) * tmask_i(:,:) ) & 522 & / MAX( 1.e-20, SUM( e1e2t(:,:) * wmask(:,:,jk) * tmask_i(:,:) ) ) 543 523 WRITE(numout,*) ' jk= ', jk,' ', ze_z * 1.e4,' cm2/s' 544 524 END DO 545 525 DO jk = 1, jpk 546 526 zkz(:,:) = az_tmx(:,:,jk) /rn_n2min 547 ze_z = SUM( e1 t(:,:) * e2t(:,:) * zkz(:,:)* tmask_i(:,:) ) &548 & / MAX( 1.e-20, SUM( e1 t(:,:) * e2t(:,:) * wmask(:,:,jk) * tmask_i(:,:) ) )527 ze_z = SUM( e1e2t(:,:) * zkz (:,:) * tmask_i(:,:) ) & 528 & / MAX( 1.e-20, SUM( e1e2t(:,:) * wmask(:,:,jk) * tmask_i(:,:) ) ) 549 529 WRITE(numout,*) 550 530 WRITE(numout,*) ' N2 min - jk= ', jk,' ', ze_z * 1.e4,' cm2/s min= ',MINVAL(zkz)*1.e4, & 551 531 & 'max= ', MAXVAL(zkz)*1.e4, ' cm2/s' 552 532 END DO 533 !!gm end bug mpp 553 534 ! 554 535 ENDIF
Note: See TracChangeset
for help on using the changeset viewer.
