- Timestamp:
- 2016-11-28T17:04:10+01:00 (7 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2016/dev_INGV_UKMO_2016/NEMOGCM/NEMO/OPA_SRC/TRA/traadv_fct.F90
r5930 r7351 39 39 40 40 !! * Substitutions 41 # include "domzgr_substitute.h90"42 41 # include "vectopt_loop_substitute.h90" 43 42 !!---------------------------------------------------------------------- … … 53 52 !! *** ROUTINE tra_adv_fct *** 54 53 !! 55 !! ** Purpose : Compute the now trend due to total advection of 56 !! tracersand add it to the general trend of tracer equations54 !! ** Purpose : Compute the now trend due to total advection of tracers 55 !! and add it to the general trend of tracer equations 57 56 !! 58 57 !! ** Method : - 2nd or 4th FCT scheme on the horizontal direction 59 58 !! (choice through the value of kn_fct) 60 !! - 4th order compact scheme on the vertical59 !! - on the vertical the 4th order is a compact scheme 61 60 !! - corrected flux (monotonic correction) 62 61 !! 63 !! ** Action : - update (pta) with the now advective tracer trends 64 !! - send the trends for further diagnostics 62 !! ** Action : - update pta with the now advective tracer trends 63 !! - send trends to trdtra module for further diagnostcs (l_trdtra=T) 64 !! - htr_adv, str_adv : poleward advective heat and salt transport (ln_diaptr=T) 65 65 !!---------------------------------------------------------------------- 66 66 INTEGER , INTENT(in ) :: kt ! ocean time-step index … … 70 70 INTEGER , INTENT(in ) :: kn_fct_h ! order of the FCT scheme (=2 or 4) 71 71 INTEGER , INTENT(in ) :: kn_fct_v ! order of the FCT scheme (=2 or 4) 72 REAL(wp) , DIMENSION( jpk ), INTENT(in ) :: p2dt ! vertical profile oftracer time-step72 REAL(wp) , INTENT(in ) :: p2dt ! tracer time-step 73 73 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ) :: pun, pvn, pwn ! 3 ocean velocity components 74 74 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: ptb, ptn ! before and now tracer fields … … 76 76 ! 77 77 INTEGER :: ji, jj, jk, jn ! dummy loop indices 78 REAL(wp) :: z 2dtt, ztra! local scalar78 REAL(wp) :: ztra ! local scalar 79 79 REAL(wp) :: zfp_ui, zfp_vj, zfp_wk, zC2t_u, zC4t_u ! - - 80 80 REAL(wp) :: zfm_ui, zfm_vj, zfm_wk, zC2t_v, zC4t_v ! - - … … 101 101 ENDIF 102 102 ! 103 ! 104 IF( lk_vvl ) zwz(:,:, 1 ) = 0._wp ! except at the surface in linear free surface case103 ! ! surface & bottom value : flux set to zero one for all 104 zwz(:,:, 1 ) = 0._wp 105 105 zwx(:,:,jpk) = 0._wp ; zwy(:,:,jpk) = 0._wp ; zwz(:,:,jpk) = 0._wp 106 106 ! 107 107 zwi(:,:,:) = 0._wp 108 ! ! =========== 109 DO jn = 1, kjpt ! tracer loop 110 ! ! =========== 108 ! 109 DO jn = 1, kjpt !== loop over the tracers ==! 111 110 ! 112 111 ! !== upstream advection with initial mass fluxes & intermediate update ==! … … 126 125 END DO 127 126 ! !* upstream tracer flux in the k direction *! 128 DO jk = 2, jpkm1 127 DO jk = 2, jpkm1 ! Interior value ( multiplied by wmask) 129 128 DO jj = 1, jpj 130 129 DO ji = 1, jpi … … 135 134 END DO 136 135 END DO 137 ! 138 IF(.NOT.lk_vvl ) THEN ! top ocean value (only in linear free surface as zwz has been w-masked) 136 IF( ln_linssh ) THEN ! top ocean value (only in linear free surface as zwz has been w-masked) 139 137 IF( ln_isfcav ) THEN ! top of the ice-shelf cavities and at the ocean surface 140 138 DO jj = 1, jpj … … 149 147 ! 150 148 DO jk = 1, jpkm1 !* trend and after field with monotonic scheme 151 z2dtt = p2dt(jk)152 149 DO jj = 2, jpjm1 153 150 DO ji = fs_2, fs_jpim1 ! vector opt. 154 ! total intermediate advective trends151 ! ! total intermediate advective trends 155 152 ztra = - ( zwx(ji,jj,jk) - zwx(ji-1,jj ,jk ) & 156 153 & + zwy(ji,jj,jk) - zwy(ji ,jj-1,jk ) & 157 & + zwz(ji,jj,jk) - zwz(ji ,jj ,jk+1) ) / ( e1e2t(ji,jj) * fse3t_n(ji,jj,jk) ) 158 ! update and guess with monotonic sheme 159 !!gm why tmask added in the two following lines ??? the mask is done in tranxt ! 160 pta(ji,jj,jk,jn) = pta(ji,jj,jk,jn) + ztra * tmask(ji,jj,jk) 161 zwi(ji,jj,jk) = ( ptb(ji,jj,jk,jn) + z2dtt * ztra ) * tmask(ji,jj,jk) 154 & + zwz(ji,jj,jk) - zwz(ji ,jj ,jk+1) ) * r1_e1e2t(ji,jj) 155 ! ! update and guess with monotonic sheme 156 pta(ji,jj,jk,jn) = pta(ji,jj,jk,jn) + ztra / e3t_n(ji,jj,jk) * tmask(ji,jj,jk) 157 zwi(ji,jj,jk) = ( e3t_b(ji,jj,jk) * ptb(ji,jj,jk,jn) + p2dt * ztra ) / e3t_a(ji,jj,jk) * tmask(ji,jj,jk) 162 158 END DO 163 159 END DO … … 166 162 ! 167 163 IF( l_trd ) THEN ! trend diagnostics (contribution of upstream fluxes) 168 ztrdx(:,:,:) = zwx(:,:,:) ; ztrdy(:,:,:) = zwy(:,:,:); ztrdz(:,:,:) = zwz(:,:,:)164 ztrdx(:,:,:) = zwx(:,:,:) ; ztrdy(:,:,:) = zwy(:,:,:) ; ztrdz(:,:,:) = zwz(:,:,:) 169 165 END IF 170 166 ! ! "Poleward" heat and salt transports (contribution of upstream fluxes) … … 174 170 ENDIF 175 171 ! 176 !177 172 ! !== anti-diffusive flux : high order minus low order ==! 178 173 ! 179 SELECT CASE( kn_fct_h ) 180 ! 181 CASE( 2 ) !2nd order centered174 SELECT CASE( kn_fct_h ) !* horizontal anti-diffusive fluxes 175 ! 176 CASE( 2 ) !- 2nd order centered 182 177 DO jk = 1, jpkm1 183 178 DO jj = 1, jpjm1 … … 189 184 END DO 190 185 ! 191 CASE( 4 ) !4th order centered192 zltu(:,:,jpk) = 0._wp 186 CASE( 4 ) !- 4th order centered 187 zltu(:,:,jpk) = 0._wp ! Bottom value : flux set to zero 193 188 zltv(:,:,jpk) = 0._wp 194 DO jk = 1, jpkm1 195 DO jj = 1, jpjm1 ! First derivative (gradient)189 DO jk = 1, jpkm1 ! Laplacian 190 DO jj = 1, jpjm1 ! 1st derivative (gradient) 196 191 DO ji = 1, fs_jpim1 ! vector opt. 197 192 ztu(ji,jj,jk) = ( ptn(ji+1,jj ,jk,jn) - ptn(ji,jj,jk,jn) ) * umask(ji,jj,jk) … … 199 194 END DO 200 195 END DO 201 DO jj = 2, jpjm1 !196 DO jj = 2, jpjm1 ! 2nd derivative * 1/ 6 202 197 DO ji = fs_2, fs_jpim1 ! vector opt. 203 198 zltu(ji,jj,jk) = ( ztu(ji,jj,jk) + ztu(ji-1,jj,jk) ) * r1_6 … … 206 201 END DO 207 202 END DO 208 !209 203 CALL lbc_lnk( zltu, 'T', 1. ) ; CALL lbc_lnk( zltv, 'T', 1. ) ! Lateral boundary cond. (unchanged sgn) 210 204 ! 211 DO jk = 1, jpkm1 205 DO jk = 1, jpkm1 ! Horizontal advective fluxes 212 206 DO jj = 1, jpjm1 213 207 DO ji = 1, fs_jpim1 ! vector opt. … … 221 215 END DO 222 216 ! 223 CASE( 41 ) !4th order centered ==>> !!gm coding attempt need to be tested224 ztu(:,:,jpk) = 0._wp 217 CASE( 41 ) !- 4th order centered ==>> !!gm coding attempt need to be tested 218 ztu(:,:,jpk) = 0._wp ! Bottom value : flux set to zero 225 219 ztv(:,:,jpk) = 0._wp 226 DO jk = 1, jpkm1 ! gradient227 DO jj = 1, jpjm1 ! First derivative (gradient)220 DO jk = 1, jpkm1 ! 1st derivative (gradient) 221 DO jj = 1, jpjm1 228 222 DO ji = 1, fs_jpim1 ! vector opt. 229 223 ztu(ji,jj,jk) = ( ptn(ji+1,jj ,jk,jn) - ptn(ji,jj,jk,jn) ) * umask(ji,jj,jk) … … 234 228 CALL lbc_lnk( ztu, 'U', -1. ) ; CALL lbc_lnk( ztv, 'V', -1. ) ! Lateral boundary cond. (unchanged sgn) 235 229 ! 236 DO jk = 1, jpkm1 230 DO jk = 1, jpkm1 ! Horizontal advective fluxes 237 231 DO jj = 2, jpjm1 238 232 DO ji = 2, fs_jpim1 ! vector opt. … … 250 244 ! 251 245 END SELECT 252 ! !* vertical anti-diffusive fluxes253 SELECT CASE( kn_fct_v ) ! Interior values (w-masked)254 ! 255 CASE( 2 ) !2nd order centered246 ! 247 SELECT CASE( kn_fct_v ) !* vertical anti-diffusive fluxes (w-masked interior values) 248 ! 249 CASE( 2 ) !- 2nd order centered 256 250 DO jk = 2, jpkm1 257 251 DO jj = 2, jpjm1 258 252 DO ji = fs_2, fs_jpim1 259 zwz(ji,jj,jk) = ( 0.5_wp * pwn(ji,jj,jk) * ( ptn(ji,jj,jk,jn) + ptn(ji,jj,jk-1,jn) ) & 260 - zwz(ji,jj,jk) ) * wmask(ji,jj,jk) 261 END DO 262 END DO 263 END DO 264 ! 265 CASE( 4 ) ! 4th order COMPACT 266 ! 267 CALL interp_4th_cpt( ptn(:,:,:,jn) , ztw ) ! COMPACT interpolation of T at w-point 268 ! 253 zwz(ji,jj,jk) = ( pwn(ji,jj,jk) * 0.5_wp * ( ptn(ji,jj,jk,jn) + ptn(ji,jj,jk-1,jn) ) & 254 & - zwz(ji,jj,jk) ) * wmask(ji,jj,jk) 255 END DO 256 END DO 257 END DO 258 ! 259 CASE( 4 ) !- 4th order COMPACT 260 CALL interp_4th_cpt( ptn(:,:,:,jn) , ztw ) ! zwt = COMPACT interpolation of T at w-point 269 261 DO jk = 2, jpkm1 270 262 DO jj = 2, jpjm1 … … 276 268 ! 277 269 END SELECT 278 ! ! top ocean value: high order = upstream ==>> zwz=0 279 zwz(:,:, 1 ) = 0._wp ! only ocean surface as interior zwz values have been w-masked 270 IF( ln_linssh ) THEN ! top ocean value: high order = upstream ==>> zwz=0 271 zwz(:,:,1) = 0._wp ! only ocean surface as interior zwz values have been w-masked 272 ENDIF 280 273 ! 281 274 CALL lbc_lnk( zwx, 'U', -1. ) ; CALL lbc_lnk( zwy, 'V', -1. ) ! Lateral bondary conditions 282 275 CALL lbc_lnk( zwz, 'W', 1. ) 283 276 ! 284 277 ! !== monotonicity algorithm ==! 285 278 ! 286 279 CALL nonosc( ptb(:,:,:,jn), zwx, zwy, zwz, zwi, p2dt ) 287 288 280 ! 289 281 ! !== final trend with corrected fluxes ==! 290 282 ! … … 295 287 & + zwy(ji,jj,jk) - zwy(ji ,jj-1,jk ) & 296 288 & + zwz(ji,jj,jk) - zwz(ji ,jj ,jk+1) ) & 297 & / ( e1e2t(ji,jj) * fse3t(ji,jj,jk))298 END DO 299 END DO 300 END DO 301 ! 302 IF( l_trd ) THEN 289 & * r1_e1e2t(ji,jj) / e3t_n(ji,jj,jk) 290 END DO 291 END DO 292 END DO 293 ! 294 IF( l_trd ) THEN ! trend diagnostics (contribution of upstream fluxes) 303 295 ztrdx(:,:,:) = ztrdx(:,:,:) + zwx(:,:,:) ! <<< Add to previously computed 304 296 ztrdy(:,:,:) = ztrdy(:,:,:) + zwy(:,:,:) ! <<< Add to previously computed … … 311 303 CALL wrk_dealloc( jpi,jpj,jpk, ztrdx, ztrdy, ztrdz ) 312 304 END IF 313 ! 305 ! ! "Poleward" heat and salt transports (contribution of upstream fluxes) 314 306 IF( cdtype == 'TRA' .AND. ln_diaptr ) THEN 315 IF( jn == jp_tem ) htr_adv(:) = ptr_sj( zwy(:,:,:) ) + htr_adv(:)316 IF( jn == jp_sal ) str_adv(:) = ptr_sj( zwy(:,:,:) ) + str_adv(:)307 IF( jn == jp_tem ) htr_adv(:) = htr_adv(:) + ptr_sj( zwy(:,:,:) ) 308 IF( jn == jp_sal ) str_adv(:) = str_adv(:) + ptr_sj( zwy(:,:,:) ) 317 309 ENDIF 318 310 ! 319 END DO 311 END DO ! end of tracer loop 320 312 ! 321 313 CALL wrk_dealloc( jpi,jpj,jpk, zwi, zwx, zwy, zwz, ztu, ztv, zltu, zltv, ztw ) … … 348 340 INTEGER , INTENT(in ) :: kjpt ! number of tracers 349 341 INTEGER , INTENT(in ) :: kn_fct_zts ! number of number of vertical sub-timesteps 350 REAL(wp) , DIMENSION( jpk ), INTENT(in ) :: p2dt ! vertical profile oftracer time-step342 REAL(wp) , INTENT(in ) :: p2dt ! tracer time-step 351 343 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ) :: pun, pvn, pwn ! 3 ocean velocity components 352 344 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: ptb, ptn ! before and now tracer fields … … 354 346 ! 355 347 REAL(wp), DIMENSION( jpk ) :: zts ! length of sub-timestep for vertical advection 356 REAL(wp) , DIMENSION( jpk ):: zr_p2dt ! reciprocal of tracer timestep348 REAL(wp) :: zr_p2dt ! reciprocal of tracer timestep 357 349 INTEGER :: ji, jj, jk, jl, jn ! dummy loop indices 358 350 INTEGER :: jtb, jtn, jta ! sub timestep pointers for leap-frog/euler forward steps 359 351 INTEGER :: jtaken ! toggle for collecting appropriate fluxes from sub timesteps 360 352 REAL(wp) :: z_rzts ! Fractional length of Euler forward sub-timestep for vertical advection 361 REAL(wp) :: z 2dtt, ztra! local scalar353 REAL(wp) :: ztra ! local scalar 362 354 REAL(wp) :: zfp_ui, zfp_vj, zfp_wk ! - - 363 355 REAL(wp) :: zfm_ui, zfm_vj, zfm_wk ! - - … … 371 363 ! 372 364 CALL wrk_alloc( jpi,jpj, zwx_sav, zwy_sav ) 373 CALL wrk_alloc( jpi,jpj, jpk,zwx, zwy, zwz, zwi, zhdiv, zwzts, zwz_sav )365 CALL wrk_alloc( jpi,jpj,jpk, zwx, zwy, zwz, zwi, zhdiv, zwzts, zwz_sav ) 374 366 CALL wrk_alloc( jpi,jpj,jpk,kjpt+1, ztrs ) 375 367 ! … … 390 382 zwi(:,:,:) = 0._wp 391 383 z_rzts = 1._wp / REAL( kn_fct_zts, wp ) 392 zr_p2dt(:) = 1._wp / p2dt(:) 384 zr_p2dt = 1._wp / p2dt 385 ! 386 ! surface & Bottom value : flux set to zero for all tracers 387 zwz(:,:, 1 ) = 0._wp 388 zwx(:,:,jpk) = 0._wp ; zwz(:,:,jpk) = 0._wp 389 zwy(:,:,jpk) = 0._wp ; zwi(:,:,jpk) = 0._wp 393 390 ! 394 391 ! ! =========== 395 392 DO jn = 1, kjpt ! tracer loop 396 393 ! ! =========== 397 ! 1. Bottom value : flux set to zero 398 ! ---------------------------------- 399 zwx(:,:,jpk) = 0._wp ; zwz(:,:,jpk) = 0._wp 400 zwy(:,:,jpk) = 0._wp ; zwi(:,:,jpk) = 0._wp 401 402 ! 2. upstream advection with initial mass fluxes & intermediate update 403 ! -------------------------------------------------------------------- 404 ! upstream tracer flux in the i and j direction 405 DO jk = 1, jpkm1 394 ! 395 ! Upstream advection with initial mass fluxes & intermediate update 396 DO jk = 1, jpkm1 ! upstream tracer flux in the i and j direction 406 397 DO jj = 1, jpjm1 407 398 DO ji = 1, fs_jpim1 ! vector opt. … … 416 407 END DO 417 408 END DO 418 419 ! upstream tracer flux in the k direction 420 DO jk = 2, jpkm1 ! Interior value 409 ! ! upstream tracer flux in the k direction 410 DO jk = 2, jpkm1 ! Interior value 421 411 DO jj = 1, jpj 422 412 DO ji = 1, jpi … … 427 417 END DO 428 418 END DO 429 ! ! top value 430 IF( lk_vvl ) THEN ! variable volume: only k=1 as zwz is multiplied by wmask 431 zwz(:,:, 1 ) = 0._wp 432 ELSE ! linear free surface 433 IF( ln_isfcav ) THEN ! ice-shelf cavities 419 IF( ln_linssh ) THEN ! top value : linear free surface case only (as zwz is multiplied by wmask) 420 IF( ln_isfcav ) THEN ! ice-shelf cavities: top value 434 421 DO jj = 1, jpj 435 422 DO ji = 1, jpi … … 437 424 END DO 438 425 END DO 439 ELSE ! standard case426 ELSE ! no cavities, surface value 440 427 zwz(:,:,1) = pwn(:,:,1) * ptb(:,:,1,jn) 441 428 ENDIF … … 443 430 ! 444 431 DO jk = 1, jpkm1 ! total advective trend 445 z2dtt = p2dt(jk)446 432 DO jj = 2, jpjm1 447 433 DO ji = fs_2, fs_jpim1 ! vector opt. 448 ! total intermediate advective trends434 ! ! total intermediate advective trends 449 435 ztra = - ( zwx(ji,jj,jk) - zwx(ji-1,jj ,jk ) & 450 436 & + zwy(ji,jj,jk) - zwy(ji ,jj-1,jk ) & 451 & + zwz(ji,jj,jk) - zwz(ji ,jj ,jk+1) ) / ( e1e2t(ji,jj) * fse3t_n(ji,jj,jk))452 ! update and guess with monotonic sheme453 pta(ji,jj,jk,jn) = pta(ji,jj,jk,jn) + ztra454 zwi(ji,jj,jk) = ( ptb(ji,jj,jk,jn) + z2dtt * ztra) * tmask(ji,jj,jk)437 & + zwz(ji,jj,jk) - zwz(ji ,jj ,jk+1) ) * r1_e1e2t(ji,jj) 438 ! ! update and guess with monotonic sheme 439 pta(ji,jj,jk,jn) = pta(ji,jj,jk,jn) + ztra / e3t_n(ji,jj,jk) * tmask(ji,jj,jk) 440 zwi(ji,jj,jk) = ( e3t_b(ji,jj,jk) * ptb(ji,jj,jk,jn) + p2dt * ztra ) / e3t_a(ji,jj,jk) * tmask(ji,jj,jk) 455 441 END DO 456 442 END DO … … 497 483 END DO 498 484 END DO 499 485 ! 500 486 ! !* vertical anti-diffusive flux 501 487 zwz_sav(:,:,:) = zwz(:,:,:) 502 488 ztrs (:,:,:,1) = ptb(:,:,:,jn) 489 ztrs (:,:,1,2) = ptb(:,:,1,jn) 490 ztrs (:,:,1,3) = ptb(:,:,1,jn) 503 491 zwzts (:,:,:) = 0._wp 504 IF( lk_vvl ) zwz(:,:, 1 ) = 0._wp ! surface value set to zero in vvl case505 492 ! 506 493 DO jl = 1, kn_fct_zts ! Start of sub timestepping loop … … 508 495 IF( jl == 1 ) THEN ! Euler forward to kick things off 509 496 jtb = 1 ; jtn = 1 ; jta = 2 510 zts(:) = p2dt (:)* z_rzts497 zts(:) = p2dt * z_rzts 511 498 jtaken = MOD( kn_fct_zts + 1 , 2) ! Toggle to collect every second flux 512 499 ! ! starting at jl =1 if kn_fct_zts is odd; … … 514 501 ELSEIF( jl == 2 ) THEN ! First leapfrog step 515 502 jtb = 1 ; jtn = 2 ; jta = 3 516 zts(:) = 2._wp * p2dt (:)* z_rzts503 zts(:) = 2._wp * p2dt * z_rzts 517 504 ELSE ! Shuffle pointers for subsequent leapfrog steps 518 505 jtb = MOD(jtb,3) + 1 … … 528 515 END DO 529 516 END DO 530 IF( .NOT.lk_vvl) THEN ! top value (only in linear free surface case)517 IF( ln_linssh ) THEN ! top value (only in linear free surface case) 531 518 IF( ln_isfcav ) THEN ! ice-shelf cavities 532 519 DO jj = 1, jpj … … 535 522 END DO 536 523 END DO 537 ELSE ! standard case524 ELSE ! no ocean cavities 538 525 zwz(:,:,1) = pwn(:,:,1) * ptb(:,:,1,jn) 539 526 ENDIF … … 547 534 ztrs(ji,jj,jk,jta) = ztrs(ji,jj,jk,jtb) & 548 535 & - zts(jk) * ( zhdiv(ji,jj,jk) + zwz(ji,jj,jk) - zwz(ji,jj,jk+1) ) & 549 & / ( e1e2t(ji,jj) * fse3t(ji,jj,jk))536 & * r1_e1e2t(ji,jj) / e3t_n(ji,jj,jk) 550 537 END DO 551 538 END DO … … 557 544 DO jj = 2, jpjm1 558 545 DO ji = fs_2, fs_jpim1 559 zwz(ji,jj,jk) = ( zwzts(ji,jj,jk) * zr_p2dt (jk)- zwz_sav(ji,jj,jk) ) * wmask(ji,jj,jk)546 zwz(ji,jj,jk) = ( zwzts(ji,jj,jk) * zr_p2dt - zwz_sav(ji,jj,jk) ) * wmask(ji,jj,jk) 560 547 END DO 561 548 END DO … … 576 563 pta(ji,jj,jk,jn) = pta(ji,jj,jk,jn) + ( zwy(ji,jj,jk) - zwy(ji ,jj-1,jk ) & 577 564 & + zwz(ji,jj,jk) - zwz(ji ,jj ,jk+1) ) & 578 & / ( e1e2t(ji,jj) * fse3t_n(ji,jj,jk))565 & * r1_e1e2t(ji,jj) / e3t_n(ji,jj,jk) 579 566 END DO 580 567 END DO … … 623 610 !! in-space based differencing for fluid 624 611 !!---------------------------------------------------------------------- 625 REAL(wp) , DIMENSION(jpk) , INTENT(in ) :: p2dt ! vertical profile oftracer time-step612 REAL(wp) , INTENT(in ) :: p2dt ! tracer time-step 626 613 REAL(wp), DIMENSION (jpi,jpj,jpk), INTENT(in ) :: pbef, paft ! before & after field 627 614 REAL(wp), DIMENSION (jpi,jpj,jpk), INTENT(inout) :: paa, pbb, pcc ! monotonic fluxes in the 3 directions … … 629 616 INTEGER :: ji, jj, jk ! dummy loop indices 630 617 INTEGER :: ikm1 ! local integer 631 REAL(wp) :: zpos, zneg, zbt, za, zb, zc, zbig, zrtrn , z2dtt! local scalars618 REAL(wp) :: zpos, zneg, zbt, za, zb, zc, zbig, zrtrn ! local scalars 632 619 REAL(wp) :: zau, zbu, zcu, zav, zbv, zcv, zup, zdo ! - - 633 620 REAL(wp), POINTER, DIMENSION(:,:,:) :: zbetup, zbetdo, zbup, zbdo … … 652 639 DO jk = 1, jpkm1 653 640 ikm1 = MAX(jk-1,1) 654 z2dtt = p2dt(jk)655 641 DO jj = 2, jpjm1 656 642 DO ji = fs_2, fs_jpim1 ! vector opt. … … 679 665 680 666 ! up & down beta terms 681 zbt = e1 t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) / z2dtt667 zbt = e1e2t(ji,jj) * e3t_n(ji,jj,jk) / p2dt 682 668 zbetup(ji,jj,jk) = ( zup - paft(ji,jj,jk) ) / ( zpos + zrtrn ) * zbt 683 669 zbetdo(ji,jj,jk) = ( paft(ji,jj,jk) - zdo ) / ( zneg + zrtrn ) * zbt
Note: See TracChangeset
for help on using the changeset viewer.