Changeset 786 for branches/dev_001_GM/NEMO/OPA_SRC/TRA/traadv_ubs.F90
- Timestamp:
- 2008-01-10T18:11:23+01:00 (16 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/dev_001_GM/NEMO/OPA_SRC/TRA/traadv_ubs.F90
r719 r786 2 2 !!============================================================================== 3 3 !! *** MODULE traadv_ubs *** 4 !! Ocean activetracers: horizontal & vertical advective trend4 !! Ocean tracers: horizontal & vertical advective trend 5 5 !!============================================================================== 6 !! History : 9.0 ! 06-08 (L. Debreu, R. Benshila) Original code 6 !! History : 1.0 ! 06-08 (L. Debreu, R. Benshila) Original code 7 !! 2.4 ! 08-01 (G. Madec) Merge TRA-TRC 7 8 !!---------------------------------------------------------------------- 8 9 … … 11 12 !! advection trends using a third order biaised scheme 12 13 !!---------------------------------------------------------------------- 13 USE oce ! ocean dynamics and active tracers14 14 USE dom_oce ! ocean space and time domain 15 15 USE trdmod … … 33 33 # include "vectopt_loop_substitute.h90" 34 34 !!---------------------------------------------------------------------- 35 !! OPA 9.0 , LOCEAN-IPSL (2006)36 !! $ Header$35 !! NEMO/OPA & TRP 2.4 , LOCEAN-IPSL (2008) 36 !! $Id:$ 37 37 !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 38 38 !!---------------------------------------------------------------------- … … 40 40 CONTAINS 41 41 42 SUBROUTINE tra_adv_ubs( kt, pun, pvn, pwn ) 42 SUBROUTINE tra_adv_ubs( kt, cdtype, ktra, pun, pvn, pwn, & 43 & ptb, ptn, pta ) 43 44 !!---------------------------------------------------------------------- 44 45 !! *** ROUTINE tra_adv_ubs *** … … 70 71 !! 71 72 !! Reference : Shchepetkin, A. F., J. C. McWilliams, 2005, Ocean Modelling, 9, 347-404. 72 !! Farrow, D.E., Stevens, D.P., 1995, J. Phys. Ocean. 25, 1731Ð1741. 73 !!---------------------------------------------------------------------- 74 USE oce, ONLY : zwx => ua ! use ua as workspace 75 USE oce, ONLY : zwy => va ! use va as workspace 76 !! 77 INTEGER , INTENT(in) :: kt ! ocean time-step index 78 REAL(wp), INTENT(in), DIMENSION(jpi,jpj,jpk) :: pun ! effective ocean velocity, u_component 79 REAL(wp), INTENT(in), DIMENSION(jpi,jpj,jpk) :: pvn ! effective ocean velocity, v_component 80 REAL(wp), INTENT(in), DIMENSION(jpi,jpj,jpk) :: pwn ! effective ocean velocity, w_component 73 !! Farrow, D.E., Stevens, D.P., 1995, J. Phys. Ocean. 25, 1731–1741. 74 !!---------------------------------------------------------------------- 75 INTEGER , INTENT(in ) :: kt ! ocean time-step index 76 CHARACTER(len=3), INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator) 77 INTEGER , INTENT(in ) :: ktra ! tracer index 78 REAL(wp) , INTENT(in ), DIMENSION(jpi,jpj,jpk) :: pun, pvn, pwn ! 3 ocean velocity components 79 REAL(wp) , INTENT(inout), DIMENSION(jpi,jpj,jpk) :: ptb, ptn ! before and now tracer fields 80 REAL(wp) , INTENT(inout), DIMENSION(jpi,jpj,jpk) :: pta ! tracer trend 81 81 !! 82 82 INTEGER :: ji, jj, jk ! dummy loop indices 83 REAL(wp) :: zta, z sa, zbtr, zcoef ! temporary scalars84 REAL(wp) :: zfui, zfp_ui, zfm_ui, zcenut , zcenus! " "85 REAL(wp) :: zfvj, zfp_vj, zfm_vj, zcenvt , zcenvs! " "86 REAL(wp) :: z_hdivn _x, z_hdivn_y, z_hdivn! " "83 REAL(wp) :: zta, zbtr, zcoef ! temporary scalars 84 REAL(wp) :: zfui, zfp_ui, zfm_ui, zcenut ! " " 85 REAL(wp) :: zfvj, zfp_vj, zfm_vj, zcenvt ! " " 86 REAL(wp) :: z_hdivn ! " " 87 87 REAL(wp), DIMENSION(jpi,jpj) :: zeeu, zeev ! temporary 2D workspace 88 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zw z , zww! temporary 3D workspace88 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zwx, zwy ! temporary 3D workspace 89 89 REAL(wp), DIMENSION(jpi,jpj,jpk) :: ztu , ztv , zltu , zltv, ztrdt ! " " 90 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zsu , zsv , zlsu , zlsv, ztrds ! " "91 90 !!---------------------------------------------------------------------- 92 91 93 92 zltu(:,:,:) = 0.e0 94 93 zltv(:,:,:) = 0.e0 95 zlsu(:,:,:) = 0.e096 zlsv(:,:,:) = 0.e097 94 98 95 IF( kt == nit000 ) THEN … … 104 101 ENDIF 105 102 106 ! Save ta and sa trends 107 ztrdt(:,:,:) = ta(:,:,:) 108 ztrds(:,:,:) = sa(:,:,:) 103 ! store pta trends 104 ztrdt(:,:,:) = pta(:,:,:) 109 105 110 106 zcoef = 1./6. … … 132 128 DO jj = 1, jpjm1 133 129 DO ji = 1, fs_jpim1 ! vector opt. 134 ztu(ji,jj,jk) = zeeu(ji,jj) * ( tb(ji+1,jj ,jk) - tb(ji,jj,jk) ) 135 zsu(ji,jj,jk) = zeeu(ji,jj) * ( sb(ji+1,jj ,jk) - sb(ji,jj,jk) ) 136 ztv(ji,jj,jk) = zeev(ji,jj) * ( tb(ji ,jj+1,jk) - tb(ji,jj,jk) ) 137 zsv(ji,jj,jk) = zeev(ji,jj) * ( sb(ji ,jj+1,jk) - sb(ji,jj,jk) ) 130 ztu(ji,jj,jk) = zeeu(ji,jj) * ( ptb(ji+1,jj ,jk) - ptb(ji,jj,jk) ) 131 ztv(ji,jj,jk) = zeev(ji,jj) * ( ptb(ji ,jj+1,jk) - ptb(ji,jj,jk) ) 138 132 END DO 139 133 END DO … … 145 139 #endif 146 140 zltu(ji,jj,jk) = ( ztu(ji,jj,jk) - ztu(ji-1,jj,jk) ) * zcoef 147 zlsu(ji,jj,jk) = ( zsu(ji,jj,jk) - zsu(ji-1,jj,jk) ) * zcoef148 141 zltv(ji,jj,jk) = ( ztv(ji,jj,jk) - ztv(ji,jj-1,jk) ) * zcoef 149 zlsv(ji,jj,jk) = ( zsv(ji,jj,jk) - zsv(ji,jj-1,jk) ) * zcoef150 142 END DO 151 143 END DO … … 155 147 156 148 ! Lateral boundary conditions on the laplacian (zlt,zls) (unchanged sgn) 157 CALL lbc_lnk( zltu, 'T', 1. ) ; CALL lbc_lnk( zlsu, 'T', 1. ) 158 CALL lbc_lnk( zltv, 'T', 1. ) ; CALL lbc_lnk( zlsv, 'T', 1. ) 149 CALL lbc_lnk( zltu, 'T', 1. ) ; CALL lbc_lnk( zltv, 'T', 1. ) 159 150 160 151 ! ! =============== … … 178 169 zfm_vj = zfvj - ABS( zfvj ) 179 170 ! centered scheme 180 zcenut = zfui * ( tn(ji,jj,jk) + tn(ji+1,jj ,jk) ) 181 zcenvt = zfvj * ( tn(ji,jj,jk) + tn(ji ,jj+1,jk) ) 182 zcenus = zfui * ( sn(ji,jj,jk) + sn(ji+1,jj ,jk) ) 183 zcenvs = zfvj * ( sn(ji,jj,jk) + sn(ji ,jj+1,jk) ) 171 zcenut = zfui * ( ptn(ji,jj,jk) + ptn(ji+1,jj ,jk) ) 172 zcenvt = zfvj * ( ptn(ji,jj,jk) + ptn(ji ,jj+1,jk) ) 184 173 ! mixed centered / upstream scheme 185 174 zwx(ji,jj,jk) = zcenut - zfp_ui * zltu(ji,jj,jk) -zfm_ui * zltu(ji+1,jj,jk) 186 175 zwy(ji,jj,jk) = zcenvt - zfp_vj * zltv(ji,jj,jk) -zfm_vj * zltv(ji,jj+1,jk) 187 zww(ji,jj,jk) = zcenus - zfp_ui * zlsu(ji,jj,jk) -zfm_ui * zlsu(ji+1,jj,jk)188 zwz(ji,jj,jk) = zcenvs - zfp_vj * zlsv(ji,jj,jk) -zfm_vj * zlsv(ji,jj+1,jk)189 176 END DO 190 177 END DO … … 201 188 zta = - zbtr * ( zwx(ji,jj,jk) - zwx(ji-1,jj ,jk) & 202 189 & + zwy(ji,jj,jk) - zwy(ji ,jj-1,jk) ) 203 zsa = - zbtr * ( zww(ji,jj,jk) - zww(ji-1,jj ,jk) &204 & + zwz(ji,jj,jk) - zwz(ji ,jj-1,jk) )205 190 ! add it to the general tracer trends 206 ta(ji,jj,jk) = ta(ji,jj,jk) + zta 207 sa(ji,jj,jk) = sa(ji,jj,jk) + zsa 191 pta(ji,jj,jk) = pta(ji,jj,jk) + zta 208 192 END DO 209 193 END DO … … 213 197 214 198 ! Horizontal trend used in tra_adv_ztvd subroutine 215 zltu(:,:,:) = ta(:,:,:) - ztrdt(:,:,:) 216 zlsu(:,:,:) = sa(:,:,:) - ztrds(:,:,:) 217 218 ! 3. Save the horizontal advective trends for diagnostic 219 ! ------------------------------------------------------ 220 IF( l_trdtra ) THEN 221 ! Recompute the hoizontal advection zta & zsa trends computed 222 ! at the step 2. above in making the difference between the new 223 ! trends and the previous one ta()/sa - ztrdt()/ztrds() and add 224 ! the term tn()/sn()*hdivn() to recover the Uh gradh(T/S) trends 225 ztrdt(:,:,:) = 0.e0 ; ztrds(:,:,:) = 0.e0 226 ! 227 ! T/S ZONAL advection trends 228 DO jk = 1, jpkm1 229 DO jj = 2, jpjm1 230 DO ji = fs_2, fs_jpim1 ! vector opt. 231 !-- Compute zonal divergence by splitting hdivn (see divcur.F90) 232 #if defined key_zco 233 zbtr = e1e2tr(ji,jj) 234 z_hdivn_x = ( e2u(ji ,jj) * pun(ji ,jj,jk) & 235 & - e2u(ji-1,jj) * pun(ji-1,jj,jk) ) * zbtr 236 #else 237 zbtr = e1e2tr(ji,jj) / fse3t(ji,jj,jk) 238 z_hdivn_x = ( e2u(ji ,jj) * fse3u(ji ,jj,jk) * pun(ji ,jj,jk) & 239 & - e2u(ji-1,jj) * fse3u(ji-1,jj,jk) * pun(ji-1,jj,jk) ) * zbtr 240 #endif 241 ztrdt(ji,jj,jk) = - ( zwx(ji,jj,jk) - zwx(ji-1,jj,jk) ) * zbtr + tn(ji,jj,jk) * z_hdivn_x 242 ztrds(ji,jj,jk) = - ( zww(ji,jj,jk) - zww(ji-1,jj,jk) ) * zbtr + sn(ji,jj,jk) * z_hdivn_x 243 END DO 244 END DO 245 END DO 246 CALL trd_mod(ztrdt, ztrds, jptra_trd_xad, 'TRA', kt) ! save the trends 247 ! 248 ! T/S MERIDIONAL advection trends 249 DO jk = 1, jpkm1 250 DO jj = 2, jpjm1 251 DO ji = fs_2, fs_jpim1 ! vector opt. 252 !-- Compute merid. divergence by splitting hdivn (see divcur.F90) 253 #if defined key_zco 254 zbtr = e1e2tr(ji,jj) 255 z_hdivn_y = ( e1v(ji, jj) * pvn(ji,jj ,jk) & 256 & - e1v(ji,jj-1) * pvn(ji,jj-1,jk) ) * zbtr 257 #else 258 zbtr = e1e2tr(ji,jj) / fse3t(ji,jj,jk) 259 z_hdivn_y = ( e1v(ji, jj) * fse3v(ji,jj ,jk) * pvn(ji,jj ,jk) & 260 & - e1v(ji,jj-1) * fse3v(ji,jj-1,jk) * pvn(ji,jj-1,jk) ) * zbtr 261 #endif 262 ztrdt(ji,jj,jk) = - ( zwy(ji,jj,jk) - zwy(ji,jj-1,jk) ) * zbtr + tn(ji,jj,jk) * z_hdivn_y 263 ztrds(ji,jj,jk) = - ( zwz(ji,jj,jk) - zwz(ji,jj-1,jk) ) * zbtr + sn(ji,jj,jk) * z_hdivn_y 264 END DO 265 END DO 266 END DO 267 CALL trd_mod(ztrdt, ztrds, jptra_trd_yad, 'TRA', kt) ! save the trends 268 ! 269 ENDIF 270 271 IF(ln_ctl) CALL prt_ctl( tab3d_1=ta, clinfo1=' ubs had - Ta: ', mask1=tmask, & 272 & tab3d_2=sa, clinfo2= ' Sa: ', mask2=tmask, clinfo3='tra' ) 273 274 ! "zonal" mean advective heat and salt transport 275 IF( ln_diaptr .AND. ( MOD( kt, nf_ptr ) == 0 ) ) THEN 199 zltu(:,:,:) = pta(:,:,:) - ztrdt(:,:,:) 200 201 ! Save the horizontal advective trends for diagnostic 202 ! ----------------------------------------------------- 203 IF( l_trdtra ) THEN 204 CALL trd_tra_adv( kt, ktra, jpt_trd_xad, cdtype, zwx, pun, ptn ) 205 CALL trd_tra_adv( kt, ktra, jpt_trd_yad, cdtype, zwy, pvn, ptn ) 206 ENDIF 207 208 ! "Poleward" heat or salt transport 209 IF( cdtype == 'TRA' .AND. ln_diaptr .AND. ( MOD( kt, nf_ptr ) == 0 ) ) THEN 276 210 IF( lk_zco ) THEN 277 211 DO jk = 1, jpkm1 278 212 DO jj = 2, jpjm1 279 213 DO ji = fs_2, fs_jpim1 ! vector opt. 280 zwy(ji,jj,jk) = zwy(ji,jj,jk) * fse3v(ji,jj,jk) 281 zwz(ji,jj,jk) = zwz(ji,jj,jk) * fse3v(ji,jj,jk) 214 zwy(ji,jj,jk) = zwy(ji,jj,jk) * fse3v(ji,jj,jk) 282 215 END DO 283 216 END DO 284 217 END DO 285 218 ENDIF 286 pht_adv(:) = ptr_vj( zwy(:,:,:) ) 287 pst_adv(:) = ptr_vj( zwz(:,:,:) ) 288 ENDIF 219 IF( ktra == jp_tem) pht_adv(:) = ptr_vj( zwy(:,:,:) ) 220 IF( ktra == jp_sal) pst_adv(:) = ptr_vj( zwy(:,:,:) ) 221 ENDIF 222 223 IF(ln_ctl) CALL prt_ctl( tab3d_1=pta, clinfo1=' ubs - had: ', mask1=tmask, clinfo3=cdtype ) 224 289 225 290 226 ! II. Vertical advection 291 227 ! ---------------------- 292 IF( l_trdtra ) THEN ! Save ta and sa trends 293 ztrdt(:,:,:) = ta(:,:,:) 294 ztrds(:,:,:) = sa(:,:,:) 295 ENDIF 228 IF( l_trdtra ) ztrdt(:,:,:) = pta(:,:,:) ! Save ta and sa trends 296 229 297 230 ! TVD scheme the vertical direction 298 CALL tra_adv_ztvd( kt, pwn, zltu, zlsu)299 300 IF( l_trdtra ) THEN ! Save the final vertical advective trends231 CALL tra_adv_ztvd( kt, pwn, zltu, ptb, ptn, pta ) 232 233 IF( l_trdtra ) THEN ! vertical advective trend diagnostics 301 234 DO jk = 1, jpkm1 302 235 DO jj = 2, jpjm1 303 236 DO ji = fs_2, fs_jpim1 ! vector opt. 304 #if defined key_zco 305 zbtr = e1e2tr(ji,jj) 306 z_hdivn_x = e2u(ji,jj)*pun(ji,jj,jk) - e2u(ji-1,jj)*pun(ji-1,jj,jk) 307 z_hdivn_y = e1v(ji,jj)*pvn(ji,jj,jk) - e1v(ji,jj-1)*pvn(ji,jj-1,jk) 308 #else 309 zbtr = e1e2tr(ji,jj) / fse3t(ji,jj,jk) 310 z_hdivn_x = e2u(ji,jj)*fse3u(ji,jj,jk)*pun(ji,jj,jk) - e2u(ji-1,jj)*fse3u(ji-1,jj,jk)*pun(ji-1,jj,jk) 311 z_hdivn_y = e1v(ji,jj)*fse3v(ji,jj,jk)*pvn(ji,jj,jk) - e1v(ji,jj-1)*fse3v(ji,jj-1,jk)*pvn(ji,jj-1,jk) 312 #endif 313 z_hdivn = (z_hdivn_x + z_hdivn_y) * zbtr 314 zbtr = e1e2tr(ji,jj) / fse3t(ji,jj,jk) 315 ztrdt(ji,jj,jk) = ta(ji,jj,jk) - ztrdt(ji,jj,jk) - tn(ji,jj,jk) * z_hdivn 316 ztrds(ji,jj,jk) = sa(ji,jj,jk) - ztrds(ji,jj,jk) - sn(ji,jj,jk) * z_hdivn 237 z_hdivn = ( pwn(ji,jj,jk) - pwn(ji,jj,jk+1) ) / fse3t(ji,jj,jk) 238 ztrdt(ji,jj,jk) = pta(ji,jj,jk) - ztrdt(ji,jj,jk+1) + ptn(ji,jj,jk) * z_hdivn 317 239 END DO 318 240 END DO 319 241 END DO 320 CALL trd_ mod(ztrdt, ztrds, jptra_trd_zad, 'TRA', kt) ! <<< ADD TO PREVIOUSLY COMPUTED242 CALL trd_tra( kt, ktra, jpt_trd_zad, cdtype, ptrd3d=ztrdt ) 321 243 ! 322 244 ENDIF 323 324 IF(ln_ctl) CALL prt_ctl( tab3d_1=ta, clinfo1=' ubs zad - Ta: ', mask1=tmask, & 325 & tab3d_2=sa, clinfo2= ' Sa: ', mask2=tmask, clinfo3='tra') 245 246 IF(ln_ctl) CALL prt_ctl( tab3d_1=pta, clinfo1=' ubs - zad: ', mask1=tmask, clinfo3=cdtype ) 326 247 ! 327 248 END SUBROUTINE tra_adv_ubs 328 249 329 250 330 SUBROUTINE tra_adv_ztvd( kt, pwn, zttrd, zstrd)251 SUBROUTINE tra_adv_ztvd( kt, pwn, zttrd, ptb, ptn, pta ) 331 252 !!---------------------------------------------------------------------- 332 253 !! *** ROUTINE tra_adv_ztvd *** … … 342 263 !! - save the trends in (ztrdt,ztrds) ('key_trdtra') 343 264 !!---------------------------------------------------------------------- 344 INTEGER , INTENT(in) :: kt ! ocean time-step 345 REAL(wp), INTENT(in), DIMENSION(jpi,jpj,jpk) :: pwn ! verical effective velocity 346 REAL(wp), INTENT(in), DIMENSION(jpi,jpj,jpk) :: zttrd, zstrd ! lateral advective trends on T & S 265 INTEGER , INTENT(in ) :: kt ! ocean time-step 266 REAL(wp), INTENT(in ), DIMENSION(jpi,jpj,jpk) :: pwn ! verical effective velocity 267 REAL(wp), INTENT(in ), DIMENSION(jpi,jpj,jpk) :: zttrd ! lateral advective trends on T & S 268 REAL(wp), INTENT(inout), DIMENSION(jpi,jpj,jpk) :: ptb, ptn ! before and now tracer fields 269 REAL(wp), INTENT(inout), DIMENSION(jpi,jpj,jpk) :: pta ! tracer trend 347 270 !! 348 271 INTEGER :: ji, jj, jk ! dummy loop indices 349 272 REAL(wp) :: z2dtt, zbtr, zew, z2 ! temporary scalar 350 REAL(wp) :: ztak, zfp_wk ! " " 351 REAL(wp) :: zsak, zfm_wk ! " " 273 REAL(wp) :: ztak, zfp_wk, zfm_wk ! " " 352 274 REAL(wp), DIMENSION (jpi,jpj,jpk) :: zti, ztw ! temporary 3D workspace 353 REAL(wp), DIMENSION (jpi,jpj,jpk) :: zsi, zsw ! " "354 275 !!---------------------------------------------------------------------- 355 276 … … 366 287 ! Bottom value : flux set to zero 367 288 ! -------------- 368 ztw(:,:,jpk) = 0.e0 ; zsw(:,:,jpk) = 0.e0 369 zti (:,:,:) = 0.e0 ; zsi (:,:,:) = 0.e0 289 ztw(:,:,jpk) = 0.e0 ; zti (:,:,:) = 0.e0 370 290 371 291 … … 375 295 IF( lk_dynspg_rl .OR. lk_vvl ) THEN ! rigid lid : flux set to zero 376 296 ztw(:,:,1) = 0.e0 377 zsw(:,:,1) = 0.e0 378 ELSE ! free surface 379 DO jj = 1, jpj 380 DO ji = 1, jpi 381 zew = e1t(ji,jj) * e2t(ji,jj) * pwn(ji,jj,1) 382 ztw(ji,jj,1) = zew * tb(ji,jj,1) 383 zsw(ji,jj,1) = zew * sb(ji,jj,1) 384 END DO 385 END DO 297 ELSE ! free surface 298 ztw(:,:,1) = e1t(:,:) * e2t(:,:) * pwn(:,:,1) * ptb(:,:,1) 386 299 ENDIF 387 300 … … 393 306 zfp_wk = zew + ABS( zew ) 394 307 zfm_wk = zew - ABS( zew ) 395 ztw(ji,jj,jk) = zfp_wk * tb(ji,jj,jk) + zfm_wk * tb(ji,jj,jk-1) 396 zsw(ji,jj,jk) = zfp_wk * sb(ji,jj,jk) + zfm_wk * sb(ji,jj,jk-1) 308 ztw(ji,jj,jk) = zfp_wk * ptb(ji,jj,jk) + zfm_wk * ptb(ji,jj,jk-1) 397 309 END DO 398 310 END DO … … 406 318 zbtr = 1./ ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) ) 407 319 ztak = - ( ztw(ji,jj,jk) - ztw(ji,jj,jk+1) ) * zbtr 408 zsak = - ( zsw(ji,jj,jk) - zsw(ji,jj,jk+1) ) * zbtr 409 ta(ji,jj,jk) = ta(ji,jj,jk) + ztak 410 sa(ji,jj,jk) = sa(ji,jj,jk) + zsak 411 zti (ji,jj,jk) = ( tb(ji,jj,jk) + z2dtt * ( ztak + zttrd(ji,jj,jk) ) ) * tmask(ji,jj,jk) 412 zsi (ji,jj,jk) = ( sb(ji,jj,jk) + z2dtt * ( zsak + zstrd(ji,jj,jk) ) ) * tmask(ji,jj,jk) 320 pta(ji,jj,jk) = pta(ji,jj,jk) + ztak 321 zti (ji,jj,jk) = ( ptb(ji,jj,jk) + z2dtt * ( ztak + zttrd(ji,jj,jk) ) ) * tmask(ji,jj,jk) 413 322 END DO 414 323 END DO … … 417 326 ! Lateral boundary conditions on zti, zsi (unchanged sign) 418 327 CALL lbc_lnk( zti, 'T', 1. ) 419 CALL lbc_lnk( zsi, 'T', 1. )420 328 421 329 422 330 ! antidiffusive flux : high order minus low order 423 331 ! ------------------------------------------------- 424 ! Surface value 425 ztw(:,:,1) = 0.e0 ; zsw(:,:,1) = 0.e0 426 427 ! Interior value 428 DO jk = 2, jpkm1 332 ztw(:,:,1) = 0.e0 ! Surface value 333 334 DO jk = 2, jpkm1 ! Interior value 429 335 DO jj = 1, jpj 430 336 DO ji = 1, jpi 431 337 zew = 0.5 * e1t(ji,jj) * e2t(ji,jj) * pwn(ji,jj,jk) 432 ztw(ji,jj,jk) = zew * ( tn(ji,jj,jk) + tn(ji,jj,jk-1) ) - ztw(ji,jj,jk) 433 zsw(ji,jj,jk) = zew * ( sn(ji,jj,jk) + sn(ji,jj,jk-1) ) - zsw(ji,jj,jk) 338 ztw(ji,jj,jk) = zew * ( ptn(ji,jj,jk) + ptn(ji,jj,jk-1) ) - ztw(ji,jj,jk) 434 339 END DO 435 340 END DO … … 438 343 ! monotonicity algorithm 439 344 ! ------------------------ 440 CALL nonosc_z( tb, ztw, zti, z2 ) 441 CALL nonosc_z( sb, zsw, zsi, z2 ) 345 CALL nonosc_z( ptb, ztw, zti, z2 ) 442 346 443 347 … … 450 354 ! k- vertical advective trends 451 355 ztak = - ( ztw(ji,jj,jk) - ztw(ji,jj,jk+1) ) * zbtr 452 zsak = - ( zsw(ji,jj,jk) - zsw(ji,jj,jk+1) ) * zbtr453 356 ! add them to the general tracer trends 454 ta(ji,jj,jk) = ta(ji,jj,jk) + ztak 455 sa(ji,jj,jk) = sa(ji,jj,jk) + zsak 357 pta(ji,jj,jk) = pta(ji,jj,jk) + ztak 456 358 END DO 457 359 END DO
Note: See TracChangeset
for help on using the changeset viewer.