- Timestamp:
- 2018-06-21T11:58:42+02:00 (6 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/UKMO/dev_r5518_nemo2cice_prints/NEMOGCM/NEMO/OPA_SRC/TRA/traadv_tvd.F90
r9816 r9817 34 34 USE timing ! Timing 35 35 USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 36 USE iom 36 37 37 38 IMPLICIT NONE … … 42 43 43 44 LOGICAL :: l_trd ! flag to compute trends 45 LOGICAL :: l_trans ! flag to output vertically integrated transports 44 46 45 47 !! * Substitutions … … 84 86 REAL(wp) :: zfp_ui, zfp_vj, zfp_wk ! - - 85 87 REAL(wp) :: zfm_ui, zfm_vj, zfm_wk ! - - 86 REAL(wp), POINTER, DIMENSION(:,:,:) :: zwi, zwz 87 REAL(wp), POINTER, DIMENSION(:,:,:) :: ztrdx, ztrdy, ztrdz 88 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zwi, zwz 89 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: ztrdx, ztrdy, ztrdz, zptry 90 REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: z2d 88 91 !!---------------------------------------------------------------------- 89 92 ! 90 93 IF( nn_timing == 1 ) CALL timing_start('tra_adv_tvd') 91 94 ! 92 CALL wrk_alloc( jpi, jpj, jpk, zwi, zwz ) 95 ALLOCATE(zwi(1:jpi, 1:jpj, 1:jpk)) 96 ALLOCATE(zwz(1:jpi, 1:jpj, 1:jpk)) 97 93 98 ! 94 99 IF( kt == kit000 ) THEN … … 97 102 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~' 98 103 ! 99 l_trd = .FALSE.100 IF( ( cdtype == 'TRA' .AND. l_trdtra ) .OR. ( cdtype == 'TRC' .AND. l_trdtrc ) ) l_trd = .TRUE.101 104 ENDIF 102 ! 103 IF( l_trd ) THEN 104 CALL wrk_alloc( jpi, jpj, jpk, ztrdx, ztrdy, ztrdz ) 105 l_trd = .FALSE. 106 l_trans = .FALSE. 107 IF( ( cdtype == 'TRA' .AND. l_trdtra ) .OR. ( cdtype == 'TRC' .AND. l_trdtrc ) ) l_trd = .TRUE. 108 IF( cdtype == 'TRA' .AND. (iom_use("uadv_heattr") .OR. iom_use("vadv_heattr") ) ) l_trans = .TRUE. 109 ! 110 IF( l_trd .OR. l_trans ) THEN 111 ALLOCATE(ztrdx(1:jpi, 1:jpj, 1:jpk)) 112 ALLOCATE(ztrdy(1:jpi, 1:jpj, 1:jpk)) 113 ALLOCATE(ztrdz(1:jpi, 1:jpj, 1:jpk)) 105 114 ztrdx(:,:,:) = 0.e0 ; ztrdy(:,:,:) = 0.e0 ; ztrdz(:,:,:) = 0.e0 115 ALLOCATE(z2d(1:jpi, 1:jpj)) 116 ENDIF 117 ! 118 IF( cdtype == 'TRA' .AND. ln_diaptr ) THEN 119 ALLOCATE(zptry(1:jpi, 1:jpj, 1:jpk)) 120 zptry(:,:,:) = 0._wp 106 121 ENDIF 107 122 ! … … 173 188 DO jj = 2, jpjm1 174 189 DO ji = fs_2, fs_jpim1 ! vector opt. 175 zbtr = 1. / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )176 190 ! total intermediate advective trends 177 ztra = - zbtr *( zwx(ji,jj,jk) - zwx(ji-1,jj ,jk ) &178 & 179 & + zwz(ji,jj,jk) - zwz(ji ,jj ,jk+1))191 ztra = - ( zwx(ji,jj,jk) - zwx(ji-1,jj ,jk ) & 192 & + zwy(ji,jj,jk) - zwy(ji ,jj-1,jk ) & 193 & + zwz(ji,jj,jk) - zwz(ji ,jj ,jk+1) ) / e1e2t(ji,jj) 180 194 ! update and guess with monotonic sheme 181 pta(ji,jj,jk,jn) = pta(ji,jj,jk,jn) + ztra* tmask(ji,jj,jk)182 zwi(ji,jj,jk) = ( ptb(ji,jj,jk,jn) + z2dtt * ztra) * tmask(ji,jj,jk)195 pta(ji,jj,jk,jn) = pta(ji,jj,jk,jn) + ztra / fse3t_n(ji,jj,jk) * tmask(ji,jj,jk) 196 zwi(ji,jj,jk) = ( fse3t_b(ji,jj,jk) * ptb(ji,jj,jk,jn) + z2dtt * ztra ) / fse3t_a(ji,jj,jk) * tmask(ji,jj,jk) 183 197 END DO 184 198 END DO … … 188 202 189 203 ! ! trend diagnostics (contribution of upstream fluxes) 190 IF( l_trd ) THEN204 IF( l_trd .OR. l_trans ) THEN 191 205 ! store intermediate advective trends 192 206 ztrdx(:,:,:) = zwx(:,:,:) ; ztrdy(:,:,:) = zwy(:,:,:) ; ztrdz(:,:,:) = zwz(:,:,:) 193 207 END IF 194 208 ! ! "Poleward" heat and salt transports (contribution of upstream fluxes) 195 IF( cdtype == 'TRA' .AND. ln_diaptr ) THEN 196 IF( jn == jp_tem ) htr_adv(:) = ptr_sj( zwy(:,:,:) ) 197 IF( jn == jp_sal ) str_adv(:) = ptr_sj( zwy(:,:,:) ) 198 ENDIF 209 IF( cdtype == 'TRA' .AND. ln_diaptr ) zptry(:,:,:) = zwy(:,:,:) 199 210 200 211 ! 3. antidiffusive flux : high order minus low order … … 254 265 255 266 ! ! trend diagnostics (contribution of upstream fluxes) 256 IF( l_trd ) THEN267 IF( l_trd .OR. l_trans ) THEN 257 268 ztrdx(:,:,:) = ztrdx(:,:,:) + zwx(:,:,:) ! <<< Add to previously computed 258 269 ztrdy(:,:,:) = ztrdy(:,:,:) + zwy(:,:,:) ! <<< Add to previously computed 259 270 ztrdz(:,:,:) = ztrdz(:,:,:) + zwz(:,:,:) ! <<< Add to previously computed 260 261 CALL trd_tra( kt, cdtype, jn, jptra_xad, ztrdx, pun, ptn(:,:,:,jn) ) 262 CALL trd_tra( kt, cdtype, jn, jptra_yad, ztrdy, pvn, ptn(:,:,:,jn) ) 263 CALL trd_tra( kt, cdtype, jn, jptra_zad, ztrdz, pwn, ptn(:,:,:,jn) ) 271 ENDIF 272 273 IF( l_trd ) THEN 274 CALL trd_tra( kt, cdtype, jn, jptra_xad, ztrdx, pun, ptn(:,:,:,jn) ) 275 CALL trd_tra( kt, cdtype, jn, jptra_yad, ztrdy, pvn, ptn(:,:,:,jn) ) 276 CALL trd_tra( kt, cdtype, jn, jptra_zad, ztrdz, pwn, ptn(:,:,:,jn) ) 264 277 END IF 265 ! ! "Poleward" heat and salt transports (contribution of upstream fluxes) 278 279 IF( l_trans .AND. jn==jp_tem ) THEN 280 z2d(:,:) = 0._wp 281 DO jk = 1, jpkm1 282 DO jj = 2, jpjm1 283 DO ji = fs_2, fs_jpim1 ! vector opt. 284 z2d(ji,jj) = z2d(ji,jj) + ztrdx(ji,jj,jk) 285 END DO 286 END DO 287 END DO 288 CALL lbc_lnk( z2d, 'U', -1. ) 289 CALL iom_put( "uadv_heattr", rau0_rcp * z2d ) ! heat transport in i-direction 290 ! 291 z2d(:,:) = 0._wp 292 DO jk = 1, jpkm1 293 DO jj = 2, jpjm1 294 DO ji = fs_2, fs_jpim1 ! vector opt. 295 z2d(ji,jj) = z2d(ji,jj) + ztrdy(ji,jj,jk) 296 END DO 297 END DO 298 END DO 299 CALL lbc_lnk( z2d, 'V', -1. ) 300 CALL iom_put( "vadv_heattr", rau0_rcp * z2d ) ! heat transport in j-direction 301 ENDIF 302 ! "Poleward" heat and salt transports (contribution of upstream fluxes) 266 303 IF( cdtype == 'TRA' .AND. ln_diaptr ) THEN 267 IF( jn == jp_tem ) htr_adv(:) = ptr_sj( zwy(:,:,:) ) + htr_adv(:)268 IF( jn == jp_sal ) str_adv(:) = ptr_sj( zwy(:,:,:) ) + str_adv(:)304 zptry(:,:,:) = zptry(:,:,:) + zwy(:,:,:) ! <<< Add to previously computed 305 CALL dia_ptr_ohst_components( jn, 'adv', zptry(:,:,:) ) 269 306 ENDIF 270 307 ! 271 308 END DO 272 309 ! 273 CALL wrk_dealloc( jpi, jpj, jpk, zwi, zwz ) 274 IF( l_trd ) CALL wrk_dealloc( jpi, jpj, jpk, ztrdx, ztrdy, ztrdz ) 310 DEALLOCATE( zwi ) 311 DEALLOCATE( zwz ) 312 IF( l_trd .OR. l_trans ) THEN 313 DEALLOCATE( ztrdx ) 314 DEALLOCATE( ztrdy ) 315 DEALLOCATE( ztrdz ) 316 DEALLOCATE( z2d ) 317 ENDIF 318 IF( cdtype == 'TRA' .AND. ln_diaptr ) DEALLOCATE( zptry ) 275 319 ! 276 320 IF( nn_timing == 1 ) CALL timing_stop('tra_adv_tvd') … … 316 360 REAL(wp) :: zfp_ui, zfp_vj, zfp_wk ! - - 317 361 REAL(wp) :: zfm_ui, zfm_vj, zfm_wk ! - - 318 REAL(wp), POINTER, DIMENSION(:,: ) :: zwx_sav , zwy_sav 319 REAL(wp), POINTER, DIMENSION(:,:,:) :: zwi, zwz, zhdiv, zwz_sav, zwzts 320 REAL(wp), POINTER, DIMENSION(:,:,:) :: ztrdx, ztrdy, ztrdz 321 REAL(wp), POINTER, DIMENSION(:,:,:,:) :: ztrs 362 REAL(wp), ALLOCATABLE, DIMENSION(:,: ) :: zwx_sav , zwy_sav 363 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zwi, zwz, zhdiv, zwz_sav, zwzts 364 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: ztrdx, ztrdy, ztrdz 365 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zptry 366 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:,:) :: ztrs 322 367 !!---------------------------------------------------------------------- 323 368 ! 324 369 IF( nn_timing == 1 ) CALL timing_start('tra_adv_tvd_zts') 325 370 ! 326 CALL wrk_alloc( jpi, jpj, zwx_sav, zwy_sav ) 327 CALL wrk_alloc( jpi, jpj, jpk, zwi, zwz , zhdiv, zwz_sav, zwzts ) 328 CALL wrk_alloc( jpi, jpj, jpk, 3, ztrs ) 371 ALLOCATE(zwx_sav(1:jpi, 1:jpj)) 372 ALLOCATE(zwy_sav(1:jpi, 1:jpj)) 373 ALLOCATE(zwi(1:jpi, 1:jpj, 1:jpk)) 374 ALLOCATE(zwz(1:jpi, 1:jpj, 1:jpk)) 375 ALLOCATE(zhdiv(1:jpi, 1:jpj, 1:jpk)) 376 ALLOCATE(zwz_sav(1:jpi, 1:jpj, 1:jpk)) 377 ALLOCATE(zwzts(1:jpi, 1:jpj, 1:jpk)) 378 ALLOCATE(ztrs(1:jpi, 1:jpj, 1:jpk, 1:kjpt+1)) 329 379 ! 330 380 IF( kt == kit000 ) THEN … … 338 388 ! 339 389 IF( l_trd ) THEN 340 CALL wrk_alloc( jpi, jpj, jpk, ztrdx, ztrdy, ztrdz ) 390 ALLOCATE(ztrdx(1:jpi, 1:jpj, 1:jpk)) 391 ALLOCATE(ztrdy(1:jpi, 1:jpj, 1:jpk)) 392 ALLOCATE(ztrdz(1:jpi, 1:jpj, 1:jpk)) 341 393 ztrdx(:,:,:) = 0._wp ; ztrdy(:,:,:) = 0._wp ; ztrdz(:,:,:) = 0._wp 394 ENDIF 395 ! 396 IF( cdtype == 'TRA' .AND. ln_diaptr ) THEN 397 ALLOCATE(zptry(1:jpi, 1:jpj, 1:jpk)) 398 zptry(:,:,:) = 0._wp 342 399 ENDIF 343 400 ! … … 410 467 DO jj = 2, jpjm1 411 468 DO ji = fs_2, fs_jpim1 ! vector opt. 412 zbtr = 1._wp / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )413 469 ! total intermediate advective trends 414 ztra = - zbtr *( zwx(ji,jj,jk) - zwx(ji-1,jj ,jk ) &415 & 416 & + zwz(ji,jj,jk) - zwz(ji ,jj ,jk+1))470 ztra = - ( zwx(ji,jj,jk) - zwx(ji-1,jj ,jk ) & 471 & + zwy(ji,jj,jk) - zwy(ji ,jj-1,jk ) & 472 & + zwz(ji,jj,jk) - zwz(ji ,jj ,jk+1) ) / e1e2t(ji,jj) 417 473 ! update and guess with monotonic sheme 418 pta(ji,jj,jk,jn) = pta(ji,jj,jk,jn) + ztra419 zwi(ji,jj,jk) = ( ptb(ji,jj,jk,jn) + z2dtt * ztra) * tmask(ji,jj,jk)474 pta(ji,jj,jk,jn) = pta(ji,jj,jk,jn) + ztra / fse3t_n(ji,jj,jk) * tmask(ji,jj,jk) 475 zwi(ji,jj,jk) = ( fse3t_b(ji,jj,jk) * ptb(ji,jj,jk,jn) + z2dtt * ztra ) / fse3t_a(ji,jj,jk) * tmask(ji,jj,jk) 420 476 END DO 421 477 END DO … … 430 486 END IF 431 487 ! ! "Poleward" heat and salt transports (contribution of upstream fluxes) 432 IF( cdtype == 'TRA' .AND. ln_diaptr ) THEN 433 IF( jn == jp_tem ) htr_adv(:) = ptr_sj( zwy(:,:,:) ) 434 IF( jn == jp_sal ) str_adv(:) = ptr_sj( zwy(:,:,:) ) 435 ENDIF 488 IF( cdtype == 'TRA' .AND. ln_diaptr ) zptry(:,:,:) = zwy(:,:,:) 436 489 437 490 ! 3. antidiffusive flux : high order minus low order 438 491 ! -------------------------------------------------- 439 492 ! antidiffusive flux on i and j 440 441 493 ! 442 494 DO jk = 1, jpkm1 443 495 ! 444 496 DO jj = 1, jpjm1 445 497 DO ji = 1, fs_jpim1 ! vector opt. … … 472 524 ! 473 525 ztrs(:,:,:,1) = ptb(:,:,:,jn) 526 ztrs(:,:,1,2) = ptb(:,:,1,jn) 527 ztrs(:,:,1,3) = ptb(:,:,1,jn) 474 528 zwzts(:,:,:) = 0._wp 475 529 … … 557 611 ! ! "Poleward" heat and salt transports (contribution of upstream fluxes) 558 612 IF( cdtype == 'TRA' .AND. ln_diaptr ) THEN 559 IF( jn == jp_tem ) htr_adv(:) = ptr_sj( zwy(:,:,:) ) + htr_adv(:)560 IF( jn == jp_sal ) str_adv(:) = ptr_sj( zwy(:,:,:) ) + str_adv(:)613 zptry(:,:,:) = zptry(:,:,:) + zwy(:,:,:) 614 CALL dia_ptr_ohst_components( jn, 'adv', zptry(:,:,:) ) 561 615 ENDIF 562 616 ! 563 617 END DO 564 618 ! 565 CALL wrk_dealloc( jpi, jpj, jpk, zwi, zwz, zhdiv, zwz_sav, zwzts ) 566 CALL wrk_dealloc( jpi, jpj, jpk, 3, ztrs ) 567 CALL wrk_dealloc( jpi, jpj, zwx_sav, zwy_sav ) 568 IF( l_trd ) CALL wrk_dealloc( jpi, jpj, jpk, ztrdx, ztrdy, ztrdz ) 619 DEALLOCATE(zwi) 620 DEALLOCATE(zwz) 621 DEALLOCATE(zhdiv) 622 DEALLOCATE(zwz_sav) 623 DEALLOCATE(zwzts) 624 DEALLOCATE(ztrs ) 625 DEALLOCATE(zwx_sav) 626 DEALLOCATE(zwy_sav ) 627 628 IF( l_trd ) THEN 629 DEALLOCATE(ztrdx) 630 DEALLOCATE(ztrdy) 631 DEALLOCATE(ztrdz) 632 END IF 633 634 IF( cdtype == 'TRA' .AND. ln_diaptr ) DEALLOCATE(zptry ) 569 635 ! 570 636 IF( nn_timing == 1 ) CALL timing_stop('tra_adv_tvd_zts') 571 637 ! 572 638 END SUBROUTINE tra_adv_tvd_zts 639 573 640 574 641 SUBROUTINE nonosc( pbef, paa, pbb, pcc, paft, p2dt ) … … 593 660 REAL(wp) :: zpos, zneg, zbt, za, zb, zc, zbig, zrtrn, z2dtt ! local scalars 594 661 REAL(wp) :: zau, zbu, zcu, zav, zbv, zcv, zup, zdo ! - - 595 REAL(wp), POINTER, DIMENSION(:,:,:) :: zbetup, zbetdo, zbup, zbdo662 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zbetup, zbetdo, zbup, zbdo 596 663 !!---------------------------------------------------------------------- 597 664 ! 598 665 IF( nn_timing == 1 ) CALL timing_start('nonosc') 599 666 ! 600 CALL wrk_alloc( jpi, jpj, jpk, zbetup, zbetdo, zbup, zbdo ) 667 ALLOCATE(zbetup(1:jpi, 1:jpj, 1:jpk)) 668 ALLOCATE(zbetdo(1:jpi, 1:jpj, 1:jpk)) 669 ALLOCATE(zbup(1:jpi, 1:jpj, 1:jpk)) 670 ALLOCATE(zbdo(1:jpi, 1:jpj, 1:jpk)) 601 671 ! 602 672 zbig = 1.e+40_wp … … 675 745 CALL lbc_lnk( paa, 'U', -1. ) ; CALL lbc_lnk( pbb, 'V', -1. ) ! lateral boundary condition (changed sign) 676 746 ! 677 CALL wrk_dealloc( jpi, jpj, jpk, zbetup, zbetdo, zbup, zbdo ) 747 DEALLOCATE(zbetup) 748 DEALLOCATE(zbetdo) 749 DEALLOCATE(zbup) 750 DEALLOCATE(zbdo) 678 751 ! 679 752 IF( nn_timing == 1 ) CALL timing_stop('nonosc')
Note: See TracChangeset
for help on using the changeset viewer.
