- Timestamp:
- 2016-11-16T12:48:46+01:00 (7 years ago)
- Location:
- branches/2016/dev_r7233_CMIP6_diags_trunk_version/NEMOGCM/NEMO/OPA_SRC
- Files:
-
- 9 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2016/dev_r7233_CMIP6_diags_trunk_version/NEMOGCM/NEMO/OPA_SRC/DIA/diaptr.F90
r6140 r7236 9 9 !! 3.3 ! 2010-10 (G. Madec) dynamical allocation 10 10 !! 3.6 ! 2014-12 (C. Ethe) use of IOM 11 !! 3.6 ! 2016-06 (T. Graham) Addition of diagnostics for CMIP6 11 12 !!---------------------------------------------------------------------- 12 13 … … 21 22 USE dom_oce ! ocean space and time domain 22 23 USE phycst ! physical constants 24 USE ldftra 23 25 ! 24 26 USE iom ! IOM library … … 38 40 PUBLIC dia_ptr_init ! call in step module 39 41 PUBLIC dia_ptr ! call in step module 42 PUBLIC dia_ptr_ohst_components ! called from tra_ldf/tra_adv routines 40 43 41 44 ! !!** namelist namptr ** 42 REAL(wp), ALLOCATABLE, SAVE, PUBLIC, DIMENSION(:) :: htr_adv, htr_ldf !: Heat TRansports (adv, diff, overturn.) 43 REAL(wp), ALLOCATABLE, SAVE, PUBLIC, DIMENSION(:) :: str_adv, str_ldf !: Salt TRansports (adv, diff, overturn.) 44 45 REAL(wp), ALLOCATABLE, SAVE, PUBLIC, DIMENSION(:,:) :: htr_adv, htr_ldf, htr_eiv, htr_vt !: Heat TRansports (adv, diff, Bolus.) 46 REAL(wp), ALLOCATABLE, SAVE, PUBLIC, DIMENSION(:,:) :: str_adv, str_ldf, str_eiv, str_vs !: Salt TRansports (adv, diff, Bolus.) 47 REAL(wp), ALLOCATABLE, SAVE, PUBLIC, DIMENSION(:,:) :: htr_ove, str_ove !: heat Salt TRansports ( overturn.) 48 REAL(wp), ALLOCATABLE, SAVE, PUBLIC, DIMENSION(:,:) :: htr_btr, str_btr !: heat Salt TRansports ( barotropic ) 45 49 46 50 LOGICAL, PUBLIC :: ln_diaptr ! Poleward transport flag (T) or not (F) 47 51 LOGICAL, PUBLIC :: ln_subbas ! Atlantic/Pacific/Indian basins calculation 48 INTEGER 52 INTEGER, PUBLIC :: nptr ! = 1 (l_subbas=F) or = 5 (glo, atl, pac, ind, ipc) (l_subbas=T) 49 53 50 54 REAL(wp) :: rc_sv = 1.e-6_wp ! conversion from m3/s to Sverdrup … … 75 79 ! 76 80 INTEGER :: ji, jj, jk, jn ! dummy loop indices 77 REAL(wp) :: z v, zsfc ! local scalar81 REAL(wp) :: zsfc,zvfc ! local scalar 78 82 REAL(wp), DIMENSION(jpi,jpj) :: z2d ! 2D workspace 79 83 REAL(wp), DIMENSION(jpi,jpj,jpk) :: z3d ! 3D workspace 80 84 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zmask ! 3D workspace 81 85 REAL(wp), DIMENSION(jpi,jpj,jpk,jpts) :: zts ! 3D workspace 82 CHARACTER( len = 10 ) :: cl1 86 REAL(wp), DIMENSION(jpj) :: vsum ! 1D workspace 87 REAL(wp), DIMENSION(jpj,jpts) :: tssum ! 1D workspace 88 89 ! 90 !overturning calculation 91 REAL(wp), DIMENSION(jpj,jpk,nptr) :: sjk , r1_sjk ! i-mean i-k-surface and its inverse 92 REAL(wp), DIMENSION(jpj,jpk,nptr) :: v_msf, sn_jk , tn_jk ! i-mean T and S, j-Stream-Function 93 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zvn ! 3D workspace 94 95 96 CHARACTER( len = 12 ) :: cl1 83 97 !!---------------------------------------------------------------------- 84 98 ! … … 109 123 END DO 110 124 ENDIF 125 IF( iom_use("sopstove") .OR. iom_use("sophtove") .OR. iom_use("sopstbtr") .OR. iom_use("sophtbtr") ) THEN 126 ! define fields multiplied by scalar 127 zmask(:,:,:) = 0._wp 128 zts(:,:,:,:) = 0._wp 129 zvn(:,:,:) = 0._wp 130 DO jk = 1, jpkm1 131 DO jj = 1, jpjm1 132 DO ji = 1, jpi 133 zvfc = e1v(ji,jj) * e3v_n(ji,jj,jk) 134 zmask(ji,jj,jk) = vmask(ji,jj,jk) * zvfc 135 zts(ji,jj,jk,jp_tem) = (tsn(ji,jj,jk,jp_tem)+tsn(ji,jj+1,jk,jp_tem)) * 0.5 * zvfc !Tracers averaged onto V grid 136 zts(ji,jj,jk,jp_sal) = (tsn(ji,jj,jk,jp_sal)+tsn(ji,jj+1,jk,jp_sal)) * 0.5 * zvfc 137 zvn(ji,jj,jk) = vn(ji,jj,jk) * zvfc 138 ENDDO 139 ENDDO 140 ENDDO 141 ENDIF 142 IF( iom_use("sopstove") .OR. iom_use("sophtove") ) THEN 143 sjk(:,:,1) = ptr_sjk( zmask(:,:,:), btmsk(:,:,1) ) 144 r1_sjk(:,:,1) = 0._wp 145 WHERE( sjk(:,:,1) /= 0._wp ) r1_sjk(:,:,1) = 1._wp / sjk(:,:,1) 146 147 ! i-mean T and S, j-Stream-Function, global 148 tn_jk(:,:,1) = ptr_sjk( zts(:,:,:,jp_tem) ) * r1_sjk(:,:,1) 149 sn_jk(:,:,1) = ptr_sjk( zts(:,:,:,jp_sal) ) * r1_sjk(:,:,1) 150 v_msf(:,:,1) = ptr_sjk( zvn(:,:,:) ) 151 152 htr_ove(:,1) = SUM( v_msf(:,:,1)*tn_jk(:,:,1) ,2 ) 153 str_ove(:,1) = SUM( v_msf(:,:,1)*sn_jk(:,:,1) ,2 ) 154 155 z2d(1,:) = htr_ove(:,1) * rc_pwatt ! (conversion in PW) 156 DO ji = 1, jpi 157 z2d(ji,:) = z2d(1,:) 158 ENDDO 159 cl1 = 'sophtove' 160 CALL iom_put( TRIM(cl1), z2d ) 161 z2d(1,:) = str_ove(:,1) * rc_ggram ! (conversion in Gg) 162 DO ji = 1, jpi 163 z2d(ji,:) = z2d(1,:) 164 ENDDO 165 cl1 = 'sopstove' 166 CALL iom_put( TRIM(cl1), z2d ) 167 IF( ln_subbas ) THEN 168 DO jn = 2, nptr 169 sjk(:,:,jn) = ptr_sjk( zmask(:,:,:), btmsk(:,:,jn) ) 170 r1_sjk(:,:,jn) = 0._wp 171 WHERE( sjk(:,:,jn) /= 0._wp ) r1_sjk(:,:,jn) = 1._wp / sjk(:,:,jn) 172 173 ! i-mean T and S, j-Stream-Function, basin 174 tn_jk(:,:,jn) = ptr_sjk( zts(:,:,:,jp_tem), btmsk(:,:,jn) ) * r1_sjk(:,:,jn) 175 sn_jk(:,:,jn) = ptr_sjk( zts(:,:,:,jp_sal), btmsk(:,:,jn) ) * r1_sjk(:,:,jn) 176 v_msf(:,:,jn) = ptr_sjk( zvn(:,:,:), btmsk(:,:,jn) ) 177 htr_ove(:,jn) = SUM( v_msf(:,:,jn)*tn_jk(:,:,jn) ,2 ) 178 str_ove(:,jn) = SUM( v_msf(:,:,jn)*sn_jk(:,:,jn) ,2 ) 179 180 z2d(1,:) = htr_ove(:,jn) * rc_pwatt ! (conversion in PW) 181 DO ji = 1, jpi 182 z2d(ji,:) = z2d(1,:) 183 ENDDO 184 cl1 = TRIM('sophtove_'//clsubb(jn)) 185 CALL iom_put( cl1, z2d ) 186 z2d(1,:) = str_ove(:,jn) * rc_ggram ! (conversion in Gg) 187 DO ji = 1, jpi 188 z2d(ji,:) = z2d(1,:) 189 ENDDO 190 cl1 = TRIM('sopstove_'//clsubb(jn)) 191 CALL iom_put( cl1, z2d ) 192 END DO 193 ENDIF 194 ENDIF 195 IF( iom_use("sopstbtr") .OR. iom_use("sophtbtr") ) THEN 196 ! Calculate barotropic heat and salt transport here 197 sjk(:,1,1) = ptr_sj( zmask(:,:,:), btmsk(:,:,1) ) 198 r1_sjk(:,1,1) = 0._wp 199 WHERE( sjk(:,1,1) /= 0._wp ) r1_sjk(:,1,1) = 1._wp / sjk(:,1,1) 200 201 vsum = ptr_sj( zvn(:,:,:), btmsk(:,:,1)) 202 tssum(:,jp_tem) = ptr_sj( zts(:,:,:,jp_tem), btmsk(:,:,1) ) 203 tssum(:,jp_sal) = ptr_sj( zts(:,:,:,jp_sal), btmsk(:,:,1) ) 204 htr_btr(:,1) = vsum * tssum(:,jp_tem) * r1_sjk(:,1,1) 205 str_btr(:,1) = vsum * tssum(:,jp_sal) * r1_sjk(:,1,1) 206 z2d(1,:) = htr_btr(:,1) * rc_pwatt ! (conversion in PW) 207 DO ji = 2, jpi 208 z2d(ji,:) = z2d(1,:) 209 ENDDO 210 cl1 = 'sophtbtr' 211 CALL iom_put( TRIM(cl1), z2d ) 212 z2d(1,:) = str_btr(:,1) * rc_ggram ! (conversion in Gg) 213 DO ji = 2, jpi 214 z2d(ji,:) = z2d(1,:) 215 ENDDO 216 cl1 = 'sopstbtr' 217 CALL iom_put( TRIM(cl1), z2d ) 218 IF( ln_subbas ) THEN 219 DO jn = 2, nptr 220 sjk(:,1,jn) = ptr_sj( zmask(:,:,:), btmsk(:,:,jn) ) 221 r1_sjk(:,1,jn) = 0._wp 222 WHERE( sjk(:,1,jn) /= 0._wp ) r1_sjk(:,1,jn) = 1._wp / sjk(:,1,jn) 223 vsum = ptr_sj( zvn(:,:,:), btmsk(:,:,jn)) 224 tssum(:,jp_tem) = ptr_sj( zts(:,:,:,jp_tem), btmsk(:,:,jn) ) 225 tssum(:,jp_sal) = ptr_sj( zts(:,:,:,jp_sal), btmsk(:,:,jn) ) 226 htr_btr(:,jn) = vsum * tssum(:,jp_tem) * r1_sjk(:,1,jn) 227 str_btr(:,jn) = vsum * tssum(:,jp_sal) * r1_sjk(:,1,jn) 228 z2d(1,:) = htr_btr(:,jn) * rc_pwatt ! (conversion in PW) 229 DO ji = 1, jpi 230 z2d(ji,:) = z2d(1,:) 231 ENDDO 232 cl1 = TRIM('sophtbtr_'//clsubb(jn)) 233 CALL iom_put( cl1, z2d ) 234 z2d(1,:) = str_btr(:,jn) * rc_ggram ! (conversion in Gg) 235 DO ji = 1, jpi 236 z2d(ji,:) = z2d(1,:) 237 ENDDO 238 cl1 = TRIM('sopstbtr_'//clsubb(jn)) 239 CALL iom_put( cl1, z2d ) 240 ENDDO 241 ENDIF !ln_subbas 242 ENDIF !iom_use("sopstbtr....) 111 243 ! 112 244 ELSE … … 148 280 ! ! Advective and diffusive heat and salt transport 149 281 IF( iom_use("sophtadv") .OR. iom_use("sopstadv") ) THEN 150 z2d(1,:) = htr_adv(: ) * rc_pwatt ! (conversion in PW)282 z2d(1,:) = htr_adv(:,1) * rc_pwatt ! (conversion in PW) 151 283 DO ji = 1, jpi 152 284 z2d(ji,:) = z2d(1,:) … … 154 286 cl1 = 'sophtadv' 155 287 CALL iom_put( TRIM(cl1), z2d ) 156 z2d(1,:) = str_adv(: ) * rc_ggram ! (conversion in Gg)288 z2d(1,:) = str_adv(:,1) * rc_ggram ! (conversion in Gg) 157 289 DO ji = 1, jpi 158 290 z2d(ji,:) = z2d(1,:) … … 160 292 cl1 = 'sopstadv' 161 293 CALL iom_put( TRIM(cl1), z2d ) 294 IF( ln_subbas ) THEN 295 DO jn=2,nptr 296 z2d(1,:) = htr_adv(:,jn) * rc_pwatt ! (conversion in PW) 297 DO ji = 1, jpi 298 z2d(ji,:) = z2d(1,:) 299 ENDDO 300 cl1 = TRIM('sophtadv_'//clsubb(jn)) 301 CALL iom_put( cl1, z2d ) 302 z2d(1,:) = str_adv(:,jn) * rc_ggram ! (conversion in Gg) 303 DO ji = 1, jpi 304 z2d(ji,:) = z2d(1,:) 305 ENDDO 306 cl1 = TRIM('sopstadv_'//clsubb(jn)) 307 CALL iom_put( cl1, z2d ) 308 ENDDO 309 ENDIF 162 310 ENDIF 163 311 ! 164 312 IF( iom_use("sophtldf") .OR. iom_use("sopstldf") ) THEN 165 z2d(1,:) = htr_ldf(: ) * rc_pwatt ! (conversion in PW)313 z2d(1,:) = htr_ldf(:,1) * rc_pwatt ! (conversion in PW) 166 314 DO ji = 1, jpi 167 315 z2d(ji,:) = z2d(1,:) … … 169 317 cl1 = 'sophtldf' 170 318 CALL iom_put( TRIM(cl1), z2d ) 171 z2d(1,:) = str_ldf(: ) * rc_ggram ! (conversion in Gg)319 z2d(1,:) = str_ldf(:,1) * rc_ggram ! (conversion in Gg) 172 320 DO ji = 1, jpi 173 321 z2d(ji,:) = z2d(1,:) … … 175 323 cl1 = 'sopstldf' 176 324 CALL iom_put( TRIM(cl1), z2d ) 325 IF( ln_subbas ) THEN 326 DO jn=2,nptr 327 z2d(1,:) = htr_ldf(:,jn) * rc_pwatt ! (conversion in PW) 328 DO ji = 1, jpi 329 z2d(ji,:) = z2d(1,:) 330 ENDDO 331 cl1 = TRIM('sophtldf_'//clsubb(jn)) 332 CALL iom_put( cl1, z2d ) 333 z2d(1,:) = str_ldf(:,jn) * rc_ggram ! (conversion in Gg) 334 DO ji = 1, jpi 335 z2d(ji,:) = z2d(1,:) 336 ENDDO 337 cl1 = TRIM('sopstldf_'//clsubb(jn)) 338 CALL iom_put( cl1, z2d ) 339 ENDDO 340 ENDIF 341 ENDIF 342 343 IF( iom_use("sopht_vt") .OR. iom_use("sopst_vs") ) THEN 344 z2d(1,:) = htr_vt(:,1) * rc_pwatt ! (conversion in PW) 345 DO ji = 1, jpi 346 z2d(ji,:) = z2d(1,:) 347 ENDDO 348 cl1 = 'sopht_vt' 349 CALL iom_put( TRIM(cl1), z2d ) 350 z2d(1,:) = str_vs(:,1) * rc_ggram ! (conversion in Gg) 351 DO ji = 1, jpi 352 z2d(ji,:) = z2d(1,:) 353 ENDDO 354 cl1 = 'sopst_vs' 355 CALL iom_put( TRIM(cl1), z2d ) 356 IF( ln_subbas ) THEN 357 DO jn=2,nptr 358 z2d(1,:) = htr_vt(:,jn) * rc_pwatt ! (conversion in PW) 359 DO ji = 1, jpi 360 z2d(ji,:) = z2d(1,:) 361 ENDDO 362 cl1 = TRIM('sopht_vt_'//clsubb(jn)) 363 CALL iom_put( cl1, z2d ) 364 z2d(1,:) = str_vs(:,jn) * rc_ggram ! (conversion in Gg) 365 DO ji = 1, jpi 366 z2d(ji,:) = z2d(1,:) 367 ENDDO 368 cl1 = TRIM('sopst_vs_'//clsubb(jn)) 369 CALL iom_put( cl1, z2d ) 370 ENDDO 371 ENDIF 372 ENDIF 373 374 IF(ln_ldfeiv) THEN 375 IF( iom_use("sophteiv") .OR. iom_use("sopsteiv") ) THEN 376 z2d(1,:) = htr_eiv(:,1) * rc_pwatt ! (conversion in PW) 377 DO ji = 1, jpi 378 z2d(ji,:) = z2d(1,:) 379 ENDDO 380 cl1 = 'sophteiv' 381 CALL iom_put( TRIM(cl1), z2d ) 382 z2d(1,:) = str_eiv(:,1) * rc_ggram ! (conversion in Gg) 383 DO ji = 1, jpi 384 z2d(ji,:) = z2d(1,:) 385 ENDDO 386 cl1 = 'sopsteiv' 387 CALL iom_put( TRIM(cl1), z2d ) 388 IF( ln_subbas ) THEN 389 DO jn=2,nptr 390 z2d(1,:) = htr_eiv(:,jn) * rc_pwatt ! (conversion in PW) 391 DO ji = 1, jpi 392 z2d(ji,:) = z2d(1,:) 393 ENDDO 394 cl1 = TRIM('sophteiv_'//clsubb(jn)) 395 CALL iom_put( cl1, z2d ) 396 z2d(1,:) = str_eiv(:,jn) * rc_ggram ! (conversion in Gg) 397 DO ji = 1, jpi 398 z2d(ji,:) = z2d(1,:) 399 ENDDO 400 cl1 = TRIM('sopsteiv_'//clsubb(jn)) 401 CALL iom_put( cl1, z2d ) 402 ENDDO 403 ENDIF 404 ENDIF 177 405 ENDIF 178 406 ! … … 254 482 ! Initialise arrays to zero because diatpr is called before they are first calculated 255 483 ! Note that this means diagnostics will not be exactly correct when model run is restarted. 256 htr_adv(:) = 0._wp ; str_adv(:) = 0._wp 257 htr_ldf(:) = 0._wp ; str_ldf(:) = 0._wp 484 htr_adv(:,:) = 0._wp ; str_adv(:,:) = 0._wp 485 htr_ldf(:,:) = 0._wp ; str_ldf(:,:) = 0._wp 486 htr_eiv(:,:) = 0._wp ; str_eiv(:,:) = 0._wp 487 htr_vt(:,:) = 0._wp ; str_vs(:,:) = 0._wp 488 htr_ove(:,:) = 0._wp ; str_ove(:,:) = 0._wp 489 htr_btr(:,:) = 0._wp ; str_btr(:,:) = 0._wp 258 490 ! 259 491 ENDIF … … 261 493 END SUBROUTINE dia_ptr_init 262 494 495 SUBROUTINE dia_ptr_ohst_components( ktra, cptr, pva ) 496 !!---------------------------------------------------------------------- 497 !! *** ROUTINE dia_ptr_ohst_components *** 498 !!---------------------------------------------------------------------- 499 !! Wrapper for heat and salt transport calculations to calculate them for each basin 500 !! Called from all advection and/or diffusion routines 501 !!---------------------------------------------------------------------- 502 INTEGER , INTENT(in ) :: ktra ! tracer index 503 CHARACTER(len=3) , INTENT(in) :: cptr ! transport type 'adv'/'ldf'/'eiv' 504 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in) :: pva ! 3D input array of advection/diffusion 505 INTEGER :: jn ! 506 507 IF( cptr == 'adv' ) THEN 508 IF( ktra == jp_tem ) htr_adv(:,1) = ptr_sj( pva(:,:,:) ) 509 IF( ktra == jp_sal ) str_adv(:,1) = ptr_sj( pva(:,:,:) ) 510 ENDIF 511 IF( cptr == 'ldf' ) THEN 512 IF( ktra == jp_tem ) htr_ldf(:,1) = ptr_sj( pva(:,:,:) ) 513 IF( ktra == jp_sal ) str_ldf(:,1) = ptr_sj( pva(:,:,:) ) 514 ENDIF 515 IF( cptr == 'eiv' ) THEN 516 IF( ktra == jp_tem ) htr_eiv(:,1) = ptr_sj( pva(:,:,:) ) 517 IF( ktra == jp_sal ) str_eiv(:,1) = ptr_sj( pva(:,:,:) ) 518 ENDIF 519 IF( cptr == 'vts' ) THEN 520 IF( ktra == jp_tem ) htr_vt(:,1) = ptr_sj( pva(:,:,:) ) 521 IF( ktra == jp_sal ) str_vs(:,1) = ptr_sj( pva(:,:,:) ) 522 ENDIF 523 ! 524 IF( ln_subbas ) THEN 525 ! 526 IF( cptr == 'adv' ) THEN 527 IF( ktra == jp_tem ) THEN 528 DO jn = 2, nptr 529 htr_adv(:,jn) = ptr_sj( pva(:,:,:), btmsk(:,:,jn) ) 530 END DO 531 ENDIF 532 IF( ktra == jp_sal ) THEN 533 DO jn = 2, nptr 534 str_adv(:,jn) = ptr_sj( pva(:,:,:), btmsk(:,:,jn) ) 535 END DO 536 ENDIF 537 ENDIF 538 IF( cptr == 'ldf' ) THEN 539 IF( ktra == jp_tem ) THEN 540 DO jn = 2, nptr 541 htr_ldf(:,jn) = ptr_sj( pva(:,:,:), btmsk(:,:,jn) ) 542 END DO 543 ENDIF 544 IF( ktra == jp_sal ) THEN 545 DO jn = 2, nptr 546 str_ldf(:,jn) = ptr_sj( pva(:,:,:), btmsk(:,:,jn) ) 547 END DO 548 ENDIF 549 ENDIF 550 IF( cptr == 'eiv' ) THEN 551 IF( ktra == jp_tem ) THEN 552 DO jn = 2, nptr 553 htr_eiv(:,jn) = ptr_sj( pva(:,:,:), btmsk(:,:,jn) ) 554 END DO 555 ENDIF 556 IF( ktra == jp_sal ) THEN 557 DO jn = 2, nptr 558 str_eiv(:,jn) = ptr_sj( pva(:,:,:), btmsk(:,:,jn) ) 559 END DO 560 ENDIF 561 ENDIF 562 IF( cptr == 'vts' ) THEN 563 IF( ktra == jp_tem ) THEN 564 DO jn = 2, nptr 565 htr_vt(:,jn) = ptr_sj( pva(:,:,:), btmsk(:,:,jn) ) 566 END DO 567 ENDIF 568 IF( ktra == jp_sal ) THEN 569 DO jn = 2, nptr 570 str_vs(:,jn) = ptr_sj( pva(:,:,:), btmsk(:,:,jn) ) 571 END DO 572 ENDIF 573 ENDIF 574 ! 575 ENDIF 576 END SUBROUTINE dia_ptr_ohst_components 577 263 578 264 579 FUNCTION dia_ptr_alloc() … … 271 586 ierr(:) = 0 272 587 ! 273 ALLOCATE( btmsk(jpi,jpj,nptr) , & 274 & htr_adv(jpj) , str_adv(jpj) , & 275 & htr_ldf(jpj) , str_ldf(jpj) , STAT=ierr(1) ) 588 ALLOCATE( btmsk(jpi,jpj,nptr) , & 589 & htr_adv(jpj,nptr) , str_adv(jpj,nptr) , & 590 & htr_eiv(jpj,nptr) , str_eiv(jpj,nptr) , & 591 & htr_vt(jpj,nptr) , str_vs(jpj,nptr) , & 592 & htr_ove(jpj,nptr) , str_ove(jpj,nptr) , & 593 & htr_btr(jpj,nptr) , str_btr(jpj,nptr) , & 594 & htr_ldf(jpj,nptr) , str_ldf(jpj,nptr) , STAT=ierr(1) ) 276 595 ! 277 596 ALLOCATE( p_fval1d(jpj), p_fval2d(jpj,jpk), Stat=ierr(2)) -
branches/2016/dev_r7233_CMIP6_diags_trunk_version/NEMOGCM/NEMO/OPA_SRC/TRA/traadv_cen.F90
r6140 r7236 189 189 CALL trd_tra( kt, cdtype, jn, jptra_zad, zwz, pwn, ptn(:,:,:,jn) ) 190 190 END IF 191 ! ! "Poleward" heat and salt transports (contribution of upstream fluxes) 192 IF( cdtype == 'TRA' .AND. ln_diaptr ) THEN 193 IF( jn == jp_tem ) htr_adv(:) = ptr_sj( zwy(:,:,:) ) 194 IF( jn == jp_sal ) str_adv(:) = ptr_sj( zwy(:,:,:) ) 195 ENDIF 191 ! ! "Poleward" heat and salt transports (contribution of upstream fluxes) 192 IF( cdtype == 'TRA' .AND. ln_diaptr ) CALL dia_ptr_ohst_components( jn, 'adv', zwy(:,:,:) ) 196 193 ! 197 194 END DO -
branches/2016/dev_r7233_CMIP6_diags_trunk_version/NEMOGCM/NEMO/OPA_SRC/TRA/traadv_fct.F90
r6771 r7236 80 80 REAL(wp) :: zfm_ui, zfm_vj, zfm_wk, zC2t_v, zC4t_v ! - - 81 81 REAL(wp), POINTER, DIMENSION(:,:,:) :: zwi, zwx, zwy, zwz, ztu, ztv, zltu, zltv, ztw 82 REAL(wp), POINTER, DIMENSION(:,:,:) :: ztrdx, ztrdy, ztrd z82 REAL(wp), POINTER, DIMENSION(:,:,:) :: ztrdx, ztrdy, ztrd, zptry 83 83 !!---------------------------------------------------------------------- 84 84 ! … … 101 101 ENDIF 102 102 ! 103 IF( cdtype == 'TRA' .AND. ln_diaptr ) THEN 104 CALL wrk_alloc( jpi, jpj, jpk, zptry ) 105 zptry(:,:,:) = 0._wp 106 ENDIF 103 107 ! ! surface & bottom value : flux set to zero one for all 104 108 zwz(:,:, 1 ) = 0._wp … … 165 169 END IF 166 170 ! ! "Poleward" heat and salt transports (contribution of upstream fluxes) 167 IF( cdtype == 'TRA' .AND. ln_diaptr ) THEN 168 IF( jn == jp_tem ) htr_adv(:) = ptr_sj( zwy(:,:,:) ) 169 IF( jn == jp_sal ) str_adv(:) = ptr_sj( zwy(:,:,:) ) 170 ENDIF 171 IF( cdtype == 'TRA' .AND. ln_diaptr ) zptry(:,:,:) = zwy(:,:,:) 171 172 ! 172 173 ! !== anti-diffusive flux : high order minus low order ==! … … 305 306 ! ! "Poleward" heat and salt transports (contribution of upstream fluxes) 306 307 IF( cdtype == 'TRA' .AND. ln_diaptr ) THEN 307 IF( jn == jp_tem ) htr_adv(:) = htr_adv(:) + ptr_sj( zwy(:,:,:) )308 IF( jn == jp_sal ) str_adv(:) = str_adv(:) + ptr_sj( zwy(:,:,:) )308 zptry(:,:,:) = zptry(:,:,:) + zwy(:,:,:) ! <<< Add to previously computed 309 CALL dia_ptr_ohst_components( jn, 'adv', zptry(:,:,:) ) 309 310 ENDIF 310 311 ! … … 312 313 ! 313 314 CALL wrk_dealloc( jpi,jpj,jpk, zwi, zwx, zwy, zwz, ztu, ztv, zltu, zltv, ztw ) 315 IF( cdtype == 'TRA' .AND. ln_diaptr ) CALL wrk_dealloc( jpi, jpj, jpk, zptry ) 314 316 ! 315 317 IF( nn_timing == 1 ) CALL timing_stop('tra_adv_fct') … … 357 359 REAL(wp), POINTER, DIMENSION(:,:,:) :: zwi, zwx, zwy, zwz, zhdiv, zwzts, zwz_sav 358 360 REAL(wp), POINTER, DIMENSION(:,:,:) :: ztrdx, ztrdy, ztrdz 361 REAL(wp), POINTER, DIMENSION(:,:,:) :: zptry 359 362 REAL(wp), POINTER, DIMENSION(:,:,:,:) :: ztrs 360 363 !!---------------------------------------------------------------------- … … 380 383 ENDIF 381 384 ! 385 IF( cdtype == 'TRA' .AND. ln_diaptr ) THEN 386 CALL wrk_alloc( jpi, jpj,jpk, zptry ) 387 zptry(:,:,:) = 0._wp 388 ENDIF 382 389 zwi(:,:,:) = 0._wp 383 390 z_rzts = 1._wp / REAL( kn_fct_zts, wp ) … … 449 456 END IF 450 457 ! ! "Poleward" heat and salt transports (contribution of upstream fluxes) 451 IF( cdtype == 'TRA' .AND. ln_diaptr ) THEN 452 IF( jn == jp_tem ) htr_adv(:) = ptr_sj( zwy(:,:,:) ) 453 IF( jn == jp_sal ) str_adv(:) = ptr_sj( zwy(:,:,:) ) 454 ENDIF 458 IF( cdtype == 'TRA' .AND. ln_diaptr ) zptry(:,:,:) = zwy(:,:,:) 455 459 456 460 ! 3. anti-diffusive flux : high order minus low order … … 582 586 ! ! "Poleward" heat and salt transports (contribution of upstream fluxes) 583 587 IF( cdtype == 'TRA' .AND. ln_diaptr ) THEN 584 IF( jn == jp_tem ) htr_adv(:) = ptr_sj( zwy(:,:,:) ) + htr_adv(:)585 IF( jn == jp_sal ) str_adv(:) = ptr_sj( zwy(:,:,:) ) + str_adv(:)588 zptry(:,:,:) = zptry(:,:,:) + zwy(:,:,:) 589 CALL dia_ptr_ohst_components( jn, 'adv', zptry(:,:,:) ) 586 590 ENDIF 587 591 ! … … 591 595 CALL wrk_alloc( jpi,jpj, jpk, zwx, zwy, zwz, zwi, zhdiv, zwzts, zwz_sav ) 592 596 CALL wrk_alloc( jpi,jpj,jpk,kjpt+1, ztrs ) 597 IF( cdtype == 'TRA' .AND. ln_diaptr ) CALL wrk_dealloc( jpi, jpj, jpk, zptry ) 593 598 ! 594 599 IF( nn_timing == 1 ) CALL timing_stop('tra_adv_fct_zts') -
branches/2016/dev_r7233_CMIP6_diags_trunk_version/NEMOGCM/NEMO/OPA_SRC/TRA/traadv_mus.F90
r6140 r7236 197 197 CALL trd_tra( kt, cdtype, jn, jptra_yad, zwy, pvn, ptb(:,:,:,jn) ) 198 198 END IF 199 ! ! "Poleward" heat and salt transports 200 IF( cdtype == 'TRA' .AND. ln_diaptr ) CALL dia_ptr_ohst_components( jn, 'adv', zwy(:,:,:) ) 199 201 ! ! "Poleward" heat and salt transports 200 IF( cdtype == 'TRA' .AND. ln_diaptr ) THEN201 IF( jn == jp_tem ) htr_adv(:) = ptr_sj( zwy(:,:,:) )202 IF( jn == jp_sal ) str_adv(:) = ptr_sj( zwy(:,:,:) )203 ENDIF204 202 ! 205 203 ! !* Vertical advective fluxes -
branches/2016/dev_r7233_CMIP6_diags_trunk_version/NEMOGCM/NEMO/OPA_SRC/TRA/traadv_qck.F90
r6140 r7236 349 349 IF( l_trd ) CALL trd_tra( kt, cdtype, jn, jptra_yad, zwy, pvn, ptn(:,:,:,jn) ) 350 350 ! ! "Poleward" heat and salt transports (contribution of upstream fluxes) 351 IF( cdtype == 'TRA' .AND. ln_diaptr ) THEN 352 IF( jn == jp_tem ) htr_adv(:) = ptr_sj( zwy(:,:,:) ) 353 IF( jn == jp_sal ) str_adv(:) = ptr_sj( zwy(:,:,:) ) 354 ENDIF 351 IF( cdtype == 'TRA' .AND. ln_diaptr ) CALL dia_ptr_ohst_components( jn, 'adv', zwy(:,:,:) ) 355 352 ! 356 353 END DO -
branches/2016/dev_r7233_CMIP6_diags_trunk_version/NEMOGCM/NEMO/OPA_SRC/TRA/traadv_ubs.F90
r6140 r7236 177 177 END IF 178 178 ! ! "Poleward" heat and salt transports (contribution of upstream fluxes) 179 IF( cdtype == 'TRA' .AND. ln_diaptr ) THEN 180 IF( jn == jp_tem ) htr_adv(:) = ptr_sj( ztv(:,:,:) ) 181 IF( jn == jp_sal ) str_adv(:) = ptr_sj( ztv(:,:,:) ) 182 ENDIF 179 IF( cdtype == 'TRA' .AND. ln_diaptr ) CALL dia_ptr_ohst_components( jn, 'adv', ztv(:,:,:) ) 183 180 ! 184 181 ! !== vertical advective trend ==! -
branches/2016/dev_r7233_CMIP6_diags_trunk_version/NEMOGCM/NEMO/OPA_SRC/TRA/traldf_iso.F90
r6140 r7236 369 369 ! 370 370 ! ! "Poleward" diffusive heat or salt transports (T-S case only) 371 IF( cdtype == 'TRA' .AND. ln_diaptr ) THEN372 371 ! note sign is reversed to give down-gradient diffusive transports (#1043) 373 IF( jn == jp_tem) htr_ldf(:) = ptr_sj( -zftv(:,:,:) ) 374 IF( jn == jp_sal) str_ldf(:) = ptr_sj( -zftv(:,:,:) ) 375 ENDIF 372 IF( cdtype == 'TRA' .AND. ln_diaptr ) CALL dia_ptr_ohst_components( jn, 'ldf', -zftv(:,:,:) ) 376 373 ! 377 374 IF( iom_use("udiff_heattr") .OR. iom_use("vdiff_heattr") ) THEN -
branches/2016/dev_r7233_CMIP6_diags_trunk_version/NEMOGCM/NEMO/OPA_SRC/TRA/traldf_lap_blp.F90
r6140 r7236 150 150 IF( ( kpass == 1 .AND. .NOT.ln_traldf_blp ) .OR. & !== first pass only ( laplacian) ==! 151 151 ( kpass == 2 .AND. ln_traldf_blp ) ) THEN !== 2nd pass only (bilaplacian) ==! 152 IF( cdtype == 'TRA' .AND. ln_diaptr ) THEN 153 IF( jn == jp_tem) htr_ldf(:) = ptr_sj( -ztv(:,:,:) ) 154 IF( jn == jp_sal) str_ldf(:) = ptr_sj( -ztv(:,:,:) ) 155 ENDIF 152 IF( cdtype == 'TRA' .AND. ln_diaptr ) CALL dia_ptr_ohst_components( jn, 'ldf', -ztv(:,:,:) ) 156 153 ENDIF 157 154 ! ! ================== -
branches/2016/dev_r7233_CMIP6_diags_trunk_version/NEMOGCM/NEMO/OPA_SRC/TRA/traldf_triad.F90
r6140 r7236 416 416 ! 417 417 ! ! "Poleward" diffusive heat or salt transports (T-S case only) 418 IF( cdtype == 'TRA' .AND. ln_diaptr ) THEN 419 IF( jn == jp_tem) htr_ldf(:) = ptr_sj( zftv(:,:,:) ) ! 3.3 names 420 IF( jn == jp_sal) str_ldf(:) = ptr_sj( zftv(:,:,:) ) 421 ENDIF 418 IF( cdtype == 'TRA' .AND. ln_diaptr ) CALL dia_ptr_ohst_components( jn, 'ldf', zftv(:,:,:) ) 422 419 ! 423 420 IF( iom_use("udiff_heattr") .OR. iom_use("vdiff_heattr") ) THEN
Note: See TracChangeset
for help on using the changeset viewer.