Changeset 1438 for trunk/NEMO/OPA_SRC/DYN/dynvor.F90
- Timestamp:
- 2009-05-11T16:34:47+02:00 (15 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/NEMO/OPA_SRC/DYN/dynvor.F90
r1152 r1438 5 5 !! planetary vorticity trends 6 6 !!====================================================================== 7 !! History : 1.0 ! 89-12 (P. Andrich) vor_ens: Original code 8 !! 5.0 ! 91-11 (G. Madec) vor_ene, vor_mix: Original code 9 !! 6.0 ! 96-01 (G. Madec) s-coord, suppress work arrays 10 !! 8.5 ! 02-08 (G. Madec) F90: Free form and module 11 !! 8.5 ! 04-02 (G. Madec) vor_een: Original code 12 !! 9.0 ! 03-08 (G. Madec) vor_ctl: Original code 13 !! 9.0 ! 05-11 (G. Madec) dyn_vor: Original code (new step architecture) 14 !! 9.0 ! 06-11 (G. Madec) flux form advection: add metric term 7 !! History : OPA ! 1989-12 (P. Andrich) vor_ens: Original code 8 !! 5.0 ! 1991-11 (G. Madec) vor_ene, vor_mix: Original code 9 !! 6.0 ! 1996-01 (G. Madec) s-coord, suppress work arrays 10 !! 8.5 ! 2002-08 (G. Madec) F90: Free form and module 11 !! NEMO 1.0 ! 2004-02 (G. Madec) vor_een: Original code 12 !! - ! 2003-08 (G. Madec) add vor_ctl 13 !! - ! 2005-11 (G. Madec) add dyn_vor (new step architecture) 14 !! 2.0 ! 2006-11 (G. Madec) flux form advection: add metric term 15 !! 3.2 ! 2009-04 (R. Benshila) vvl: correction of een scheme 15 16 !!---------------------------------------------------------------------- 16 17 … … 37 38 PUBLIC dyn_vor ! routine called by step.F90 38 39 39 ! !* Namelist nam_dynvor: vorticity term40 ! !!* Namelist nam_dynvor: vorticity term 40 41 LOGICAL, PUBLIC :: ln_dynvor_ene = .FALSE. !: energy conserving scheme 41 42 LOGICAL, PUBLIC :: ln_dynvor_ens = .TRUE. !: enstrophy conserving scheme … … 52 53 # include "vectopt_loop_substitute.h90" 53 54 !!---------------------------------------------------------------------- 54 !! OPA 9.0 , LOCEAN-IPSL (2006)55 !! NEMO/OPA 3,2 , LOCEAN-IPSL (2009) 55 56 !! $Id$ 56 57 !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) … … 174 175 IF(ln_ctl) CALL prt_ctl( tab3d_1=ua, clinfo1=' vor - Ua: ', mask1=umask, & 175 176 & tab3d_2=va, clinfo2= ' Va: ', mask2=vmask, clinfo3='dyn' ) 176 177 ! 177 178 END SUBROUTINE dyn_vor 178 179 … … 206 207 INTEGER , INTENT(in ) :: kt ! ocean time-step index 207 208 INTEGER , INTENT(in ) :: kvor ! =ncor (planetary) ; =ntot (total) ; 208 !! =nrvm (relative vorticity or metric)209 ! ! =nrvm (relative vorticity or metric) 209 210 REAL(wp), INTENT(inout), DIMENSION(jpi,jpj,jpk) :: pua ! total u-trend 210 211 REAL(wp), INTENT(inout), DIMENSION(jpi,jpj,jpk) :: pva ! total v-trend … … 222 223 ENDIF 223 224 224 ! Local constant initialization 225 zfact2 = 0.5 * 0.5 225 zfact2 = 0.5 * 0.5 ! Local constant initialization 226 226 227 227 !CDIR PARALLEL DO PRIVATE( zwx, zwy, zwz ) … … 229 229 DO jk = 1, jpkm1 ! Horizontal slab 230 230 ! ! =============== 231 ! 231 232 ! Potential vorticity and horizontal fluxes 232 233 ! ----------------------------------------- … … 315 316 INTEGER, INTENT(in) :: kt ! ocean timestep index 316 317 !! 317 INTEGER :: ji, jj, jk ! dummy loop indices318 INTEGER :: ji, jj, jk ! dummy loop indices 318 319 REAL(wp) :: zfact1, zua, zcua, zx1, zy1 ! temporary scalars 319 320 REAL(wp) :: zfact2, zva, zcva, zx2, zy2 ! " " … … 327 328 ENDIF 328 329 329 ! Local constant initialization 330 zfact1 = 0.5 * 0.25 330 zfact1 = 0.5 * 0.25 ! Local constant initialization 331 331 zfact2 = 0.5 * 0.5 332 332 … … 335 335 DO jk = 1, jpkm1 ! Horizontal slab 336 336 ! ! =============== 337 337 ! 338 338 ! Relative and planetary potential vorticity and horizontal fluxes 339 339 ! ---------------------------------------------------------------- … … 438 438 ENDIF 439 439 440 ! Local constant initialization 441 zfact1 = 0.5 * 0.25 440 zfact1 = 0.5 * 0.25 ! Local constant initialization 442 441 443 442 !CDIR PARALLEL DO PRIVATE( zwx, zwy, zwz ) … … 445 444 DO jk = 1, jpkm1 ! Horizontal slab 446 445 ! ! =============== 446 ! 447 447 ! Potential vorticity and horizontal fluxes 448 448 ! ----------------------------------------- … … 465 465 & + ( ( vn(ji+1,jj ,jk) + vn (ji,jj,jk) ) * ( e2v(ji+1,jj ) - e2v(ji,jj) ) & 466 466 & - ( un(ji ,jj+1,jk) + un (ji,jj,jk) ) * ( e1u(ji ,jj+1) - e1u(ji,jj) ) ) & 467 & * 0.5 / ( e1f(ji,jj) * e2f(ji,jj) ) &467 & * 0.5 / ( e1f(ji,jj) * e2f(ji,jj) ) & 468 468 & ) 469 469 END DO 470 470 END DO 471 471 END SELECT 472 472 ! 473 473 IF( ln_sco ) THEN 474 474 DO jj = 1, jpj ! caution: don't use (:,:) for this loop … … 487 487 END DO 488 488 ENDIF 489 489 ! 490 490 ! Compute and add the vorticity term trend 491 491 ! ---------------------------------------- … … 514 514 !! 515 515 !! ** Method : Trend evaluated using now fields (centered in time) 516 !! and the Arakawa and Lamb (19 XX) flux form formulation : conserves516 !! and the Arakawa and Lamb (1980) flux form formulation : conserves 517 517 !! both the horizontal kinetic energy and the potential enstrophy 518 !! when horizontal divergence is zero. 519 !! The trend of the vorticity term is given by: 520 !! * s-coordinate (ln_sco=T), the e3. are inside the derivatives: 521 !! * z-coordinate (default key), e3t=e3u=e3v, the trend becomes: 522 !! Add this trend to the general momentum trend (ua,va): 523 !! (ua,va) = (ua,va) + ( voru , vorv ) 518 !! when horizontal divergence is zero (see the NEMO documentation) 519 !! Add this trend to the general momentum trend (ua,va). 524 520 !! 525 521 !! ** Action : - Update (ua,va) with the now vorticity term trend … … 531 527 INTEGER , INTENT(in ) :: kt ! ocean time-step index 532 528 INTEGER , INTENT(in ) :: kvor ! =ncor (planetary) ; =ntot (total) ; 533 !! =nrvm (relative vorticity or metric)529 ! ! =nrvm (relative vorticity or metric) 534 530 REAL(wp), INTENT(inout), DIMENSION(jpi,jpj,jpk) :: pua ! total u-trend 535 531 REAL(wp), INTENT(inout), DIMENSION(jpi,jpj,jpk) :: pva ! total v-trend 536 532 !! 537 INTEGER :: ji, jj, jk! dummy loop indices533 INTEGER :: ji, jj, jk ! dummy loop indices 538 534 REAL(wp) :: zfac12, zua, zva ! temporary scalars 539 535 REAL(wp), DIMENSION(jpi,jpj) :: zwx, zwy, zwz ! temporary 2D workspace 540 536 REAL(wp), DIMENSION(jpi,jpj) :: ztnw, ztne, ztsw, ztse ! temporary 3D workspace 537 #if defined key_vvl 538 REAL(wp), DIMENSION(jpi,jpj,jpk) :: ze3f 539 #else 541 540 REAL(wp), DIMENSION(jpi,jpj,jpk), SAVE :: ze3f 541 #endif 542 542 !!---------------------------------------------------------------------- 543 543 … … 546 546 IF(lwp) WRITE(numout,*) 'dyn:vor_een : vorticity term: energy and enstrophy conserving scheme' 547 547 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~' 548 548 ENDIF 549 550 IF( kt == nit000 .OR. lk_vvl ) THEN ! reciprocal of e3 at F-point (masked averaging of e3t) 549 551 DO jk = 1, jpk 550 552 DO jj = 1, jpjm1 … … 559 561 ENDIF 560 562 561 ! Local constant initialization 562 zfac12 = 1.e0 / 12.e0 563 zfac12 = 1.e0 / 12.e0 ! Local constant initialization 563 564 564 565 … … 571 572 ! ----------------------------------------- 572 573 SELECT CASE( kvor ) ! vorticity considered 573 CASE ( 1 ) ; zwz(:,:) = ff(:,:) * ze3f(:,:,jk) ! planetary vorticity (Coriolis) 574 CASE ( 2 ) ; zwz(:,:) = rotn(:,:,jk) * ze3f(:,:,jk) ! relative vorticity 574 CASE ( 1 ) ! planetary vorticity (Coriolis) 575 zwz(:,:) = ff(:,:) * ze3f(:,:,jk) 576 CASE ( 2 ) ! relative vorticity 577 zwz(:,:) = rotn(:,:,jk) * ze3f(:,:,jk) 575 578 CASE ( 3 ) ! metric term 576 579 DO jj = 1, jpjm1 … … 581 584 END DO 582 585 END DO 583 CASE ( 4 ) ; zwz(:,:) = ( rotn(:,:,jk) + ff(:,:) ) * ze3f(:,:,jk) ! total (relative + planetary vorticity) 586 CASE ( 4 ) ! total (relative + planetary vorticity) 587 zwz(:,:) = ( rotn(:,:,jk) + ff(:,:) ) * ze3f(:,:,jk) 584 588 CASE ( 5 ) ! total (coriolis + metric) 585 589 DO jj = 1, jpjm1 … … 588 592 & + ( ( vn(ji+1,jj ,jk) + vn (ji,jj,jk) ) * ( e2v(ji+1,jj ) - e2v(ji,jj) ) & 589 593 & - ( un(ji ,jj+1,jk) + un (ji,jj,jk) ) * ( e1u(ji ,jj+1) - e1u(ji,jj) ) ) & 590 & * 0.5 / ( e1f(ji,jj) * e2f(ji,jj) ) &594 & * 0.5 / ( e1f(ji,jj) * e2f(ji,jj) ) & 591 595 & ) * ze3f(ji,jj,jk) 592 596 END DO … … 599 603 ! Compute and add the vorticity term trend 600 604 ! ---------------------------------------- 601 jj =2602 ztne(1,:) = 0 ; ztnw(1,:) = 0 ; ztse(1,:) = 0 ;ztsw(1,:) = 0605 jj = 2 606 ztne(1,:) = 0 ; ztnw(1,:) = 0 ; ztse(1,:) = 0 ; ztsw(1,:) = 0 603 607 DO ji = 2, jpi 604 608 ztne(ji,jj) = zwz(ji-1,jj ) + zwz(ji ,jj ) + zwz(ji ,jj-1) … … 636 640 !! 637 641 !! ** Purpose : Control the consistency between cpp options for 638 !! tracer advection schemes642 !! tracer advection schemes 639 643 !!---------------------------------------------------------------------- 640 644 INTEGER :: ioptio ! temporary integer
Note: See TracChangeset
for help on using the changeset viewer.