Changeset 811 for branches/dev_001_SBC/NEMO/OPA_SRC/TRA/traadv_tvd.F90
- Timestamp:
- 2008-02-07T17:00:12+01:00 (16 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/dev_001_SBC/NEMO/OPA_SRC/TRA/traadv_tvd.F90
r699 r811 75 75 INTEGER :: ji, jj, jk ! dummy loop indices 76 76 REAL(wp) :: & ! temporary scalar 77 ztai, ztaj, ztak, & ! " " 78 zsai, zsaj, zsak, & ! " " 77 ztat, zsat, & ! " " 79 78 z_hdivn_x, z_hdivn_y, z_hdivn 80 79 REAL(wp) :: & … … 158 157 ! total advective trend 159 158 DO jk = 1, jpkm1 159 z2dtt = z2 * rdttra(jk) 160 160 DO jj = 2, jpjm1 161 161 DO ji = fs_2, fs_jpim1 ! vector opt. 162 162 zbtr = 1./ ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) ) 163 ! i- j- horizontal & k- vertical advective trends164 ztai = - ( ztu(ji,jj,jk) - ztu(ji-1,jj ,jk ) ) * zbtr165 ztaj = - ( ztv(ji,jj,jk) - ztv(ji ,jj-1,jk ) ) * zbtr166 ztak = - ( ztw(ji,jj,jk) - ztw(ji ,jj ,jk+1) ) * zbtr167 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 163 ! total intermediate advective trends 171 zti(ji,jj,jk) = ztai + ztaj + ztak 172 zsi(ji,jj,jk) = zsai + zsaj + zsak 164 ztat = - ( ztu(ji,jj,jk) - ztu(ji-1,jj ,jk ) & 165 & + ztv(ji,jj,jk) - ztv(ji ,jj-1,jk ) & 166 & + ztw(ji,jj,jk) - ztw(ji ,jj ,jk+1) ) * zbtr 167 zsat = - ( zsu(ji,jj,jk) - zsu(ji-1,jj ,jk ) & 168 & + zsv(ji,jj,jk) - zsv(ji ,jj-1,jk ) & 169 & + zsw(ji,jj,jk) - zsw(ji ,jj ,jk+1) ) * zbtr 170 ! update and guess with monotonic sheme 171 ta(ji,jj,jk) = ta(ji,jj,jk) + ztat 172 sa(ji,jj,jk) = sa(ji,jj,jk) + zsat 173 zti (ji,jj,jk) = ( tb(ji,jj,jk) + z2dtt * ztat ) * tmask(ji,jj,jk) 174 zsi (ji,jj,jk) = ( sb(ji,jj,jk) + z2dtt * zsat ) * tmask(ji,jj,jk) 173 175 END DO 174 176 END DO … … 222 224 ! 223 225 ENDIF 224 225 ! update and guess with monotonic sheme226 DO jk = 1, jpkm1227 z2dtt = z2 * rdttra(jk)228 DO jj = 2, jpjm1229 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 DO235 END DO236 END DO237 226 238 227 ! Lateral boundary conditions on zti, zsi (unchanged sign) … … 290 279 DO ji = fs_2, fs_jpim1 ! vector opt. 291 280 zbtr = 1. / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) ) 292 ! i- j- horizontal & k- vertical advective trends 293 ztai = - ( ztu(ji,jj,jk) - ztu(ji-1,jj ,jk )) * zbtr 294 ztaj = - ( ztv(ji,jj,jk) - ztv(ji ,jj-1,jk )) * zbtr 295 ztak = - ( ztw(ji,jj,jk) - ztw(ji ,jj ,jk+1)) * zbtr 296 zsai = - ( zsu(ji,jj,jk) - zsu(ji-1,jj ,jk )) * zbtr 297 zsaj = - ( zsv(ji,jj,jk) - zsv(ji ,jj-1,jk )) * zbtr 298 zsak = - ( zsw(ji,jj,jk) - zsw(ji ,jj ,jk+1)) * zbtr 299 281 ! total advective trends 282 ztat = - ( ztu(ji,jj,jk) - ztu(ji-1,jj ,jk ) & 283 & + ztv(ji,jj,jk) - ztv(ji ,jj-1,jk ) & 284 & + ztw(ji,jj,jk) - ztw(ji ,jj ,jk+1) ) * zbtr 285 zsat = - ( zsu(ji,jj,jk) - zsu(ji-1,jj ,jk ) & 286 & + zsv(ji,jj,jk) - zsv(ji ,jj-1,jk ) & 287 & + zsw(ji,jj,jk) - zsw(ji ,jj ,jk+1) ) * zbtr 300 288 ! add them to the general tracer trends 301 ta(ji,jj,jk) = ta(ji,jj,jk) + zta i + ztaj + ztak302 sa(ji,jj,jk) = sa(ji,jj,jk) + zsa i + zsaj + zsak289 ta(ji,jj,jk) = ta(ji,jj,jk) + ztat 290 sa(ji,jj,jk) = sa(ji,jj,jk) + zsat 303 291 END DO 304 292 END DO … … 396 384 !! in-space based differencing for fluid 397 385 !!---------------------------------------------------------------------- 398 REAL(wp), INTENT( in ) :: prdt 386 REAL(wp), INTENT( in ) :: prdt ! ??? 399 387 REAL(wp), DIMENSION (jpi,jpj,jpk), INTENT( inout ) :: & 400 388 pbef, & ! before field … … 407 395 INTEGER :: ikm1 408 396 REAL(wp), DIMENSION (jpi,jpj,jpk) :: zbetup, zbetdo 397 REAL(wp), DIMENSION (jpi,jpj,jpk) :: zbup, zbdo 409 398 REAL(wp) :: zpos, zneg, zbt, za, zb, zc, zbig, zrtrn, z2dtt 399 REAL(wp) :: zau, zbu, zcu, zav, zbv, zcv 400 REAL(wp) :: zup, zdo 410 401 !!---------------------------------------------------------------------- 411 402 412 403 zbig = 1.e+40 413 404 zrtrn = 1.e-15 414 zbetup(:,:,:) = 0.e0 ; zbetdo(:,:,:) = 0.e0 405 zbetup(:,:,jpk) = 0.e0 ; zbetdo(:,:,jpk) = 0.e0 406 415 407 416 408 ! Search local extrema 417 409 ! -------------------- 418 ! large negative value (-zbig) inside land 419 pbef(:,:,:) = pbef(:,:,:) * tmask(:,:,:) - zbig * ( 1.e0 - tmask(:,:,:) ) 420 paft(:,:,:) = paft(:,:,:) * tmask(:,:,:) - zbig * ( 1.e0 - tmask(:,:,:) ) 421 ! search maximum in neighbourhood 410 ! max/min of pbef & paft with large negative/positive value (-/+zbig) inside land 411 zbup = MAX( pbef * tmask - zbig * ( 1.e0 - tmask ), & 412 & paft * tmask - zbig * ( 1.e0 - tmask ) ) 413 zbdo = MIN( pbef * tmask + zbig * ( 1.e0 - tmask ), & 414 & paft * tmask + zbig * ( 1.e0 - tmask ) ) 415 422 416 DO jk = 1, jpkm1 423 417 ikm1 = MAX(jk-1,1) 424 DO jj = 2, jpjm1425 DO ji = fs_2, fs_jpim1 ! vector opt.426 zbetup(ji,jj,jk) = MAX( pbef(ji ,jj ,jk ), paft(ji ,jj ,jk ), &427 & pbef(ji-1,jj ,jk ), pbef(ji+1,jj ,jk ), &428 & paft(ji-1,jj ,jk ), paft(ji+1,jj ,jk ), &429 & pbef(ji ,jj-1,jk ), pbef(ji ,jj+1,jk ), &430 & paft(ji ,jj-1,jk ), paft(ji ,jj+1,jk ), &431 & pbef(ji ,jj ,ikm1), pbef(ji ,jj ,jk+1), &432 & paft(ji ,jj ,ikm1), paft(ji ,jj ,jk+1) )433 END DO434 END DO435 END DO436 ! large positive value (+zbig) inside land437 pbef(:,:,:) = pbef(:,:,:) * tmask(:,:,:) + zbig * ( 1.e0 - tmask(:,:,:) )438 paft(:,:,:) = paft(:,:,:) * tmask(:,:,:) + zbig * ( 1.e0 - tmask(:,:,:) )439 ! search minimum in neighbourhood440 DO jk = 1, jpkm1441 ikm1 = MAX(jk-1,1)442 DO jj = 2, jpjm1443 DO ji = fs_2, fs_jpim1 ! vector opt.444 zbetdo(ji,jj,jk) = MIN( pbef(ji ,jj ,jk ), paft(ji ,jj ,jk ), &445 & pbef(ji-1,jj ,jk ), pbef(ji+1,jj ,jk ), &446 & paft(ji-1,jj ,jk ), paft(ji+1,jj ,jk ), &447 & pbef(ji ,jj-1,jk ), pbef(ji ,jj+1,jk ), &448 & paft(ji ,jj-1,jk ), paft(ji ,jj+1,jk ), &449 & pbef(ji ,jj ,ikm1), pbef(ji ,jj ,jk+1), &450 & paft(ji ,jj ,ikm1), paft(ji ,jj ,jk+1) )451 END DO452 END DO453 END DO454 455 ! restore masked values to zero456 pbef(:,:,:) = pbef(:,:,:) * tmask(:,:,:)457 paft(:,:,:) = paft(:,:,:) * tmask(:,:,:)458 459 460 ! 2. Positive and negative part of fluxes and beta terms461 ! ------------------------------------------------------462 463 DO jk = 1, jpkm1464 418 z2dtt = prdt * rdttra(jk) 465 419 DO jj = 2, jpjm1 466 420 DO ji = fs_2, fs_jpim1 ! vector opt. 467 ! positive & negative part of the flux 421 422 ! search maximum in neighbourhood 423 zup = MAX( zbup(ji ,jj ,jk ), & 424 & zbup(ji-1,jj ,jk ), zbup(ji+1,jj ,jk ), & 425 & zbup(ji ,jj-1,jk ), zbup(ji ,jj+1,jk ), & 426 & zbup(ji ,jj ,ikm1), zbup(ji ,jj ,jk+1) ) 427 428 ! search minimum in neighbourhood 429 zdo = MIN( zbdo(ji ,jj ,jk ), & 430 & zbdo(ji-1,jj ,jk ), zbdo(ji+1,jj ,jk ), & 431 & zbdo(ji ,jj-1,jk ), zbdo(ji ,jj+1,jk ), & 432 & zbdo(ji ,jj ,ikm1), zbdo(ji ,jj ,jk+1) ) 433 434 ! positive part of the flux 468 435 zpos = MAX( 0., paa(ji-1,jj ,jk ) ) - MIN( 0., paa(ji ,jj ,jk ) ) & 469 436 & + MAX( 0., pbb(ji ,jj-1,jk ) ) - MIN( 0., pbb(ji ,jj ,jk ) ) & 470 437 & + MAX( 0., pcc(ji ,jj ,jk+1) ) - MIN( 0., pcc(ji ,jj ,jk ) ) 438 439 ! negative part of the flux 471 440 zneg = MAX( 0., paa(ji ,jj ,jk ) ) - MIN( 0., paa(ji-1,jj ,jk ) ) & 472 441 & + MAX( 0., pbb(ji ,jj ,jk ) ) - MIN( 0., pbb(ji ,jj-1,jk ) ) & 473 442 & + MAX( 0., pcc(ji ,jj ,jk ) ) - MIN( 0., pcc(ji ,jj ,jk+1) ) 443 474 444 ! up & down beta terms 475 445 zbt = e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) / z2dtt 476 zbetup(ji,jj,jk) = ( z betup(ji,jj,jk) - paft(ji,jj,jk) ) / (zpos+zrtrn) * zbt477 zbetdo(ji,jj,jk) = ( paft(ji,jj,jk) - z betdo(ji,jj,jk) ) / (zneg+zrtrn) * zbt446 zbetup(ji,jj,jk) = ( zup - paft(ji,jj,jk) ) / ( zpos + zrtrn ) * zbt 447 zbetdo(ji,jj,jk) = ( paft(ji,jj,jk) - zdo ) / ( zneg + zrtrn ) * zbt 478 448 END DO 479 449 END DO … … 490 460 DO jj = 2, jpjm1 491 461 DO ji = fs_2, fs_jpim1 ! vector opt. 492 za = MIN( 1.e0, zbetdo(ji,jj,jk), zbetup(ji+1,jj,jk) ) 493 zb = MIN( 1.e0, zbetup(ji,jj,jk), zbetdo(ji+1,jj,jk) ) 494 zc = 0.5 * ( 1.e0 + SIGN( 1.e0, paa(ji,jj,jk) ) ) 495 paa(ji,jj,jk) = paa(ji,jj,jk) * ( zc * za + ( 1.e0 - zc) * zb ) 496 497 za = MIN( 1.e0, zbetdo(ji,jj,jk), zbetup(ji,jj+1,jk) ) 498 zb = MIN( 1.e0, zbetup(ji,jj,jk), zbetdo(ji,jj+1,jk) ) 499 zc = 0.5 * ( 1.e0 + SIGN( 1.e0, pbb(ji,jj,jk) ) ) 500 pbb(ji,jj,jk) = pbb(ji,jj,jk) * ( zc * za + ( 1.e0 - zc) * zb ) 501 END DO 502 END DO 503 END DO 504 462 zau = MIN( 1.e0, zbetdo(ji,jj,jk), zbetup(ji+1,jj,jk) ) 463 zbu = MIN( 1.e0, zbetup(ji,jj,jk), zbetdo(ji+1,jj,jk) ) 464 zcu = ( 0.5 + SIGN( 0.5 , paa(ji,jj,jk) ) ) 465 paa(ji,jj,jk) = paa(ji,jj,jk) * ( zcu * zau + ( 1.e0 - zcu) * zbu ) 466 467 zav = MIN( 1.e0, zbetdo(ji,jj,jk), zbetup(ji,jj+1,jk) ) 468 zbv = MIN( 1.e0, zbetup(ji,jj,jk), zbetdo(ji,jj+1,jk) ) 469 zcv = ( 0.5 + SIGN( 0.5 , pbb(ji,jj,jk) ) ) 470 pbb(ji,jj,jk) = pbb(ji,jj,jk) * ( zcv * zav + ( 1.e0 - zcv) * zbv ) 505 471 506 472 ! monotonic flux in the k direction, i.e. pcc 507 473 ! ------------------------------------------- 508 DO jk = 2, jpkm1 509 DO jj = 2, jpjm1 510 DO ji = fs_2, fs_jpim1 ! vector opt. 511 512 za = MIN( 1., zbetdo(ji,jj,jk), zbetup(ji,jj,jk-1) ) 513 zb = MIN( 1., zbetup(ji,jj,jk), zbetdo(ji,jj,jk-1) ) 514 zc = 0.5 * ( 1.e0 + SIGN( 1.e0, pcc(ji,jj,jk) ) ) 515 pcc(ji,jj,jk) = pcc(ji,jj,jk) * ( zc * za + ( 1.e0 - zc) * zb ) 474 za = MIN( 1., zbetdo(ji,jj,jk+1), zbetup(ji,jj,jk) ) 475 zb = MIN( 1., zbetup(ji,jj,jk+1), zbetdo(ji,jj,jk) ) 476 zc = ( 0.5 + SIGN( 0.5 , pcc(ji,jj,jk+1) ) ) 477 pcc(ji,jj,jk+1) = pcc(ji,jj,jk+1) * ( zc * za + ( 1.e0 - zc) * zb ) 516 478 END DO 517 479 END DO … … 521 483 CALL lbc_lnk( paa, 'U', -1. ) ! changed sign 522 484 CALL lbc_lnk( pbb, 'V', -1. ) ! changed sign 523 CALL lbc_lnk( pcc, 'W', 1. ) ! NO changed sign524 485 ! 525 486 END SUBROUTINE nonosc
Note: See TracChangeset
for help on using the changeset viewer.