Changeset 503 for trunk/NEMO/OPA_SRC/TRA/traadv_muscl2.F90
- Timestamp:
- 2006-09-27T10:52:29+02:00 (18 years ago)
- File:
-
- 1 edited
-
trunk/NEMO/OPA_SRC/TRA/traadv_muscl2.F90 (modified) (10 diffs)
Legend:
- Unmodified
- Added
- Removed
-
trunk/NEMO/OPA_SRC/TRA/traadv_muscl2.F90
r457 r503 4 4 !! Ocean active tracers: horizontal & vertical advective trend 5 5 !!============================================================================== 6 !! History : 9.0 ! 02-06 (G. Madec) from traadv_muscl 7 !!---------------------------------------------------------------------- 6 8 7 9 !!---------------------------------------------------------------------- … … 9 11 !! and vertical advection trends using MUSCL2 scheme 10 12 !!---------------------------------------------------------------------- 11 !! * Modules used12 13 USE oce ! ocean dynamics and active tracers 13 14 USE dom_oce ! ocean space and time domain … … 32 33 # include "vectopt_loop_substitute.h90" 33 34 !!---------------------------------------------------------------------- 34 !! OPA 9.0 , LOCEAN-IPSL (200 5)35 !! OPA 9.0 , LOCEAN-IPSL (2006) 35 36 !! $Header$ 36 !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt37 !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 37 38 !!---------------------------------------------------------------------- 38 39 … … 50 51 !! 51 52 !! ** Action : - update (ta,sa) with the now advective tracer trends 52 !! - save trends in (ztdta,ztdsa) ('key_trdtra') 53 !! 54 !! References : 55 !! Estubier, A., and M. Levy, Notes Techn. Pole de Modelisation 56 !! IPSL, Sept. 2000 (http://www.lodyc.jussieu.fr/opa) 57 !! 58 !! History : 59 !! ! 06-00 (A.Estublier) for passive tracers 60 !! ! 01-08 (E.Durand G.Madec) adapted for T & S 61 !! 8.5 ! 02-06 (G. Madec) F90: Free form and module 62 !! 9.0 ! 04-08 (C. Talandier) New trends organization 53 !! - save trends in (ztrdt,ztrds) ('key_trdtra') 54 !! 55 !! References : Estubier, A., and M. Levy, Notes Techn. Pole de Modelisation 56 !! IPSL, Sept. 2000 (http://www.lodyc.jussieu.fr/opa) 63 57 !!---------------------------------------------------------------------- 64 !! * Arguments 65 INTEGER , INTENT( in ) :: kt ! ocean time-step index 66 REAL(wp), INTENT( in ), DIMENSION(jpi,jpj,jpk) :: & 67 pun, pvn, pwn ! now ocean velocity fields 68 69 70 !! * Local declarations 71 INTEGER :: ji, jj, jk ! dummy loop indices 58 USE oce , ztrdt => ua ! use ua as workspace 59 USE oce , ztrds => va ! use va as workspace 60 !! 61 INTEGER , INTENT(in) :: kt ! ocean time-step index 62 REAL(wp), INTENT(in), DIMENSION(jpi,jpj,jpk) :: pun ! ocean velocity u-component 63 REAL(wp), INTENT(in), DIMENSION(jpi,jpj,jpk) :: pvn ! ocean velocity v-component 64 REAL(wp), INTENT(in), DIMENSION(jpi,jpj,jpk) :: pwn ! ocean velocity w-component 65 !! 66 INTEGER :: ji, jj, jk ! dummy loop indices 72 67 REAL(wp) :: & 73 68 zu, zv, zw, zeu, zev, & … … 76 71 zzt1, zzt2, zalpha, & 77 72 zzs1, zzs2, z2, & 78 zta, zsa 79 REAL(wp), DIMENSION (jpi,jpj,jpk) :: & 80 zt1, zt2, ztp1, ztp2, & 81 zs1, zs2, zsp1, zsp2, & 82 ztdta, ztdsa 73 zta, zsa, & 74 z_hdivn_x, z_hdivn_y, z_hdivn 75 REAL(wp), DIMENSION (jpi,jpj,jpk) :: zt1, zt2, ztp1, ztp2 ! 3D workspace 76 REAL(wp), DIMENSION (jpi,jpj,jpk) :: zs1, zs2, zsp1, zsp2 ! " " 83 77 !!---------------------------------------------------------------------- 84 78 … … 89 83 ENDIF 90 84 91 IF( neuler == 0 .AND. kt == nit000 ) THEN 92 z2=1. 93 ELSE 94 z2=2. 95 ENDIF 96 97 ! Save ta and sa trends 98 IF( l_trdtra ) THEN 99 ztdta(:,:,:) = ta(:,:,:) 100 ztdsa(:,:,:) = sa(:,:,:) 101 l_adv = 'mu2' 102 ENDIF 103 85 IF( neuler == 0 .AND. kt == nit000 ) THEN ; z2 = 1. 86 ELSE ; z2 = 2. 87 ENDIF 104 88 105 89 ! I. Horizontal advective fluxes … … 194 178 zt1(ji,jj,jk) = zeu * ( zalpha * zzt1 + (1.-zalpha) * zzt2 ) 195 179 zs1(ji,jj,jk) = zeu * ( zalpha * zzs1 + (1.-zalpha) * zzs2 ) 196 180 ! 197 181 z0v = SIGN( 0.5, pvn(ji,jj,jk) ) 198 182 zalpha = 0.5 - z0v … … 271 255 CALL lbc_lnk( zt2, 'V', -1. ) ; CALL lbc_lnk( zs2, 'V', -1. ) 272 256 273 ! Save MUSCL fluxes to compute i- & j- horizontal274 ! advection trends in the MLD275 IF( l_trdtra ) THEN276 ! save i- terms277 tladi(:,:,:) = zt1(:,:,:)278 sladi(:,:,:) = zs1(:,:,:)279 ! save j- terms280 tladj(:,:,:) = zt2(:,:,:)281 sladj(:,:,:) = zs2(:,:,:)282 ENDIF283 284 257 ! Compute & add the horizontal advective trend 285 258 … … 305 278 306 279 ! Save the horizontal advective trends for diagnostic 307 308 IF( l_trdtra ) THEN 309 ! Recompute the horizontal advection zta & zsa trends computed 310 ! at the step 2. above in making the difference between the new 311 ! trends and the previous one ta()/sa - ztdta()/ztdsa() and add 312 ! the term tn()/sn()*hdivn() to recover the Uh gradh(T/S) trends 313 ztdta(:,:,:) = ta(:,:,:) - ztdta(:,:,:) + tn(:,:,:) * hdivn(:,:,:) 314 ztdsa(:,:,:) = sa(:,:,:) - ztdsa(:,:,:) + sn(:,:,:) * hdivn(:,:,:) 315 316 CALL trd_mod(ztdta, ztdsa, jpttdlad, 'TRA', kt) 317 318 ! Save the new ta and sa trends 319 ztdta(:,:,:) = ta(:,:,:) 320 ztdsa(:,:,:) = sa(:,:,:) 321 322 ENDIF 323 324 IF(ln_ctl) THEN 325 CALL prt_ctl(tab3d_1=ta, clinfo1=' muscl2 had - Ta: ', mask1=tmask, & 326 & tab3d_2=sa, clinfo2=' Sa: ', mask2=tmask, clinfo3='tra') 327 ENDIF 280 IF( l_trdtra ) THEN 281 ztrdt(:,:,:) = 0.e0 ; ztrds(:,:,:) = 0.e0 282 ! 283 ! T/S ZONAL advection trends 284 DO jk = 1, jpkm1 285 DO jj = 2, jpjm1 286 DO ji = fs_2, fs_jpim1 ! vector opt. 287 !-- Compute zonal divergence by splitting hdivn (see divcur.F90) 288 ! N.B. This computation is not valid along OBCs (if any) 289 #if defined key_zco 290 zbtr = 1. / ( e1t(ji,jj) * e2t(ji,jj) ) 291 z_hdivn_x = ( e2u(ji ,jj) * pun(ji ,jj,jk) & 292 & - e2u(ji-1,jj) * pun(ji-1,jj,jk) ) * zbtr 293 #else 294 zbtr = 1. / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) ) 295 z_hdivn_x = ( e2u(ji ,jj) * fse3u(ji ,jj,jk) * pun(ji ,jj,jk) & 296 & - e2u(ji-1,jj) * fse3u(ji-1,jj,jk) * pun(ji-1,jj,jk) ) * zbtr 297 #endif 298 ztrdt(ji,jj,jk) = - zbtr * ( zt1(ji,jj,jk) - zt1(ji-1,jj,jk) ) + tn(ji,jj,jk) * z_hdivn_x 299 ztrds(ji,jj,jk) = - zbtr * ( zs1(ji,jj,jk) - zs1(ji-1,jj,jk) ) + sn(ji,jj,jk) * z_hdivn_x 300 END DO 301 END DO 302 END DO 303 CALL trd_mod(ztrdt, ztrds, jptra_trd_xad, 'TRA', kt) 304 305 ! T/S MERIDIONAL advection trends 306 DO jk = 1, jpkm1 307 DO jj = 2, jpjm1 308 DO ji = fs_2, fs_jpim1 ! vector opt. 309 !-- Compute merid. divergence by splitting hdivn (see divcur.F90) 310 ! N.B. This computation is not valid along OBCs (if any) 311 #if defined key_zco 312 zbtr = 1. / ( e1t(ji,jj) * e2t(ji,jj) ) 313 z_hdivn_y = ( e1v(ji,jj ) * pvn(ji,jj ,jk) & 314 & - e1v(ji,jj-1) * pvn(ji,jj-1,jk) ) * zbtr 315 #else 316 zbtr = 1. / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) ) 317 z_hdivn_y = ( e1v(ji, jj) * fse3v(ji,jj ,jk) * pvn(ji,jj ,jk) & 318 & - e1v(ji,jj-1) * fse3v(ji,jj-1,jk) * pvn(ji,jj-1,jk) ) * zbtr 319 #endif 320 ztrdt(ji,jj,jk) = - zbtr * ( zt2(ji,jj,jk) - zt2(ji,jj-1,jk) ) + tn(ji,jj,jk) * z_hdivn_y 321 ztrds(ji,jj,jk) = - zbtr * ( zs2(ji,jj,jk) - zs2(ji,jj-1,jk) ) + sn(ji,jj,jk) * z_hdivn_y 322 END DO 323 END DO 324 END DO 325 CALL trd_mod(ztrdt, ztrds, jptra_trd_yad, 'TRA', kt) 326 327 ! Save the up-to-date ta and sa trends 328 ztrdt(:,:,:) = ta(:,:,:) 329 ztrds(:,:,:) = sa(:,:,:) 330 ! 331 ENDIF 332 333 IF(ln_ctl) CALL prt_ctl( tab3d_1=ta, clinfo1=' muscl2 had - Ta: ', mask1=tmask, & 334 & tab3d_2=sa, clinfo2= ' Sa: ', mask2=tmask, clinfo3='tra') 328 335 329 336 ! "zonal" mean advective heat and salt transport … … 451 458 452 459 ! Save the vertical advective trends for diagnostic 453 454 460 IF( l_trdtra ) THEN 455 461 ! Recompute the vertical advection zta & zsa trends computed 456 462 ! at the step 2. above in making the difference between the new 457 ! trends and the previous one: ta()/sa - zt dta()/ztdsa() and substract463 ! trends and the previous one: ta()/sa - ztrdt()/ztrds() and substract 458 464 ! the term tn()/sn()*hdivn() to recover the W gradz(T/S) trends 459 ztdta(:,:,:) = ta(:,:,:) - ztdta(:,:,:) - tn(:,:,:) * hdivn(:,:,:) 460 ztdsa(:,:,:) = sa(:,:,:) - ztdsa(:,:,:) - sn(:,:,:) * hdivn(:,:,:) 461 462 CALL trd_mod(ztdta, ztdsa, jpttdzad, 'TRA', kt) 463 ENDIF 464 465 IF(ln_ctl) THEN 466 CALL prt_ctl(tab3d_1=ta, clinfo1=' muscl2 zad - Ta: ', mask1=tmask, & 467 & tab3d_2=sa, clinfo2=' Sa: ', mask2=tmask, clinfo3='tra') 468 469 ENDIF 470 465 466 DO jk = 1, jpkm1 467 DO jj = 2, jpjm1 468 DO ji = fs_2, fs_jpim1 ! vector opt. 469 #if defined key_zco 470 zbtr = 1. / ( e1t(ji,jj) * e2t(ji,jj) ) 471 z_hdivn_x = e2u(ji,jj)*pun(ji,jj,jk) - e2u(ji-1,jj)*pun(ji-1,jj,jk) 472 z_hdivn_y = e1v(ji,jj)*pvn(ji,jj,jk) - e1v(ji,jj-1)*pvn(ji,jj-1,jk) 473 #else 474 zbtr = 1. / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) ) 475 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) 476 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) 477 #endif 478 z_hdivn = (z_hdivn_x + z_hdivn_y) * zbtr 479 ztrdt(ji,jj,jk) = ta(ji,jj,jk) - ztrdt(ji,jj,jk) - tn(ji,jj,jk) * z_hdivn 480 ztrds(ji,jj,jk) = sa(ji,jj,jk) - ztrds(ji,jj,jk) - sn(ji,jj,jk) * z_hdivn 481 END DO 482 END DO 483 END DO 484 CALL trd_mod(ztrdt, ztrds, jptra_trd_zad, 'TRA', kt) 485 ! 486 ENDIF 487 488 IF(ln_ctl) CALL prt_ctl( tab3d_1=ta, clinfo1=' muscl2 zad - Ta: ', mask1=tmask, & 489 & tab3d_2=sa, clinfo2= ' Sa: ', mask2=tmask, clinfo3='tra' ) 490 ! 471 491 END SUBROUTINE tra_adv_muscl2 472 492
Note: See TracChangeset
for help on using the changeset viewer.
