Changeset 10874 for NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traadv_fct.F90
- Timestamp:
- 2019-04-15T15:57:37+02:00 (5 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traadv_fct.F90
r10806 r10874 52 52 CONTAINS 53 53 54 SUBROUTINE tra_adv_fct( kt, kit000, ktlev1, ktlev2, ktlev3, cdtype, p2dt, pu, pv, pw, &55 & pt_lev1, pt_lev2, pt_rhs, kjpt, kn_fct_h, kn_fct_v )54 SUBROUTINE tra_adv_fct( kt, kit000, cdtype, p2dt, pun, pvn, pwn, & 55 & ptb, ptn, pta, kjpt, kn_fct_h, kn_fct_v ) 56 56 !!---------------------------------------------------------------------- 57 57 !! *** ROUTINE tra_adv_fct *** … … 65 65 !! - corrected flux (monotonic correction) 66 66 !! 67 !! ** Action : - update pt _rhswith the now advective tracer trends67 !! ** Action : - update pta with the now advective tracer trends 68 68 !! - send trends to trdtra module for further diagnostics (l_trdtra=T) 69 69 !! - htr_adv, str_adv : poleward advective heat and salt transport (ln_diaptr=T) 70 70 !!---------------------------------------------------------------------- 71 71 INTEGER , INTENT(in ) :: kt ! ocean time-step index 72 INTEGER , INTENT(in ) :: ktlev1, ktlev2, ktlev3 ! time level indices for source terms73 72 INTEGER , INTENT(in ) :: kit000 ! first time step index 74 73 CHARACTER(len=3) , INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator) … … 77 76 INTEGER , INTENT(in ) :: kn_fct_v ! order of the FCT scheme (=2 or 4) 78 77 REAL(wp) , INTENT(in ) :: p2dt ! tracer time-step 79 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ) :: pu , pv, pw! 3 ocean velocity components80 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: pt _lev1, pt_lev2! before and now tracer fields81 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pt _rhs! tracer trend78 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ) :: pun, pvn, pwn ! 3 ocean velocity components 79 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: ptb, ptn ! before and now tracer fields 80 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pta ! tracer trend 82 81 ! 83 82 INTEGER :: ji, jj, jk, jn ! dummy loop indices … … 126 125 DO ji = 1, fs_jpim1 ! vector opt. 127 126 ! upstream scheme 128 zfp_ui = pu (ji,jj,jk) + ABS( pu(ji,jj,jk) )129 zfm_ui = pu (ji,jj,jk) - ABS( pu(ji,jj,jk) )130 zfp_vj = pv (ji,jj,jk) + ABS( pv(ji,jj,jk) )131 zfm_vj = pv (ji,jj,jk) - ABS( pv(ji,jj,jk) )132 zwx(ji,jj,jk) = 0.5 * ( zfp_ui * pt _lev1(ji,jj,jk,jn) + zfm_ui * pt_lev1(ji+1,jj ,jk,jn) )133 zwy(ji,jj,jk) = 0.5 * ( zfp_vj * pt _lev1(ji,jj,jk,jn) + zfm_vj * pt_lev1(ji ,jj+1,jk,jn) )127 zfp_ui = pun(ji,jj,jk) + ABS( pun(ji,jj,jk) ) 128 zfm_ui = pun(ji,jj,jk) - ABS( pun(ji,jj,jk) ) 129 zfp_vj = pvn(ji,jj,jk) + ABS( pvn(ji,jj,jk) ) 130 zfm_vj = pvn(ji,jj,jk) - ABS( pvn(ji,jj,jk) ) 131 zwx(ji,jj,jk) = 0.5 * ( zfp_ui * ptb(ji,jj,jk,jn) + zfm_ui * ptb(ji+1,jj ,jk,jn) ) 132 zwy(ji,jj,jk) = 0.5 * ( zfp_vj * ptb(ji,jj,jk,jn) + zfm_vj * ptb(ji ,jj+1,jk,jn) ) 134 133 END DO 135 134 END DO … … 139 138 DO jj = 1, jpj 140 139 DO ji = 1, jpi 141 zfp_wk = pw (ji,jj,jk) + ABS( pw(ji,jj,jk) )142 zfm_wk = pw (ji,jj,jk) - ABS( pw(ji,jj,jk) )143 zwz(ji,jj,jk) = 0.5 * ( zfp_wk * pt _lev1(ji,jj,jk,jn) + zfm_wk * pt_lev1(ji,jj,jk-1,jn) ) * wmask(ji,jj,jk)140 zfp_wk = pwn(ji,jj,jk) + ABS( pwn(ji,jj,jk) ) 141 zfm_wk = pwn(ji,jj,jk) - ABS( pwn(ji,jj,jk) ) 142 zwz(ji,jj,jk) = 0.5 * ( zfp_wk * ptb(ji,jj,jk,jn) + zfm_wk * ptb(ji,jj,jk-1,jn) ) * wmask(ji,jj,jk) 144 143 END DO 145 144 END DO … … 149 148 DO jj = 1, jpj 150 149 DO ji = 1, jpi 151 zwz(ji,jj, mikt(ji,jj) ) = pw (ji,jj,mikt(ji,jj)) * pt_lev1(ji,jj,mikt(ji,jj),jn) ! linear free surface150 zwz(ji,jj, mikt(ji,jj) ) = pwn(ji,jj,mikt(ji,jj)) * ptb(ji,jj,mikt(ji,jj),jn) ! linear free surface 152 151 END DO 153 152 END DO 154 153 ELSE ! no cavities: only at the ocean surface 155 zwz(:,:,1) = pw (:,:,1) * pt_lev1(:,:,1,jn)154 zwz(:,:,1) = pwn(:,:,1) * ptb(:,:,1,jn) 156 155 ENDIF 157 156 ENDIF … … 165 164 & + zwz(ji,jj,jk) - zwz(ji ,jj ,jk+1) ) * r1_e1e2t(ji,jj) 166 165 ! ! update and guess with monotonic sheme 167 pt _rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn) + ztra / e3t(ji,jj,jk,ktlev2) * tmask(ji,jj,jk)168 zwi(ji,jj,jk) = ( e3t (ji,jj,jk,ktlev1) * pt_lev1(ji,jj,jk,jn) + p2dt * ztra ) / e3t(ji,jj,jk,ktlev3) * tmask(ji,jj,jk)166 pta(ji,jj,jk,jn) = pta(ji,jj,jk,jn) + ztra / e3t_n(ji,jj,jk) * tmask(ji,jj,jk) 167 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) 169 168 END DO 170 169 END DO … … 185 184 DO jj = 1, jpjm1 186 185 DO ji = 1, fs_jpim1 ! vector opt. 187 zwx(ji,jj,jk) = 0.5_wp * pu (ji,jj,jk) * ( pt_lev2(ji,jj,jk,jn) + pt_lev2(ji+1,jj,jk,jn) ) - zwx(ji,jj,jk)188 zwy(ji,jj,jk) = 0.5_wp * pv (ji,jj,jk) * ( pt_lev2(ji,jj,jk,jn) + pt_lev2(ji,jj+1,jk,jn) ) - zwy(ji,jj,jk)186 zwx(ji,jj,jk) = 0.5_wp * pun(ji,jj,jk) * ( ptn(ji,jj,jk,jn) + ptn(ji+1,jj,jk,jn) ) - zwx(ji,jj,jk) 187 zwy(ji,jj,jk) = 0.5_wp * pvn(ji,jj,jk) * ( ptn(ji,jj,jk,jn) + ptn(ji,jj+1,jk,jn) ) - zwy(ji,jj,jk) 189 188 END DO 190 189 END DO … … 197 196 DO jj = 1, jpjm1 ! 1st derivative (gradient) 198 197 DO ji = 1, fs_jpim1 ! vector opt. 199 ztu(ji,jj,jk) = ( pt _lev2(ji+1,jj ,jk,jn) - pt_lev2(ji,jj,jk,jn) ) * umask(ji,jj,jk)200 ztv(ji,jj,jk) = ( pt _lev2(ji ,jj+1,jk,jn) - pt_lev2(ji,jj,jk,jn) ) * vmask(ji,jj,jk)198 ztu(ji,jj,jk) = ( ptn(ji+1,jj ,jk,jn) - ptn(ji,jj,jk,jn) ) * umask(ji,jj,jk) 199 ztv(ji,jj,jk) = ( ptn(ji ,jj+1,jk,jn) - ptn(ji,jj,jk,jn) ) * vmask(ji,jj,jk) 201 200 END DO 202 201 END DO … … 213 212 DO jj = 1, jpjm1 214 213 DO ji = 1, fs_jpim1 ! vector opt. 215 zC2t_u = pt _lev2(ji,jj,jk,jn) + pt_lev2(ji+1,jj ,jk,jn) ! 2 x C2 interpolation of T at u- & v-points216 zC2t_v = pt _lev2(ji,jj,jk,jn) + pt_lev2(ji ,jj+1,jk,jn)214 zC2t_u = ptn(ji,jj,jk,jn) + ptn(ji+1,jj ,jk,jn) ! 2 x C2 interpolation of T at u- & v-points 215 zC2t_v = ptn(ji,jj,jk,jn) + ptn(ji ,jj+1,jk,jn) 217 216 ! ! C4 minus upstream advective fluxes 218 zwx(ji,jj,jk) = 0.5_wp * pu (ji,jj,jk) * ( zC2t_u + zltu(ji,jj,jk) - zltu(ji+1,jj,jk) ) - zwx(ji,jj,jk)219 zwy(ji,jj,jk) = 0.5_wp * pv (ji,jj,jk) * ( zC2t_v + zltv(ji,jj,jk) - zltv(ji,jj+1,jk) ) - zwy(ji,jj,jk)217 zwx(ji,jj,jk) = 0.5_wp * pun(ji,jj,jk) * ( zC2t_u + zltu(ji,jj,jk) - zltu(ji+1,jj,jk) ) - zwx(ji,jj,jk) 218 zwy(ji,jj,jk) = 0.5_wp * pvn(ji,jj,jk) * ( zC2t_v + zltv(ji,jj,jk) - zltv(ji,jj+1,jk) ) - zwy(ji,jj,jk) 220 219 END DO 221 220 END DO … … 228 227 DO jj = 1, jpjm1 229 228 DO ji = 1, fs_jpim1 ! vector opt. 230 ztu(ji,jj,jk) = ( pt _lev2(ji+1,jj ,jk,jn) - pt_lev2(ji,jj,jk,jn) ) * umask(ji,jj,jk)231 ztv(ji,jj,jk) = ( pt _lev2(ji ,jj+1,jk,jn) - pt_lev2(ji,jj,jk,jn) ) * vmask(ji,jj,jk)229 ztu(ji,jj,jk) = ( ptn(ji+1,jj ,jk,jn) - ptn(ji,jj,jk,jn) ) * umask(ji,jj,jk) 230 ztv(ji,jj,jk) = ( ptn(ji ,jj+1,jk,jn) - ptn(ji,jj,jk,jn) ) * vmask(ji,jj,jk) 232 231 END DO 233 232 END DO … … 238 237 DO jj = 2, jpjm1 239 238 DO ji = 2, fs_jpim1 ! vector opt. 240 zC2t_u = pt _lev2(ji,jj,jk,jn) + pt_lev2(ji+1,jj ,jk,jn) ! 2 x C2 interpolation of T at u- & v-points (x2)241 zC2t_v = pt _lev2(ji,jj,jk,jn) + pt_lev2(ji ,jj+1,jk,jn)239 zC2t_u = ptn(ji,jj,jk,jn) + ptn(ji+1,jj ,jk,jn) ! 2 x C2 interpolation of T at u- & v-points (x2) 240 zC2t_v = ptn(ji,jj,jk,jn) + ptn(ji ,jj+1,jk,jn) 242 241 ! ! C4 interpolation of T at u- & v-points (x2) 243 242 zC4t_u = zC2t_u + r1_6 * ( ztu(ji-1,jj ,jk) - ztu(ji+1,jj ,jk) ) 244 243 zC4t_v = zC2t_v + r1_6 * ( ztv(ji ,jj-1,jk) - ztv(ji ,jj+1,jk) ) 245 244 ! ! C4 minus upstream advective fluxes 246 zwx(ji,jj,jk) = 0.5_wp * pu (ji,jj,jk) * zC4t_u - zwx(ji,jj,jk)247 zwy(ji,jj,jk) = 0.5_wp * pv (ji,jj,jk) * zC4t_v - zwy(ji,jj,jk)245 zwx(ji,jj,jk) = 0.5_wp * pun(ji,jj,jk) * zC4t_u - zwx(ji,jj,jk) 246 zwy(ji,jj,jk) = 0.5_wp * pvn(ji,jj,jk) * zC4t_v - zwy(ji,jj,jk) 248 247 END DO 249 248 END DO … … 258 257 DO jj = 2, jpjm1 259 258 DO ji = fs_2, fs_jpim1 260 zwz(ji,jj,jk) = ( pw (ji,jj,jk) * 0.5_wp * ( pt_lev2(ji,jj,jk,jn) + pt_lev2(ji,jj,jk-1,jn) ) &259 zwz(ji,jj,jk) = ( pwn(ji,jj,jk) * 0.5_wp * ( ptn(ji,jj,jk,jn) + ptn(ji,jj,jk-1,jn) ) & 261 260 & - zwz(ji,jj,jk) ) * wmask(ji,jj,jk) 262 261 END DO … … 265 264 ! 266 265 CASE( 4 ) !- 4th order COMPACT 267 CALL interp_4th_cpt( pt _lev2(:,:,:,jn) , ztw ) ! zwt = COMPACT interpolation of T at w-point266 CALL interp_4th_cpt( ptn(:,:,:,jn) , ztw ) ! zwt = COMPACT interpolation of T at w-point 268 267 DO jk = 2, jpkm1 269 268 DO jj = 2, jpjm1 270 269 DO ji = fs_2, fs_jpim1 271 zwz(ji,jj,jk) = ( pw (ji,jj,jk) * ztw(ji,jj,jk) - zwz(ji,jj,jk) ) * wmask(ji,jj,jk)270 zwz(ji,jj,jk) = ( pwn(ji,jj,jk) * ztw(ji,jj,jk) - zwz(ji,jj,jk) ) * wmask(ji,jj,jk) 272 271 END DO 273 272 END DO … … 283 282 ! !== monotonicity algorithm ==! 284 283 ! 285 CALL nonosc( pt _lev1(:,:,:,jn), zwx, zwy, zwz, zwi, p2dt, e3t(:,:,:,ktlev2))284 CALL nonosc( ptb(:,:,:,jn), zwx, zwy, zwz, zwi, p2dt ) 286 285 ! 287 286 ! !== final trend with corrected fluxes ==! … … 290 289 DO jj = 2, jpjm1 291 290 DO ji = fs_2, fs_jpim1 ! vector opt. 292 pt _rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn) - ( zwx(ji,jj,jk) - zwx(ji-1,jj ,jk ) &291 pta(ji,jj,jk,jn) = pta(ji,jj,jk,jn) - ( zwx(ji,jj,jk) - zwx(ji-1,jj ,jk ) & 293 292 & + zwy(ji,jj,jk) - zwy(ji ,jj-1,jk ) & 294 293 & + zwz(ji,jj,jk) - zwz(ji ,jj ,jk+1) ) & 295 & * r1_e1e2t(ji,jj) / e3t (ji,jj,jk,ktlev2)294 & * r1_e1e2t(ji,jj) / e3t_n(ji,jj,jk) 296 295 END DO 297 296 END DO … … 304 303 ! 305 304 IF( l_trd ) THEN ! trend diagnostics 306 CALL trd_tra( kt, cdtype, jn, jptra_xad, ztrdx, pu , pt_lev2(:,:,:,jn) )307 CALL trd_tra( kt, cdtype, jn, jptra_yad, ztrdy, pv , pt_lev2(:,:,:,jn) )308 CALL trd_tra( kt, cdtype, jn, jptra_zad, ztrdz, pw , pt_lev2(:,:,:,jn) )305 CALL trd_tra( kt, cdtype, jn, jptra_xad, ztrdx, pun, ptn(:,:,:,jn) ) 306 CALL trd_tra( kt, cdtype, jn, jptra_yad, ztrdy, pvn, ptn(:,:,:,jn) ) 307 CALL trd_tra( kt, cdtype, jn, jptra_zad, ztrdz, pwn, ptn(:,:,:,jn) ) 309 308 ENDIF 310 309 ! ! heat/salt transport … … 329 328 330 329 331 SUBROUTINE nonosc( pbef, paa, pbb, pcc, paft, p2dt , pe3t)330 SUBROUTINE nonosc( pbef, paa, pbb, pcc, paft, p2dt ) 332 331 !!--------------------------------------------------------------------- 333 332 !! *** ROUTINE nonosc *** … … 342 341 !! in-space based differencing for fluid 343 342 !!---------------------------------------------------------------------- 344 REAL(wp) , INTENT(in ) :: p2dt 345 REAL(wp), DIMENSION (jpi,jpj,jpk), INTENT(in ) :: pbef, paft , pe3t ! before & after field, now e3tfield346 REAL(wp), DIMENSION (jpi,jpj,jpk), INTENT(inout) :: paa, pbb, pcc 343 REAL(wp) , INTENT(in ) :: p2dt ! tracer time-step 344 REAL(wp), DIMENSION (jpi,jpj,jpk), INTENT(in ) :: pbef, paft ! before & after field 345 REAL(wp), DIMENSION (jpi,jpj,jpk), INTENT(inout) :: paa, pbb, pcc ! monotonic fluxes in the 3 directions 347 346 ! 348 347 INTEGER :: ji, jj, jk ! dummy loop indices … … 393 392 394 393 ! up & down beta terms 395 zbt = e1e2t(ji,jj) * pe3t(ji,jj,jk) / p2dt394 zbt = e1e2t(ji,jj) * e3t_n(ji,jj,jk) / p2dt 396 395 zbetup(ji,jj,jk) = ( zup - paft(ji,jj,jk) ) / ( zpos + zrtrn ) * zbt 397 396 zbetdo(ji,jj,jk) = ( paft(ji,jj,jk) - zdo ) / ( zneg + zrtrn ) * zbt … … 635 634 !! The tri-diagonals matrix is given as input 3D arrays: pD, pU, pL 636 635 !! (i.e. the Diagonal, the Upper diagonal, and the Lower diagonal). 637 !! The solution is pt _rhs.636 !! The solution is pta. 638 637 !! The 3d array zwt is used as a work space array. 639 638 !!----------------------------------------------------------------------
Note: See TracChangeset
for help on using the changeset viewer.
