Changeset 5883 for branches/2015/dev_r5836_NOC3_vvl_by_default/NEMOGCM/NEMO/OPA_SRC/TRA/traadv_fct.F90
- Timestamp:
- 2015-11-13T08:01:08+01:00 (8 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2015/dev_r5836_NOC3_vvl_by_default/NEMOGCM/NEMO/OPA_SRC/TRA/traadv_fct.F90
r5866 r5883 53 53 !! *** ROUTINE tra_adv_fct *** 54 54 !! 55 !! ** Purpose : Compute the now trend due to total advection of 56 !! tracersand add it to the general trend of tracer equations55 !! ** Purpose : Compute the now trend due to total advection of tracers 56 !! and add it to the general trend of tracer equations 57 57 !! 58 58 !! ** Method : - 2nd or 4th FCT scheme on the horizontal direction 59 59 !! (choice through the value of kn_fct) 60 !! - 4th order compact scheme on the vertical60 !! - on the vertical the 4th order is a compact scheme 61 61 !! - corrected flux (monotonic correction) 62 62 !! 63 !! ** Action : - update (pta) with the now advective tracer trends 64 !! - send the trends for further diagnostics 63 !! ** Action : - update pta with the now advective tracer trends 64 !! - send trends to trdtra module for further diagnostcs (l_trdtra=T) 65 !! - htr_adv, str_adv : poleward advective heat and salt transport (ln_diaptr=T) 65 66 !!---------------------------------------------------------------------- 66 67 INTEGER , INTENT(in ) :: kt ! ocean time-step index … … 101 102 ENDIF 102 103 ! 103 ! 104 IF( .NOT.ln_linssh ) zwz(:,:, 1 ) = 0._wp ! except at the surface in linear free surface case104 ! ! surface & bottom value : flux set to zero one for all 105 zwz(:,:, 1 ) = 0._wp 105 106 zwx(:,:,jpk) = 0._wp ; zwy(:,:,jpk) = 0._wp ; zwz(:,:,jpk) = 0._wp 106 107 ! 107 108 zwi(:,:,:) = 0._wp 108 ! ! =========== 109 DO jn = 1, kjpt ! tracer loop 110 ! ! =========== 109 ! 110 DO jn = 1, kjpt !== loop over the tracers ==! 111 111 ! 112 112 ! !== upstream advection with initial mass fluxes & intermediate update ==! … … 126 126 END DO 127 127 ! !* upstream tracer flux in the k direction *! 128 DO jk = 2, jpkm1 128 DO jk = 2, jpkm1 ! Interior value ( multiplied by wmask) 129 129 DO jj = 1, jpj 130 130 DO ji = 1, jpi … … 135 135 END DO 136 136 END DO 137 !138 137 IF( ln_linssh ) THEN ! top ocean value (only in linear free surface as zwz has been w-masked) 139 138 IF( ln_isfcav ) THEN ! top of the ice-shelf cavities and at the ocean surface … … 155 154 ztra = - ( zwx(ji,jj,jk) - zwx(ji-1,jj ,jk ) & 156 155 & + zwy(ji,jj,jk) - zwy(ji ,jj-1,jk ) & 157 & + zwz(ji,jj,jk) - zwz(ji ,jj ,jk+1) ) / ( e1e2t(ji,jj) * e3t_n(ji,jj,jk))156 & + zwz(ji,jj,jk) - zwz(ji ,jj ,jk+1) ) * r1_e1e2t(ji,jj) / e3t_n(ji,jj,jk) 158 157 ! update and guess with monotonic sheme 159 158 !!gm why tmask added in the two following lines ??? the mask is done in tranxt ! … … 174 173 ENDIF 175 174 ! 176 !177 175 ! !== anti-diffusive flux : high order minus low order ==! 178 176 ! 179 SELECT CASE( kn_fct_h ) 180 ! 181 CASE( 2 ) !2nd order centered177 SELECT CASE( kn_fct_h ) !* horizontal anti-diffusive fluxes 178 ! 179 CASE( 2 ) !- 2nd order centered 182 180 DO jk = 1, jpkm1 183 181 DO jj = 1, jpjm1 … … 189 187 END DO 190 188 ! 191 CASE( 4 ) !4th order centered192 zltu(:,:,jpk) = 0._wp 189 CASE( 4 ) !- 4th order centered 190 zltu(:,:,jpk) = 0._wp ! Bottom value : flux set to zero 193 191 zltv(:,:,jpk) = 0._wp 194 DO jk = 1, jpkm1 195 DO jj = 1, jpjm1 ! First derivative (gradient)192 DO jk = 1, jpkm1 ! Laplacian 193 DO jj = 1, jpjm1 ! 1st derivative (gradient) 196 194 DO ji = 1, fs_jpim1 ! vector opt. 197 195 ztu(ji,jj,jk) = ( ptn(ji+1,jj ,jk,jn) - ptn(ji,jj,jk,jn) ) * umask(ji,jj,jk) … … 199 197 END DO 200 198 END DO 201 DO jj = 2, jpjm1 !199 DO jj = 2, jpjm1 ! 2nd derivative * 1/ 6 202 200 DO ji = fs_2, fs_jpim1 ! vector opt. 203 201 zltu(ji,jj,jk) = ( ztu(ji,jj,jk) + ztu(ji-1,jj,jk) ) * r1_6 … … 206 204 END DO 207 205 END DO 208 !209 206 CALL lbc_lnk( zltu, 'T', 1. ) ; CALL lbc_lnk( zltv, 'T', 1. ) ! Lateral boundary cond. (unchanged sgn) 210 207 ! 211 DO jk = 1, jpkm1 208 DO jk = 1, jpkm1 ! Horizontal advective fluxes 212 209 DO jj = 1, jpjm1 213 210 DO ji = 1, fs_jpim1 ! vector opt. … … 221 218 END DO 222 219 ! 223 CASE( 41 ) !4th order centered ==>> !!gm coding attempt need to be tested224 ztu(:,:,jpk) = 0._wp 220 CASE( 41 ) !- 4th order centered ==>> !!gm coding attempt need to be tested 221 ztu(:,:,jpk) = 0._wp ! Bottom value : flux set to zero 225 222 ztv(:,:,jpk) = 0._wp 226 DO jk = 1, jpkm1 ! gradient227 DO jj = 1, jpjm1 ! First derivative (gradient)223 DO jk = 1, jpkm1 ! 1st derivative (gradient) 224 DO jj = 1, jpjm1 228 225 DO ji = 1, fs_jpim1 ! vector opt. 229 226 ztu(ji,jj,jk) = ( ptn(ji+1,jj ,jk,jn) - ptn(ji,jj,jk,jn) ) * umask(ji,jj,jk) … … 234 231 CALL lbc_lnk( ztu, 'U', -1. ) ; CALL lbc_lnk( ztv, 'V', -1. ) ! Lateral boundary cond. (unchanged sgn) 235 232 ! 236 DO jk = 1, jpkm1 233 DO jk = 1, jpkm1 ! Horizontal advective fluxes 237 234 DO jj = 2, jpjm1 238 235 DO ji = 2, fs_jpim1 ! vector opt. … … 250 247 ! 251 248 END SELECT 252 ! !* vertical anti-diffusive fluxes253 SELECT CASE( kn_fct_v ) ! Interior values (w-masked)254 ! 255 CASE( 2 ) !2nd order centered249 ! 250 SELECT CASE( kn_fct_v ) !* vertical anti-diffusive fluxes (w-masked interior values) 251 ! 252 CASE( 2 ) !- 2nd order centered 256 253 DO jk = 2, jpkm1 257 254 DO jj = 2, jpjm1 258 255 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 ! 256 zwz(ji,jj,jk) = ( pwn(ji,jj,jk) * 0.5_wp * ( ptn(ji,jj,jk,jn) + ptn(ji,jj,jk-1,jn) ) & 257 & - zwz(ji,jj,jk) ) * wmask(ji,jj,jk) 258 END DO 259 END DO 260 END DO 261 ! 262 CASE( 4 ) !- 4th order COMPACT 263 CALL interp_4th_cpt( ptn(:,:,:,jn) , ztw ) ! zwt = COMPACT interpolation of T at w-point 269 264 DO jk = 2, jpkm1 270 265 DO jj = 2, jpjm1 … … 276 271 ! 277 272 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 273 IF( ln_linssh ) THEN ! top ocean value: high order = upstream ==>> zwz=0 274 zwz(:,:,1) = 0._wp ! only ocean surface as interior zwz values have been w-masked 275 ENDIF 280 276 ! 281 277 CALL lbc_lnk( zwx, 'U', -1. ) ; CALL lbc_lnk( zwy, 'V', -1. ) ! Lateral bondary conditions 282 278 CALL lbc_lnk( zwz, 'W', 1. ) 283 279 ! 284 280 ! !== monotonicity algorithm ==! 285 281 ! 286 282 CALL nonosc( ptb(:,:,:,jn), zwx, zwy, zwz, zwi, p2dt ) 287 288 283 ! 289 284 ! !== final trend with corrected fluxes ==! 290 285 ! … … 300 295 END DO 301 296 ! 302 IF( l_trd ) THEN 297 IF( l_trd ) THEN ! trend diagnostics (contribution of upstream fluxes) 303 298 ztrdx(:,:,:) = ztrdx(:,:,:) + zwx(:,:,:) ! <<< Add to previously computed 304 299 ztrdy(:,:,:) = ztrdy(:,:,:) + zwy(:,:,:) ! <<< Add to previously computed … … 311 306 CALL wrk_dealloc( jpi,jpj,jpk, ztrdx, ztrdy, ztrdz ) 312 307 END IF 313 ! 308 ! ! "Poleward" heat and salt transports (contribution of upstream fluxes) 314 309 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(:)310 IF( jn == jp_tem ) htr_adv(:) = htr_adv(:) + ptr_sj( zwy(:,:,:) ) 311 IF( jn == jp_sal ) str_adv(:) = str_adv(:) + ptr_sj( zwy(:,:,:) ) 317 312 ENDIF 318 313 ! 319 END DO 314 END DO ! end of tracer loop 320 315 ! 321 316 CALL wrk_dealloc( jpi,jpj,jpk, zwi, zwx, zwy, zwz, ztu, ztv, zltu, zltv, ztw ) … … 392 387 zr_p2dt(:) = 1._wp / p2dt(:) 393 388 ! 389 ! surface & Bottom value : flux set to zero for all tracers 390 zwz(:,:, 1 ) = 0._wp 391 zwx(:,:,jpk) = 0._wp ; zwz(:,:,jpk) = 0._wp 392 zwy(:,:,jpk) = 0._wp ; zwi(:,:,jpk) = 0._wp 393 ! 394 394 ! ! =========== 395 395 DO jn = 1, kjpt ! tracer loop 396 396 ! ! =========== 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 397 ! 398 ! Upstream advection with initial mass fluxes & intermediate update 399 DO jk = 1, jpkm1 ! upstream tracer flux in the i and j direction 406 400 DO jj = 1, jpjm1 407 401 DO ji = 1, fs_jpim1 ! vector opt. … … 416 410 END DO 417 411 END DO 418 419 ! upstream tracer flux in the k direction 420 DO jk = 2, jpkm1 ! Interior value 412 ! ! upstream tracer flux in the k direction 413 DO jk = 2, jpkm1 ! Interior value 421 414 DO jj = 1, jpj 422 415 DO ji = 1, jpi … … 427 420 END DO 428 421 END DO 429 ! ! top value 430 IF( .NOT.ln_linssh ) 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 422 IF( ln_linssh ) THEN ! top value : linear free surface case only (as zwz is multiplied by wmask) 423 IF( ln_isfcav ) THEN ! ice-shelf cavities: top value 434 424 DO jj = 1, jpj 435 425 DO ji = 1, jpi … … 437 427 END DO 438 428 END DO 439 ELSE ! standard case429 ELSE ! no cavities, surface value 440 430 zwz(:,:,1) = pwn(:,:,1) * ptb(:,:,1,jn) 441 431 ENDIF … … 446 436 DO jj = 2, jpjm1 447 437 DO ji = fs_2, fs_jpim1 ! vector opt. 448 ! total intermediate advective trends438 ! ! total intermediate advective trends 449 439 ztra = - ( zwx(ji,jj,jk) - zwx(ji-1,jj ,jk ) & 450 440 & + zwy(ji,jj,jk) - zwy(ji ,jj-1,jk ) & 451 & + zwz(ji,jj,jk) - zwz(ji ,jj ,jk+1) ) / ( e1e2t(ji,jj) * e3t_n(ji,jj,jk))452 ! update and guess with monotonic sheme441 & + zwz(ji,jj,jk) - zwz(ji ,jj ,jk+1) ) * r1_e1e2t(ji,jj) / e3t_n(ji,jj,jk) 442 ! ! update and guess with monotonic sheme 453 443 pta(ji,jj,jk,jn) = pta(ji,jj,jk,jn) + ztra 454 444 zwi(ji,jj,jk) = ( ptb(ji,jj,jk,jn) + z2dtt * ztra ) * tmask(ji,jj,jk) … … 497 487 END DO 498 488 END DO 499 489 ! 500 490 ! !* vertical anti-diffusive flux 501 491 zwz_sav(:,:,:) = zwz(:,:,:) 502 492 ztrs (:,:,:,1) = ptb(:,:,:,jn) 503 493 zwzts (:,:,:) = 0._wp 504 IF( .NOT.ln_linssh ) zwz(:,:, 1 ) = 0._wp ! surface value set to zero in vvl case505 494 ! 506 495 DO jl = 1, kn_fct_zts ! Start of sub timestepping loop … … 535 524 END DO 536 525 END DO 537 ELSE ! standard case526 ELSE ! no ocean cavities 538 527 zwz(:,:,1) = pwn(:,:,1) * ptb(:,:,1,jn) 539 528 ENDIF
Note: See TracChangeset
for help on using the changeset viewer.