Changeset 6225 for branches/2014/dev_r4704_NOC5_MPP_BDY_UPDATE/NEMOGCM/NEMO/OPA_SRC/TRA/traadv_qck.F90
- Timestamp:
- 2016-01-08T10:35:19+01:00 (8 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2014/dev_r4704_NOC5_MPP_BDY_UPDATE/NEMOGCM/NEMO/OPA_SRC/TRA/traadv_qck.F90
r4499 r6225 17 17 USE oce ! ocean dynamics and active tracers 18 18 USE dom_oce ! ocean space and time domain 19 USE trdmod_oce ! ocean space and time domain 20 USE trdtra ! ocean tracers trends 21 USE trabbl ! advective term in the BBL 19 USE trc_oce ! share passive tracers/Ocean variables 20 USE trd_oce ! trends: ocean variables 21 USE trdtra ! trends manager: tracers 22 USE diaptr ! poleward transport diagnostics 23 ! 22 24 USE lib_mpp ! distribued memory computing 23 25 USE lbclnk ! ocean lateral boundary condition (or mpp link) 24 USE dynspg_oce ! surface pressure gradient variables25 26 USE in_out_manager ! I/O manager 26 USE diaptr ! poleward transport diagnostics27 USE trc_oce ! share passive tracers/Ocean variables28 27 USE wrk_nemo ! Memory Allocation 29 28 USE timing ! Timing … … 39 38 40 39 !! * Substitutions 41 # include "domzgr_substitute.h90"42 40 # include "vectopt_loop_substitute.h90" 43 41 !!---------------------------------------------------------------------- … … 79 77 !! prevent the appearance of spurious numerical oscillations 80 78 !! 81 !! ** Action : - update (pta) with the now advective tracer trends 82 !! - save the trends 79 !! ** Action : - update pta with the now advective tracer trends 80 !! - send trends to trdtra module for further diagnostcs (l_trdtra=T) 81 !! - htr_adv, str_adv : poleward advective heat and salt transport (ln_diaptr=T) 83 82 !! 84 83 !! ** Reference : Leonard (1979, 1991) … … 88 87 CHARACTER(len=3) , INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator) 89 88 INTEGER , INTENT(in ) :: kjpt ! number of tracers 90 REAL(wp) , DIMENSION( jpk ), INTENT(in ) :: p2dt ! vertical profile oftracer time-step89 REAL(wp) , INTENT(in ) :: p2dt ! tracer time-step 91 90 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ) :: pun, pvn, pwn ! 3 ocean velocity components 92 91 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: ptb, ptn ! before and now tracer fields 93 92 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pta ! tracer trend 94 93 !!---------------------------------------------------------------------- 95 96 94 ! 97 95 IF( nn_timing == 1 ) CALL timing_start('tra_adv_qck') … … 105 103 ! 106 104 l_trd = .FALSE. 107 IF( ( cdtype == 'TRA' .AND. l_trdtra ) .OR. ( cdtype == 'TRC' .AND. l_trdtrc ) ) l_trd = .TRUE.108 109 ! I. Thehorizontal fluxes are computed with the QUICKEST + ULTIMATE scheme105 IF( ( cdtype == 'TRA' .AND. l_trdtra ) .OR. ( cdtype == 'TRC' .AND. l_trdtrc ) ) l_trd = .TRUE. 106 ! 107 ! ! horizontal fluxes are computed with the QUICKEST + ULTIMATE scheme 110 108 CALL tra_adv_qck_i( kt, cdtype, p2dt, pun, ptb, ptn, pta, kjpt ) 111 109 CALL tra_adv_qck_j( kt, cdtype, p2dt, pvn, ptb, ptn, pta, kjpt ) 112 110 113 ! II. Thevertical fluxes are computed with the 2nd order centered scheme111 ! ! vertical fluxes are computed with the 2nd order centered scheme 114 112 CALL tra_adv_cen2_k( kt, cdtype, pwn, ptn, pta, kjpt ) 115 113 ! … … 124 122 !! 125 123 !!---------------------------------------------------------------------- 126 USE oce , ONLY: zwx => ua ! ua used as workspace127 !128 124 INTEGER , INTENT(in ) :: kt ! ocean time-step index 129 125 CHARACTER(len=3) , INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator) 130 126 INTEGER , INTENT(in ) :: kjpt ! number of tracers 131 REAL(wp) , DIMENSION( jpk ), INTENT(in ) :: p2dt ! vertical profile oftracer time-step127 REAL(wp) , INTENT(in ) :: p2dt ! tracer time-step 132 128 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ) :: pun ! i-velocity components 133 129 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: ptb, ptn ! before and now tracer fields 134 130 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pta ! tracer trend 135 131 !! 136 INTEGER :: ji, jj, jk, jn ! dummy loop indices137 REAL(wp) :: ztra, zbtr, zdir, zdx, zdt, zmsk ! local scalars138 REAL(wp), POINTER, DIMENSION(:,:,:) :: zfu, zfc, zfd132 INTEGER :: ji, jj, jk, jn ! dummy loop indices 133 REAL(wp) :: ztra, zbtr, zdir, zdx, zmsk ! local scalars 134 REAL(wp), POINTER, DIMENSION(:,:,:) :: zwx, zfu, zfc, zfd 139 135 !---------------------------------------------------------------------- 140 136 ! 141 CALL wrk_alloc( jpi, jpj, jpk, z fu, zfc, zfd )137 CALL wrk_alloc( jpi, jpj, jpk, zwx, zfu, zfc, zfd ) 142 138 ! ! =========== 143 139 DO jn = 1, kjpt ! tracer loop 144 140 ! ! =========== 145 zfu(:,:,:) = 0.0 ; zfc(:,:,:) = 0.0 146 zfd(:,:,:) = 0.0 ; zwx(:,:,:) = 0.0 147 ! 148 DO jk = 1, jpkm1 149 ! 150 !--- Computation of the ustream and downstream value of the tracer and the mask 141 zfu(:,:,:) = 0._wp ; zfc(:,:,:) = 0._wp 142 zfd(:,:,:) = 0._wp ; zwx(:,:,:) = 0._wp 143 ! 144 !!gm why not using a SHIFT instruction... 145 DO jk = 1, jpkm1 !--- Computation of the ustream and downstream value of the tracer and the mask 151 146 DO jj = 2, jpjm1 152 147 DO ji = fs_2, fs_jpim1 ! vector opt. 153 ! Upstream in the x-direction for the tracer 154 zfc(ji,jj,jk) = ptb(ji-1,jj,jk,jn) 155 ! Downstream in the x-direction for the tracer 156 zfd(ji,jj,jk) = ptb(ji+1,jj,jk,jn) 148 zfc(ji,jj,jk) = ptb(ji-1,jj,jk,jn) ! Upstream in the x-direction for the tracer 149 zfd(ji,jj,jk) = ptb(ji+1,jj,jk,jn) ! Downstream in the x-direction for the tracer 157 150 END DO 158 151 END DO … … 163 156 ! Horizontal advective fluxes 164 157 ! --------------------------- 165 !166 158 DO jk = 1, jpkm1 167 159 DO jj = 2, jpjm1 … … 174 166 ! 175 167 DO jk = 1, jpkm1 176 zdt = p2dt(jk)177 168 DO jj = 2, jpjm1 178 169 DO ji = fs_2, fs_jpim1 ! vector opt. 179 170 zdir = 0.5 + SIGN( 0.5, pun(ji,jj,jk) ) ! if pun > 0 : zdir = 1 otherwise zdir = 0 180 zdx = ( zdir * e1t(ji,jj) + ( 1. - zdir ) * e1t(ji+1,jj) ) * e2u(ji,jj) * fse3u(ji,jj,jk)181 zwx(ji,jj,jk) = ABS( pun(ji,jj,jk) ) * zdt / zdx ! (0<zc_cfl<1 : Courant number on x-direction)171 zdx = ( zdir * e1t(ji,jj) + ( 1. - zdir ) * e1t(ji+1,jj) ) * e2u(ji,jj) * e3u_n(ji,jj,jk) 172 zwx(ji,jj,jk) = ABS( pun(ji,jj,jk) ) * p2dt / zdx ! (0<zc_cfl<1 : Courant number on x-direction) 182 173 zfc(ji,jj,jk) = zdir * ptb(ji ,jj,jk,jn) + ( 1. - zdir ) * ptb(ji+1,jj,jk,jn) ! FC in the x-direction for T 183 174 zfd(ji,jj,jk) = zdir * ptb(ji+1,jj,jk,jn) + ( 1. - zdir ) * ptb(ji ,jj,jk,jn) ! FD in the x-direction for T … … 224 215 DO jj = 2, jpjm1 225 216 DO ji = fs_2, fs_jpim1 ! vector opt. 226 zbtr = 1. / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk))217 zbtr = r1_e1e2t(ji,jj) / e3t_n(ji,jj,jk) 227 218 ! horizontal advective trends 228 219 ztra = - zbtr * ( zwx(ji,jj,jk) - zwx(ji-1,jj,jk) ) … … 232 223 END DO 233 224 END DO 234 ! ! trend diagnostics (contribution of upstream fluxes)235 IF( l_trd ) CALL trd_tra( kt, cdtype, jn, jptra_trd_xad, zwx, pun, ptn(:,:,:,jn) )225 ! ! trend diagnostics 226 IF( l_trd ) CALL trd_tra( kt, cdtype, jn, jptra_xad, zwx, pun, ptn(:,:,:,jn) ) 236 227 ! 237 228 END DO 238 229 ! 239 CALL wrk_dealloc( jpi, jpj, jpk, z fu, zfc, zfd )230 CALL wrk_dealloc( jpi, jpj, jpk, zwx, zfu, zfc, zfd ) 240 231 ! 241 232 END SUBROUTINE tra_adv_qck_i … … 247 238 !! 248 239 !!---------------------------------------------------------------------- 249 USE oce , ONLY: zwy => ua ! ua used as workspace250 !251 240 INTEGER , INTENT(in ) :: kt ! ocean time-step index 252 241 CHARACTER(len=3) , INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator) 253 242 INTEGER , INTENT(in ) :: kjpt ! number of tracers 254 REAL(wp) , DIMENSION( jpk ), INTENT(in ) :: p2dt ! vertical profile oftracer time-step243 REAL(wp) , INTENT(in ) :: p2dt ! tracer time-step 255 244 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ) :: pvn ! j-velocity components 256 245 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: ptb, ptn ! before and now tracer fields … … 258 247 !! 259 248 INTEGER :: ji, jj, jk, jn ! dummy loop indices 260 REAL(wp) :: ztra, zbtr, zdir, zdx, z dt, zmsk ! local scalars261 REAL(wp), POINTER, DIMENSION(:,:,:) :: z fu, zfc, zfd249 REAL(wp) :: ztra, zbtr, zdir, zdx, zmsk ! local scalars 250 REAL(wp), POINTER, DIMENSION(:,:,:) :: zwy, zfu, zfc, zfd 262 251 !---------------------------------------------------------------------- 263 252 ! 264 CALL wrk_alloc( jpi, jpj, jpk, z fu, zfc, zfd )253 CALL wrk_alloc( jpi, jpj, jpk, zwy, zfu, zfc, zfd ) 265 254 ! 266 255 ! ! =========== … … 299 288 ! 300 289 DO jk = 1, jpkm1 301 zdt = p2dt(jk)302 290 DO jj = 2, jpjm1 303 291 DO ji = fs_2, fs_jpim1 ! vector opt. 304 292 zdir = 0.5 + SIGN( 0.5, pvn(ji,jj,jk) ) ! if pun > 0 : zdir = 1 otherwise zdir = 0 305 zdx = ( zdir * e2t(ji,jj) + ( 1. - zdir ) * e2t(ji,jj+1) ) * e1v(ji,jj) * fse3v(ji,jj,jk)306 zwy(ji,jj,jk) = ABS( pvn(ji,jj,jk) ) * zdt / zdx ! (0<zc_cfl<1 : Courant number on x-direction)293 zdx = ( zdir * e2t(ji,jj) + ( 1. - zdir ) * e2t(ji,jj+1) ) * e1v(ji,jj) * e3v_n(ji,jj,jk) 294 zwy(ji,jj,jk) = ABS( pvn(ji,jj,jk) ) * p2dt / zdx ! (0<zc_cfl<1 : Courant number on x-direction) 307 295 zfc(ji,jj,jk) = zdir * ptb(ji,jj ,jk,jn) + ( 1. - zdir ) * ptb(ji,jj+1,jk,jn) ! FC in the x-direction for T 308 296 zfd(ji,jj,jk) = zdir * ptb(ji,jj+1,jk,jn) + ( 1. - zdir ) * ptb(ji,jj ,jk,jn) ! FD in the x-direction for T … … 350 338 DO jj = 2, jpjm1 351 339 DO ji = fs_2, fs_jpim1 ! vector opt. 352 zbtr = 1. / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk))340 zbtr = r1_e1e2t(ji,jj) / e3t_n(ji,jj,jk) 353 341 ! horizontal advective trends 354 342 ztra = - zbtr * ( zwy(ji,jj,jk) - zwy(ji,jj-1,jk) ) … … 358 346 END DO 359 347 END DO 360 ! ! trend diagnostics (contribution of upstream fluxes)361 IF( l_trd ) CALL trd_tra( kt, cdtype, jn, jptra_trd_yad, zwy, pvn, ptn(:,:,:,jn) )348 ! ! trend diagnostics 349 IF( l_trd ) CALL trd_tra( kt, cdtype, jn, jptra_yad, zwy, pvn, ptn(:,:,:,jn) ) 362 350 ! ! "Poleward" heat and salt transports (contribution of upstream fluxes) 363 IF( cdtype == 'TRA' .AND. ln_diaptr .AND. ( MOD( kt, nn_fptr ) == 0 )) THEN364 IF( jn == jp_tem ) htr_adv(:) = ptr_ vj( zwy(:,:,:) )365 IF( jn == jp_sal ) str_adv(:) = ptr_ vj( zwy(:,:,:) )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(:,:,:) ) 366 354 ENDIF 367 355 ! 368 356 END DO 369 357 ! 370 CALL wrk_dealloc( jpi, jpj, jpk, z fu, zfc, zfd )358 CALL wrk_dealloc( jpi, jpj, jpk, zwy, zfu, zfc, zfd ) 371 359 ! 372 360 END SUBROUTINE tra_adv_qck_j … … 378 366 !! 379 367 !!---------------------------------------------------------------------- 380 USE oce, ONLY: zwz => ua ! ua used as workspace381 !382 368 INTEGER , INTENT(in ) :: kt ! ocean time-step index 383 369 CHARACTER(len=3) , INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator) … … 388 374 ! 389 375 INTEGER :: ji, jj, jk, jn ! dummy loop indices 390 REAL(wp) :: zbtr , ztra ! local scalars 391 !!---------------------------------------------------------------------- 392 376 REAL(wp), POINTER, DIMENSION(:,:,:) :: zwz 377 !!---------------------------------------------------------------------- 378 ! 379 CALL wrk_alloc( jpi,jpj,jpk, zwz ) 380 ! 381 zwz(:,:, 1 ) = 0._wp ! surface & bottom values set to zero for all tracers 382 zwz(:,:,jpk) = 0._wp 383 ! 393 384 ! ! =========== 394 385 DO jn = 1, kjpt ! tracer loop 395 386 ! ! =========== 396 ! 1. Bottom value : flux set to zero 397 zwz(:,:,jpk) = 0.e0 ! Bottom value : flux set to zero 398 ! 399 ! ! Surface value 400 IF( lk_vvl ) THEN ; zwz(:,:, 1 ) = 0.e0 ! Variable volume : flux set to zero 401 ELSE ; zwz(:,:, 1 ) = pwn(:,:,1) * ptn(:,:,1,jn) ! Constant volume : advective flux through the surface 387 ! 388 DO jk = 2, jpkm1 !* Interior point (w-masked 2nd order centered flux) 389 DO jj = 2, jpjm1 390 DO ji = fs_2, fs_jpim1 ! vector opt. 391 zwz(ji,jj,jk) = 0.5 * pwn(ji,jj,jk) * ( ptn(ji,jj,jk-1,jn) + ptn(ji,jj,jk,jn) ) * wmask(ji,jj,jk) 392 END DO 393 END DO 394 END DO 395 IF( ln_linssh ) THEN !* top value (only in linear free surf. as zwz is multiplied by wmask) 396 IF( ln_isfcav ) THEN ! ice-shelf cavities (top of the ocean) 397 DO jj = 1, jpj 398 DO ji = 1, jpi 399 zwz(ji,jj, mikt(ji,jj) ) = pwn(ji,jj,mikt(ji,jj)) * ptn(ji,jj,mikt(ji,jj),jn) ! linear free surface 400 END DO 401 END DO 402 ELSE ! no ocean cavities (only ocean surface) 403 zwz(:,:,1) = pwn(:,:,1) * ptn(:,:,1,jn) 404 ENDIF 402 405 ENDIF 403 406 ! 404 DO jk = 2, jpkm1 ! Interior point: second order centered tracer flux at w-point407 DO jk = 1, jpkm1 !== Tracer flux divergence added to the general trend ==! 405 408 DO jj = 2, jpjm1 406 409 DO ji = fs_2, fs_jpim1 ! vector opt. 407 zwz(ji,jj,jk) = 0.5 * pwn(ji,jj,jk) * ( ptn(ji,jj,jk-1,jn) + ptn(ji,jj,jk,jn) ) 408 END DO 409 END DO 410 END DO 411 ! 412 DO jk = 1, jpkm1 !== Tracer flux divergence added to the general trend ==! 413 DO jj = 2, jpjm1 414 DO ji = fs_2, fs_jpim1 ! vector opt. 415 zbtr = 1. / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) ) 416 ! k- vertical advective trends 417 ztra = - zbtr * ( zwz(ji,jj,jk) - zwz(ji,jj,jk+1) ) 418 ! added to the general tracer trends 419 pta(ji,jj,jk,jn) = pta(ji,jj,jk,jn) + ztra 420 END DO 421 END DO 422 END DO 423 ! ! Save the vertical advective trends for diagnostic 424 IF( l_trd ) CALL trd_tra( kt, cdtype, jn, jptra_trd_zad, zwz, pwn, ptn(:,:,:,jn) ) 410 pta(ji,jj,jk,jn) = pta(ji,jj,jk,jn) - ( zwz(ji,jj,jk) - zwz(ji,jj,jk+1) ) & 411 & * r1_e1e2t(ji,jj) / e3t_n(ji,jj,jk) 412 END DO 413 END DO 414 END DO 415 ! ! Send trends for diagnostic 416 IF( l_trd ) CALL trd_tra( kt, cdtype, jn, jptra_zad, zwz, pwn, ptn(:,:,:,jn) ) 425 417 ! 426 418 END DO 419 ! 420 CALL wrk_dealloc( jpi,jpj,jpk, zwz ) 427 421 ! 428 422 END SUBROUTINE tra_adv_cen2_k
Note: See TracChangeset
for help on using the changeset viewer.