Changeset 786 for branches/dev_001_GM/NEMO/OPA_SRC/TRA/traadv_tvd.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_tvd.F90
r719 r786 4 4 !! Ocean active tracers: horizontal & vertical advective trend 5 5 !!============================================================================== 6 !! History : 7 !! 8 !! 9 !! 10 !! 6 !! History : 7.0 ! 95-12 (L. Mortier) Original code 7 !! 8.0 ! 00-01 (H. Loukos) adapted to ORCA 8 !! - ! 00-10 (MA Foujols E.Kestenare) include file not routine 9 !! - ! 00-12 (E. Kestenare M. Levy) fix bug in trtrd indexes 10 !! - ! 01-07 (E. Durand G. Madec) adaptation to ORCA config 11 11 !! 8.5 ! 02-06 (G. Madec) F90: Free form and module 12 !! 9.0 ! 04-01 (A. de Miranda, G. Madec, J.M. Molines ): advective bbl 13 !! 9.0 ! 08-04 (S. Cravatte) add the i-, j- & k- trends computation 14 !! " " ! 05-11 (V. Garnier) Surface pressure gradient organization 15 !!---------------------------------------------------------------------- 16 17 18 !!---------------------------------------------------------------------- 19 !! tra_adv_tvd : update the tracer trend with the horizontal 20 !! and vertical advection trends using a TVD scheme 21 !! nonosc : compute monotonic tracer fluxes by a nonoscillatory 22 !! algorithm 23 !!---------------------------------------------------------------------- 24 USE oce ! ocean dynamics and active tracers 12 !! NEMO 1.0 ! 04-01 (A. de Miranda, G. Madec, J.M. Molines ): advective bbl 13 !! - ! 08-04 (S. Cravatte) add the i-, j- & k- trends computation 14 !! - ! 05-11 (V. Garnier) Surface pressure gradient organization 15 !! 2.4 ! 08-01 (G. Madec) Merge TRA-TRC 16 !!---------------------------------------------------------------------- 17 18 !!---------------------------------------------------------------------- 19 !! tra_adv_tvd : update the tracer trend with the horizontal and 20 !! vertical advection trends using a TVD scheme 21 !! nonosc : compute monotonic tracer fluxes by a nonoscillatory algorithm 22 !!---------------------------------------------------------------------- 25 23 USE dom_oce ! ocean space and time domain 26 24 USE trdmod ! ocean active tracers trends … … 29 27 USE dynspg_oce ! choice/control of key cpp for surface pressure gradient 30 28 USE trabbl ! Advective term of BBL 31 USE lib_mpp 29 USE lib_mpp ! 32 30 USE lbclnk ! ocean lateral boundary condition (or mpp link) 33 31 USE diaptr ! poleward transport diagnostics 34 32 USE prtctl ! Print control 35 33 36 37 34 IMPLICIT NONE 38 35 PRIVATE 39 36 40 PUBLIC tra_adv_tvd ! routine called by step.F9037 PUBLIC tra_adv_tvd ! routine called by traadv.F90 41 38 42 39 !! * Substitutions … … 44 41 # include "vectopt_loop_substitute.h90" 45 42 !!---------------------------------------------------------------------- 46 !! OPA 9.0 , LOCEAN-IPSL (2006)47 !! $ Header$43 !! NEMO/OPA 2.4 , LOCEAN-IPSL (2008) 44 !! $Id:$ 48 45 !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 49 46 !!---------------------------------------------------------------------- … … 51 48 CONTAINS 52 49 53 SUBROUTINE tra_adv_tvd( kt, pun, pvn, pwn ) 50 SUBROUTINE tra_adv_tvd( kt, cdtype, ktra, pun, pvn, pwn, & 51 & ptb, ptn, pta ) 54 52 !!---------------------------------------------------------------------- 55 53 !! *** ROUTINE tra_adv_tvd *** … … 62 60 !! note: - this advection scheme needs a leap-frog time scheme 63 61 !! 64 !! ** Action : - update (ta,sa)with the now advective tracer trends62 !! ** Action : - update pta with the now advective tracer trends 65 63 !! - save the trends in (ztrdt,ztrds) ('key_trdtra') 66 64 !!---------------------------------------------------------------------- 67 USE oce , ztrdt => ua ! use ua as workspace 68 USE oce , ztrds => va ! use va as workspace 69 !! 70 INTEGER , INTENT(in) :: kt ! ocean time-step index 71 REAL(wp), INTENT(in), DIMENSION(jpi,jpj,jpk) :: pun ! ocean velocity u-component 72 REAL(wp), INTENT(in), DIMENSION(jpi,jpj,jpk) :: pvn ! ocean velocity v-component 73 REAL(wp), INTENT(in), DIMENSION(jpi,jpj,jpk) :: pwn ! ocean velocity w-component 65 INTEGER , INTENT(in ) :: kt ! ocean time-step index 66 CHARACTER(len=3), INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator) 67 INTEGER , INTENT(in ) :: ktra ! tracer index 68 REAL(wp) , INTENT(in ), DIMENSION(jpi,jpj,jpk) :: pun, pvn, pwn ! 3 ocean velocity components 69 REAL(wp) , INTENT(inout), DIMENSION(jpi,jpj,jpk) :: ptb, ptn ! before and now tracer fields 70 REAL(wp) , INTENT(inout), DIMENSION(jpi,jpj,jpk) :: pta ! tracer trend 74 71 !! 75 72 INTEGER :: ji, jj, jk ! dummy loop indices 76 REAL(wp) :: & ! temporary scalar 77 ztai, ztaj, ztak, & ! " " 78 zsai, zsaj, zsak, & ! " " 79 z_hdivn_x, z_hdivn_y, z_hdivn 80 REAL(wp) :: & 81 z2dtt, zbtr, zeu, zev, & ! temporary scalar 82 zew, z2, zbtr1, & ! temporary scalar 83 zfp_ui, zfp_vj, zfp_wk, & ! " " 84 zfm_ui, zfm_vj, zfm_wk ! " " 73 REAL(wp) :: ztai, ztaj, ztak 74 REAL(wp) :: z2dtt, zbtr, zeu, zev ! temporary scalar 75 REAL(wp) :: zew, z2 ! temporary scalar 76 REAL(wp) :: zfp_ui, zfp_vj, zfp_wk ! " " 77 REAL(wp) :: zfm_ui, zfm_vj, zfm_wk ! " " 85 78 REAL(wp), DIMENSION (jpi,jpj,jpk) :: zti, ztu, ztv, ztw ! temporary workspace 86 REAL(wp), DIMENSION (jpi,jpj,jpk) :: zsi, zsu, zsv, zsw ! " " 87 !!---------------------------------------------------------------------- 88 89 zti(:,:,:) = 0.e0 ; zsi(:,:,:) = 0.e0 79 !!---------------------------------------------------------------------- 80 81 zti(:,:,:) = 0.e0 90 82 91 83 IF( kt == nit000 .AND. lwp ) THEN … … 101 93 ! 1. Bottom value : flux set to zero 102 94 ! --------------- 103 ztu(:,:,jpk) = 0.e0 ; zsu(:,:,jpk) = 0.e0 104 ztv(:,:,jpk) = 0.e0 ; zsv(:,:,jpk) = 0.e0 105 ztw(:,:,jpk) = 0.e0 ; zsw(:,:,jpk) = 0.e0 106 zti(:,:,jpk) = 0.e0 ; zsi(:,:,jpk) = 0.e0 95 ztu(:,:,jpk) = 0.e0 ; ztv(:,:,jpk) = 0.e0 96 ztw(:,:,jpk) = 0.e0 ; zti(:,:,jpk) = 0.e0 107 97 108 98 … … 120 110 zfp_vj = zev + ABS( zev ) 121 111 zfm_vj = zev - ABS( zev ) 122 ztu(ji,jj,jk) = zfp_ui * tb(ji,jj,jk) + zfm_ui * tb(ji+1,jj ,jk) 123 ztv(ji,jj,jk) = zfp_vj * tb(ji,jj,jk) + zfm_vj * tb(ji ,jj+1,jk) 124 zsu(ji,jj,jk) = zfp_ui * sb(ji,jj,jk) + zfm_ui * sb(ji+1,jj ,jk) 125 zsv(ji,jj,jk) = zfp_vj * sb(ji,jj,jk) + zfm_vj * sb(ji ,jj+1,jk) 112 ztu(ji,jj,jk) = zfp_ui * ptb(ji,jj,jk) + zfm_ui * ptb(ji+1,jj ,jk) 113 ztv(ji,jj,jk) = zfp_vj * ptb(ji,jj,jk) + zfm_vj * ptb(ji ,jj+1,jk) 126 114 END DO 127 115 END DO … … 132 120 IF( lk_dynspg_rl .OR. lk_vvl ) THEN ! rigid lid or variable volume: flux set to zero 133 121 ztw(:,:,1) = 0.e0 134 zsw(:,:,1) = 0.e0135 122 ELSE ! free surface 136 DO jj = 1, jpj 137 DO ji = 1, jpi 138 zew = e1t(ji,jj) * e2t(ji,jj) * pwn(ji,jj,1) 139 ztw(ji,jj,1) = zew * tb(ji,jj,1) 140 zsw(ji,jj,1) = zew * sb(ji,jj,1) 141 END DO 142 END DO 123 ztw(:,:,1) = e1t(:,:) * e2t(:,:) * pwn(:,:,1) * ptb(:,:,1) 143 124 ENDIF 144 125 … … 150 131 zfp_wk = zew + ABS( zew ) 151 132 zfm_wk = zew - ABS( zew ) 152 ztw(ji,jj,jk) = zfp_wk * tb(ji,jj,jk) + zfm_wk * tb(ji,jj,jk-1) 153 zsw(ji,jj,jk) = zfp_wk * sb(ji,jj,jk) + zfm_wk * sb(ji,jj,jk-1) 133 ztw(ji,jj,jk) = zfp_wk * ptb(ji,jj,jk) + zfm_wk * ptb(ji,jj,jk-1) 154 134 END DO 155 135 END DO … … 165 145 ztaj = - ( ztv(ji,jj,jk) - ztv(ji ,jj-1,jk ) ) * zbtr 166 146 ztak = - ( ztw(ji,jj,jk) - ztw(ji ,jj ,jk+1) ) * zbtr 167 zsai = - ( zsu(ji,jj,jk) - zsu(ji-1,jj ,jk ) ) * zbtr168 zsaj = - ( zsv(ji,jj,jk) - zsv(ji ,jj-1,jk ) ) * zbtr169 zsak = - ( zsw(ji,jj,jk) - zsw(ji ,jj ,jk+1) ) * zbtr170 147 ! total intermediate advective trends 171 148 zti(ji,jj,jk) = ztai + ztaj + ztak 172 zsi(ji,jj,jk) = zsai + zsaj + zsak 173 END DO 174 END DO 175 END DO 176 177 178 ! Save the intermediate i / j / k advective trends for diagnostics 179 ! ------------------------------------------------------------------- 180 ! Warning : We should use zun instead of un in the computations below, but we 181 ! also use hdivn which is computed with un, vn (check ???). So we use un, vn 182 ! for consistency. Results are therefore approximate with key_trabbl_adv. 183 149 END DO 150 END DO 151 END DO 152 153 ! Save the horizontal advective trends for diagnostic 154 ! ----------------------------------------------------- 184 155 IF( l_trdtra ) THEN 185 ztrdt(:,:,:) = 0.e0 ; ztrds(:,:,:) = 0.e0 186 ! 187 ! T/S ZONAL advection trends 188 DO jk = 1, jpkm1 189 DO jj = 2, jpjm1 190 DO ji = fs_2, fs_jpim1 ! vector opt. 191 zbtr = 1. / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) ) 192 ztrdt(ji,jj,jk) = - ( ztu(ji,jj,jk) - ztu(ji-1,jj,jk) ) * zbtr 193 ztrds(ji,jj,jk) = - ( zsu(ji,jj,jk) - zsu(ji-1,jj,jk) ) * zbtr 194 END DO 195 END DO 196 END DO 197 CALL trd_mod(ztrdt, ztrds, jptra_trd_xad, 'TRA', kt) ! save the trends 198 ! 199 ! T/S MERIDIONAL advection trends 200 DO jk = 1, jpkm1 201 DO jj = 2, jpjm1 202 DO ji = fs_2, fs_jpim1 ! vector opt. 203 zbtr = 1. / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) ) 204 ztrdt(ji,jj,jk) = - ( ztv(ji,jj,jk) - ztv(ji,jj-1,jk) ) * zbtr 205 ztrds(ji,jj,jk) = - ( zsv(ji,jj,jk) - zsv(ji,jj-1,jk) ) * zbtr 206 END DO 207 END DO 208 END DO 209 CALL trd_mod(ztrdt, ztrds, jptra_trd_yad, 'TRA', kt) ! save the trends 210 ! 211 ! T/S VERTICAL advection trends 212 DO jk = 1, jpkm1 213 DO jj = 2, jpjm1 214 DO ji = fs_2, fs_jpim1 ! vector opt. 215 zbtr = 1. / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) ) 216 ztrdt(ji,jj,jk) = - ( ztw(ji,jj,jk) - ztw(ji,jj,jk+1) ) * zbtr 217 ztrds(ji,jj,jk) = - ( zsw(ji,jj,jk) - zsw(ji,jj,jk+1) ) * zbtr 218 END DO 219 END DO 220 END DO 221 CALL trd_mod(ztrdt, ztrds, jptra_trd_zad, 'TRA', kt) ! save the trends 222 ! 156 CALL trd_tra_adv( kt, ktra, jpt_trd_xad, cdtype, ztu, pun, ptn ) 157 CALL trd_tra_adv( kt, ktra, jpt_trd_yad, cdtype, ztv, pvn, ptn ) 158 CALL trd_tra_adv( kt, ktra, jpt_trd_zad, cdtype, ztw, pwn, ptn ) 223 159 ENDIF 224 160 … … 228 164 DO jj = 2, jpjm1 229 165 DO ji = fs_2, fs_jpim1 ! vector opt. 230 ta(ji,jj,jk) = ta(ji,jj,jk) + zti(ji,jj,jk) 231 sa(ji,jj,jk) = sa(ji,jj,jk) + zsi(ji,jj,jk) 232 zti (ji,jj,jk) = ( tb(ji,jj,jk) + z2dtt * zti(ji,jj,jk) ) * tmask(ji,jj,jk) 233 zsi (ji,jj,jk) = ( sb(ji,jj,jk) + z2dtt * zsi(ji,jj,jk) ) * tmask(ji,jj,jk) 234 END DO 235 END DO 236 END DO 237 238 ! Lateral boundary conditions on zti, zsi (unchanged sign) 166 pta(ji,jj,jk) = pta(ji,jj,jk) + zti(ji,jj,jk) 167 zti (ji,jj,jk) = ( ptb(ji,jj,jk) + z2dtt * zti(ji,jj,jk) ) * tmask(ji,jj,jk) 168 END DO 169 END DO 170 END DO 171 172 ! Lateral boundary conditions on zti (unchanged sign) 239 173 CALL lbc_lnk( zti, 'T', 1. ) 240 CALL lbc_lnk( zsi, 'T', 1. )241 174 242 175 … … 249 182 zeu = 0.5 * e2u(ji,jj) * fse3u(ji,jj,jk) * pun(ji,jj,jk) 250 183 zev = 0.5 * e1v(ji,jj) * fse3v(ji,jj,jk) * pvn(ji,jj,jk) 251 ztu(ji,jj,jk) = zeu * ( tn(ji,jj,jk) + tn(ji+1,jj,jk) ) - ztu(ji,jj,jk) 252 zsu(ji,jj,jk) = zeu * ( sn(ji,jj,jk) + sn(ji+1,jj,jk) ) - zsu(ji,jj,jk) 253 ztv(ji,jj,jk) = zev * ( tn(ji,jj,jk) + tn(ji,jj+1,jk) ) - ztv(ji,jj,jk) 254 zsv(ji,jj,jk) = zev * ( sn(ji,jj,jk) + sn(ji,jj+1,jk) ) - zsv(ji,jj,jk) 184 ztu(ji,jj,jk) = zeu * ( ptn(ji,jj,jk) + ptn(ji+1,jj,jk) ) - ztu(ji,jj,jk) 185 ztv(ji,jj,jk) = zev * ( ptn(ji,jj,jk) + ptn(ji,jj+1,jk) ) - ztv(ji,jj,jk) 255 186 END DO 256 187 END DO … … 260 191 ! Surface value 261 192 ztw(:,:,1) = 0.e0 262 zsw(:,:,1) = 0.e0263 193 264 194 ! Interior value … … 267 197 DO ji = 1, jpi 268 198 zew = 0.5 * e1t(ji,jj) * e2t(ji,jj) * pwn(ji,jj,jk) 269 ztw(ji,jj,jk) = zew * ( tn(ji,jj,jk) + tn(ji,jj,jk-1) ) - ztw(ji,jj,jk) 270 zsw(ji,jj,jk) = zew * ( sn(ji,jj,jk) + sn(ji,jj,jk-1) ) - zsw(ji,jj,jk) 199 ztw(ji,jj,jk) = zew * ( ptn(ji,jj,jk) + ptn(ji,jj,jk-1) ) - ztw(ji,jj,jk) 271 200 END DO 272 201 END DO … … 274 203 275 204 ! Lateral bondary conditions 276 CALL lbc_lnk( ztu, 'U', -1. ) ; CALL lbc_lnk( zsu, 'U', -1. )277 CALL lbc_lnk( ztv, 'V', -1. ) ; CALL lbc_lnk( zsv, 'V', -1. )278 CALL lbc_lnk( ztw, 'W', 1. ) ; CALL lbc_lnk( zsw, 'W', 1. )205 CALL lbc_lnk( ztu, 'U', -1. ) 206 CALL lbc_lnk( ztv, 'V', -1. ) 207 CALL lbc_lnk( ztw, 'W', 1. ) 279 208 280 209 ! 4. monotonicity algorithm 281 210 ! ------------------------- 282 CALL nonosc( tb, ztu, ztv, ztw, zti, z2 ) 283 CALL nonosc( sb, zsu, zsv, zsw, zsi, z2 ) 211 CALL nonosc( ptb, ztu, ztv, ztw, zti, z2 ) 284 212 285 213 … … 294 222 ztaj = - ( ztv(ji,jj,jk) - ztv(ji ,jj-1,jk )) * zbtr 295 223 ztak = - ( ztw(ji,jj,jk) - ztw(ji ,jj ,jk+1)) * zbtr 296 zsai = - ( zsu(ji,jj,jk) - zsu(ji-1,jj ,jk )) * zbtr297 zsaj = - ( zsv(ji,jj,jk) - zsv(ji ,jj-1,jk )) * zbtr298 zsak = - ( zsw(ji,jj,jk) - zsw(ji ,jj ,jk+1)) * zbtr299 224 300 225 ! add them to the general tracer trends 301 ta(ji,jj,jk) = ta(ji,jj,jk) + ztai + ztaj + ztak 302 sa(ji,jj,jk) = sa(ji,jj,jk) + zsai + zsaj + zsak 303 END DO 304 END DO 305 END DO 306 307 308 ! Save the advective trends for diagnostics 309 ! -------------------------------------------- 310 226 pta(ji,jj,jk) = pta(ji,jj,jk) + ztai + ztaj + ztak 227 END DO 228 END DO 229 END DO 230 231 !!gm the transport computation is wrong, the upstream part is missing ! 232 ! "zonal" mean advective heat and salt transport 233 IF( cdtype == 'TRA' .AND. ln_diaptr .AND. ( MOD( kt, nf_ptr ) == 0 ) ) THEN 234 IF( ktra == jp_tem) pht_adv(:) = ptr_vj( ztv(:,:,:) ) 235 IF( ktra == jp_sal) pst_adv(:) = ptr_vj( ztv(:,:,:) ) 236 ENDIF 237 238 ! Save the horizontal advective trends for diagnostic 239 ! ----------------------------------------------------- 311 240 IF( l_trdtra ) THEN 312 ztrdt(:,:,:) = 0.e0 ; ztrds(:,:,:) = 0.e0 313 ! 314 ! T/S ZONAL advection trends 315 DO jk = 1, jpkm1 316 DO jj = 2, jpjm1 317 DO ji = fs_2, fs_jpim1 ! vector opt. 318 !-- Compute zonal divergence by splitting hdivn (see divcur.F90) 319 ! N.B. This computation is not valid along OBCs (if any) 320 zbtr = 1. / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) ) 321 z_hdivn_x = ( e2u(ji ,jj) * fse3u(ji ,jj,jk) * pun(ji ,jj,jk) & 322 & - e2u(ji-1,jj) * fse3u(ji-1,jj,jk) * pun(ji-1,jj,jk) ) * zbtr 323 !-- Compute T/S zonal advection trends 324 ztrdt(ji,jj,jk) = - ( ztu(ji,jj,jk) - ztu(ji-1,jj,jk) ) * zbtr + tn(ji,jj,jk) * z_hdivn_x 325 ztrds(ji,jj,jk) = - ( zsu(ji,jj,jk) - zsu(ji-1,jj,jk) ) * zbtr + sn(ji,jj,jk) * z_hdivn_x 326 END DO 327 END DO 328 END DO 329 CALL trd_mod(ztrdt, ztrds, jptra_trd_xad, 'TRA', kt, cnbpas='bis') ! <<< ADD TO PREVIOUSLY COMPUTED 330 ! 331 ! T/S MERIDIONAL advection trends 332 DO jk = 1, jpkm1 333 DO jj = 2, jpjm1 334 DO ji = fs_2, fs_jpim1 ! vector opt. 335 !-- Compute merid. divergence by splitting hdivn (see divcur.F90) 336 ! N.B. This computation is not valid along OBCs (if any) 337 zbtr = 1. / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) ) 338 z_hdivn_y = ( e1v(ji, jj) * fse3v(ji,jj ,jk) * pvn(ji,jj ,jk) & 339 & - e1v(ji,jj-1) * fse3v(ji,jj-1,jk) * pvn(ji,jj-1,jk) ) * zbtr 340 !-- Compute T/S meridional advection trends 341 ztrdt(ji,jj,jk) = - ( ztv(ji,jj,jk) - ztv(ji,jj-1,jk) ) * zbtr + tn(ji,jj,jk) * z_hdivn_y 342 ztrds(ji,jj,jk) = - ( zsv(ji,jj,jk) - zsv(ji,jj-1,jk) ) * zbtr + sn(ji,jj,jk) * z_hdivn_y 343 END DO 344 END DO 345 END DO 346 CALL trd_mod(ztrdt, ztrds, jptra_trd_yad, 'TRA', kt, cnbpas='bis') ! <<< ADD TO PREVIOUSLY COMPUTED 347 ! 348 ! T/S VERTICAL advection trends 349 DO jk = 1, jpkm1 350 DO jj = 2, jpjm1 351 DO ji = fs_2, fs_jpim1 ! vector opt. 352 zbtr1 = 1. / ( e1t(ji,jj) * e2t(ji,jj) ) 353 #if defined key_zco 354 zbtr = zbtr1 355 z_hdivn_x = e2u(ji,jj)*pun(ji,jj,jk) - e2u(ji-1,jj)*pun(ji-1,jj,jk) 356 z_hdivn_y = e1v(ji,jj)*pvn(ji,jj,jk) - e1v(ji,jj-1)*pvn(ji,jj-1,jk) 357 #else 358 zbtr = zbtr1 / fse3t(ji,jj,jk) 359 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) 360 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) 361 #endif 362 z_hdivn = (z_hdivn_x + z_hdivn_y) * zbtr 363 zbtr = zbtr1 / fse3t(ji,jj,jk) 364 ztrdt(ji,jj,jk) = - ( ztw(ji,jj,jk) - ztw(ji,jj,jk+1) ) * zbtr - tn(ji,jj,jk) * z_hdivn 365 ztrds(ji,jj,jk) = - ( zsw(ji,jj,jk) - zsw(ji,jj,jk+1) ) * zbtr - sn(ji,jj,jk) * z_hdivn 366 END DO 367 END DO 368 END DO 369 CALL trd_mod(ztrdt, ztrds, jptra_trd_zad, 'TRA', kt, cnbpas='bis') ! <<< ADD TO PREVIOUSLY COMPUTED 370 ! 371 ENDIF 372 373 IF(ln_ctl) CALL prt_ctl( tab3d_1=ta, clinfo1=' tvd adv - Ta: ', mask1=tmask, & 374 & tab3d_2=sa, clinfo2= ' Sa: ', mask2=tmask, clinfo3='tra' ) 375 376 ! "zonal" mean advective heat and salt transport 377 IF( ln_diaptr .AND. ( MOD( kt, nf_ptr ) == 0 ) ) THEN 378 pht_adv(:) = ptr_vj( ztv(:,:,:) ) 379 pst_adv(:) = ptr_vj( zsv(:,:,:) ) 380 ENDIF 241 CALL trd_tra_adv( kt, ktra, jpt_trd_xad, cdtype, ztu, pun, ptn, cnbpas='bis' ) ! <<< Add to iad trend 242 CALL trd_tra_adv( kt, ktra, jpt_trd_yad, cdtype, ztv, pvn, ptn, cnbpas='bis' ) ! <<< Add to jad trend 243 CALL trd_tra_adv( kt, ktra, jpt_trd_zad, cdtype, ztw, pwn, ptn, cnbpas='bis' ) ! <<< Add to zad trend 244 ENDIF 245 246 IF(ln_ctl) CALL prt_ctl( tab3d_1=pta, clinfo1=' tvd - adv: ', mask1=tmask, clinfo3=cdtype ) 381 247 ! 382 248 END SUBROUTINE tra_adv_tvd … … 396 262 !! in-space based differencing for fluid 397 263 !!---------------------------------------------------------------------- 398 REAL(wp), INTENT( in ) :: prdt ! ??? 399 REAL(wp), DIMENSION (jpi,jpj,jpk), INTENT( inout ) :: & 400 pbef, & ! before field 401 paft, & ! after field 402 paa, & ! monotonic flux in the i direction 403 pbb, & ! monotonic flux in the j direction 404 pcc ! monotonic flux in the k direction 264 REAL(wp), INTENT(in ) :: prdt ! ??? 265 REAL(wp), INTENT(inout), DIMENSION (jpi,jpj,jpk) :: pbef, paft ! before & after field 266 REAL(wp), INTENT(inout), DIMENSION (jpi,jpj,jpk) :: paa, pbb, pcc ! monotonic flux in the 3 directions 405 267 !! 406 268 INTEGER :: ji, jj, jk ! dummy loop indices
Note: See TracChangeset
for help on using the changeset viewer.