New URL for NEMO forge! http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.

dynvor.F90 in trunk/NEMOGCM/NEMO/OPA_SRC/DYN – NEMO

Context Navigation

source: trunk/NEMOGCM/NEMO/OPA_SRC/DYN/dynvor.F90 @ 2715

Last change on this file since 2715 was 2715, checked in by rblod, 13 years ago
First attempt to put dynamic allocation on the trunk
Property svn:keywords set to `Id`
File size: 39.3 KB

Rev	Line
[3]	1	MODULE dynvor
	2	!!======================================================================
	3	!! * MODULE dynvor *
	4	!! Ocean dynamics: Update the momentum trend with the relative and
	5	!! planetary vorticity trends
	6	!!======================================================================
[2715]	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	!! NEMO 0.5 ! 2002-08 (G. Madec) F90: Free form and module
	11	!! 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
	16	!! 3.3 ! 2010-10 (C. Ethe, G. Madec) reorganisation of initialisation phase
[503]	17	!!----------------------------------------------------------------------
[3]	18
	19	!!----------------------------------------------------------------------
[2528]	20	!! dyn_vor : Update the momentum trend with the vorticity trend
	21	!! vor_ens : enstrophy conserving scheme (ln_dynvor_ens=T)
	22	!! vor_ene : energy conserving scheme (ln_dynvor_ene=T)
	23	!! vor_mix : mixed enstrophy/energy conserving (ln_dynvor_mix=T)
	24	!! vor_een : energy and enstrophy conserving (ln_dynvor_een=T)
	25	!! dyn_vor_init : set and control of the different vorticity option
[3]	26	!!----------------------------------------------------------------------
[503]	27	USE oce ! ocean dynamics and tracers
	28	USE dom_oce ! ocean space and time domain
[643]	29	USE dynadv ! momentum advection (use ln_dynadv_vec value)
[503]	30	USE trdmod ! ocean dynamics trends
	31	USE trdmod_oce ! ocean variables trends
	32	USE lbclnk ! ocean lateral boundary conditions (or mpp link)
	33	USE prtctl ! Print control
	34	USE in_out_manager ! I/O manager
[2715]	35	USE lib_mpp
[3]	36
	37	IMPLICIT NONE
	38	PRIVATE
	39
[2528]	40	PUBLIC dyn_vor ! routine called by step.F90
	41	PUBLIC dyn_vor_init ! routine called by opa.F90
[3]	42
[1601]	43	! !!* Namelist namdyn_vor: vorticity term
[32]	44	LOGICAL, PUBLIC :: ln_dynvor_ene = .FALSE. !: energy conserving scheme
	45	LOGICAL, PUBLIC :: ln_dynvor_ens = .TRUE. !: enstrophy conserving scheme
	46	LOGICAL, PUBLIC :: ln_dynvor_mix = .FALSE. !: mixed scheme
[503]	47	LOGICAL, PUBLIC :: ln_dynvor_een = .FALSE. !: energy and enstrophy conserving scheme
[3]	48
[503]	49	INTEGER :: nvor = 0 ! type of vorticity trend used
[643]	50	INTEGER :: ncor = 1 ! coriolis
	51	INTEGER :: nrvm = 2 ! =2 relative vorticity ; =3 metric term
	52	INTEGER :: ntot = 4 ! =4 total vorticity (relative + planetary) ; =5 coriolis + metric term
[455]	53
[3]	54	!! * Substitutions
	55	# include "domzgr_substitute.h90"
	56	# include "vectopt_loop_substitute.h90"
	57	!!----------------------------------------------------------------------
[2528]	58	!! NEMO/OPA 3.3 , NEMO Consortium (2010)
[1152]	59	!! $Id$
[2715]	60	!! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
[3]	61	!!----------------------------------------------------------------------
	62	CONTAINS
	63
[455]	64	SUBROUTINE dyn_vor( kt )
[3]	65	!!----------------------------------------------------------------------
	66	!!
[455]	67	!! ** Purpose : compute the lateral ocean tracer physics.
	68	!!
	69	!! ** Action : - Update (ua,va) with the now vorticity term trend
[503]	70	!! - save the trends in (ztrdu,ztrdv) in 2 parts (relative
[455]	71	!! and planetary vorticity trends) ('key_trddyn')
[503]	72	!!----------------------------------------------------------------------
[2715]	73	USE oce, ONLY: ztrdu => ta , ztrdv => sa ! (ta,sa) used as 3D workspace
	74	!
[503]	75	INTEGER, INTENT( in ) :: kt ! ocean time-step index
[455]	76	!!----------------------------------------------------------------------
[2715]	77	!
[643]	78	! ! vorticity term
[455]	79	SELECT CASE ( nvor ) ! compute the vorticity trend and add it to the general trend
[643]	80	!
[455]	81	CASE ( -1 ) ! esopa: test all possibility with control print
[643]	82	CALL vor_ene( kt, ntot, ua, va )
[503]	83	CALL prt_ctl( tab3d_1=ua, clinfo1=' vor0 - Ua: ', mask1=umask, &
	84	& tab3d_2=va, clinfo2= ' Va: ', mask2=vmask, clinfo3='dyn' )
[643]	85	CALL vor_ens( kt, ntot, ua, va )
[503]	86	CALL prt_ctl( tab3d_1=ua, clinfo1=' vor1 - Ua: ', mask1=umask, &
	87	& tab3d_2=va, clinfo2= ' Va: ', mask2=vmask, clinfo3='dyn' )
[455]	88	CALL vor_mix( kt )
[503]	89	CALL prt_ctl( tab3d_1=ua, clinfo1=' vor2 - Ua: ', mask1=umask, &
	90	& tab3d_2=va, clinfo2= ' Va: ', mask2=vmask, clinfo3='dyn' )
[643]	91	CALL vor_een( kt, ntot, ua, va )
[503]	92	CALL prt_ctl( tab3d_1=ua, clinfo1=' vor3 - Ua: ', mask1=umask, &
	93	& tab3d_2=va, clinfo2= ' Va: ', mask2=vmask, clinfo3='dyn' )
[643]	94	!
[455]	95	CASE ( 0 ) ! energy conserving scheme
	96	IF( l_trddyn ) THEN
	97	ztrdu(:,:,:) = ua(:,:,:)
	98	ztrdv(:,:,:) = va(:,:,:)
[643]	99	CALL vor_ene( kt, nrvm, ua, va ) ! relative vorticity or metric trend
[455]	100	ztrdu(:,:,:) = ua(:,:,:) - ztrdu(:,:,:)
	101	ztrdv(:,:,:) = va(:,:,:) - ztrdv(:,:,:)
[503]	102	CALL trd_mod( ztrdu, ztrdv, jpdyn_trd_rvo, 'DYN', kt )
[455]	103	ztrdu(:,:,:) = ua(:,:,:)
	104	ztrdv(:,:,:) = va(:,:,:)
[643]	105	CALL vor_ene( kt, ncor, ua, va ) ! planetary vorticity trend
[455]	106	ztrdu(:,:,:) = ua(:,:,:) - ztrdu(:,:,:)
	107	ztrdv(:,:,:) = va(:,:,:) - ztrdv(:,:,:)
[1104]	108	CALL trd_mod( ztrdu, ztrdv, jpdyn_trd_pvo, 'DYN', kt )
[503]	109	CALL trd_mod( ztrdu, ztrdv, jpdyn_trd_dat, 'DYN', kt )
[455]	110	ELSE
[643]	111	CALL vor_ene( kt, ntot, ua, va ) ! total vorticity
[455]	112	ENDIF
[643]	113	!
[455]	114	CASE ( 1 ) ! enstrophy conserving scheme
	115	IF( l_trddyn ) THEN
	116	ztrdu(:,:,:) = ua(:,:,:)
	117	ztrdv(:,:,:) = va(:,:,:)
[643]	118	CALL vor_ens( kt, nrvm, ua, va ) ! relative vorticity or metric trend
[455]	119	ztrdu(:,:,:) = ua(:,:,:) - ztrdu(:,:,:)
	120	ztrdv(:,:,:) = va(:,:,:) - ztrdv(:,:,:)
[503]	121	CALL trd_mod( ztrdu, ztrdv, jpdyn_trd_rvo, 'DYN', kt )
[455]	122	ztrdu(:,:,:) = ua(:,:,:)
	123	ztrdv(:,:,:) = va(:,:,:)
[643]	124	CALL vor_ens( kt, ncor, ua, va ) ! planetary vorticity trend
[455]	125	ztrdu(:,:,:) = ua(:,:,:) - ztrdu(:,:,:)
	126	ztrdv(:,:,:) = va(:,:,:) - ztrdv(:,:,:)
[1104]	127	CALL trd_mod( ztrdu, ztrdv, jpdyn_trd_pvo, 'DYN', kt )
[503]	128	CALL trd_mod( ztrdu, ztrdv, jpdyn_trd_dat, 'DYN', kt )
[455]	129	ELSE
[643]	130	CALL vor_ens( kt, ntot, ua, va ) ! total vorticity
[455]	131	ENDIF
[643]	132	!
[455]	133	CASE ( 2 ) ! mixed ene-ens scheme
	134	IF( l_trddyn ) THEN
	135	ztrdu(:,:,:) = ua(:,:,:)
	136	ztrdv(:,:,:) = va(:,:,:)
[643]	137	CALL vor_ens( kt, nrvm, ua, va ) ! relative vorticity or metric trend (ens)
[455]	138	ztrdu(:,:,:) = ua(:,:,:) - ztrdu(:,:,:)
	139	ztrdv(:,:,:) = va(:,:,:) - ztrdv(:,:,:)
[503]	140	CALL trd_mod( ztrdu, ztrdv, jpdyn_trd_rvo, 'DYN', kt )
[455]	141	ztrdu(:,:,:) = ua(:,:,:)
	142	ztrdv(:,:,:) = va(:,:,:)
[643]	143	CALL vor_ene( kt, ncor, ua, va ) ! planetary vorticity trend (ene)
[455]	144	ztrdu(:,:,:) = ua(:,:,:) - ztrdu(:,:,:)
	145	ztrdv(:,:,:) = va(:,:,:) - ztrdv(:,:,:)
[1104]	146	CALL trd_mod( ztrdu, ztrdv, jpdyn_trd_pvo, 'DYN', kt )
[503]	147	CALL trd_mod( ztrdu, ztrdv, jpdyn_trd_dat, 'DYN', kt )
[455]	148	ELSE
	149	CALL vor_mix( kt ) ! total vorticity (mix=ens-ene)
	150	ENDIF
[643]	151	!
[455]	152	CASE ( 3 ) ! energy and enstrophy conserving scheme
	153	IF( l_trddyn ) THEN
	154	ztrdu(:,:,:) = ua(:,:,:)
	155	ztrdv(:,:,:) = va(:,:,:)
[643]	156	CALL vor_een( kt, nrvm, ua, va ) ! relative vorticity or metric trend
[455]	157	ztrdu(:,:,:) = ua(:,:,:) - ztrdu(:,:,:)
	158	ztrdv(:,:,:) = va(:,:,:) - ztrdv(:,:,:)
[503]	159	CALL trd_mod( ztrdu, ztrdv, jpdyn_trd_rvo, 'DYN', kt )
[455]	160	ztrdu(:,:,:) = ua(:,:,:)
	161	ztrdv(:,:,:) = va(:,:,:)
[643]	162	CALL vor_een( kt, ncor, ua, va ) ! planetary vorticity trend
[455]	163	ztrdu(:,:,:) = ua(:,:,:) - ztrdu(:,:,:)
	164	ztrdv(:,:,:) = va(:,:,:) - ztrdv(:,:,:)
[1104]	165	CALL trd_mod( ztrdu, ztrdv, jpdyn_trd_pvo, 'DYN', kt )
[503]	166	CALL trd_mod( ztrdu, ztrdv, jpdyn_trd_dat, 'DYN', kt )
[455]	167	ELSE
[643]	168	CALL vor_een( kt, ntot, ua, va ) ! total vorticity
[455]	169	ENDIF
[643]	170	!
[455]	171	END SELECT
[2715]	172	!
[455]	173	! ! print sum trends (used for debugging)
[2715]	174	IF(ln_ctl) CALL prt_ctl( tab3d_1=ua, clinfo1=' vor - Ua: ', mask1=umask, &
[455]	175	& tab3d_2=va, clinfo2= ' Va: ', mask2=vmask, clinfo3='dyn' )
[1438]	176	!
[455]	177	END SUBROUTINE dyn_vor
	178
	179
[643]	180	SUBROUTINE vor_ene( kt, kvor, pua, pva )
[455]	181	!!----------------------------------------------------------------------
	182	!! * ROUTINE vor_ene *
	183	!!
[3]	184	!! ** Purpose : Compute the now total vorticity trend and add it to
	185	!! the general trend of the momentum equation.
	186	!!
	187	!! ** Method : Trend evaluated using now fields (centered in time)
	188	!! and the Sadourny (1975) flux form formulation : conserves the
	189	!! horizontal kinetic energy.
	190	!! The trend of the vorticity term is given by:
[455]	191	!! * s-coordinate (ln_sco=T), the e3. are inside the derivatives:
[3]	192	!! voru = 1/e1u mj-1[ (rotn+f)/e3f mi(e1v*e3v vn) ]
	193	!! vorv = 1/e2v mi-1[ (rotn+f)/e3f mj(e2u*e3u un) ]
	194	!! * z-coordinate (default key), e3t=e3u=e3v, the trend becomes:
	195	!! voru = 1/e1u mj-1[ (rotn+f) mi(e1v vn) ]
	196	!! vorv = 1/e2v mi-1[ (rotn+f) mj(e2u un) ]
	197	!! Add this trend to the general momentum trend (ua,va):
	198	!! (ua,va) = (ua,va) + ( voru , vorv )
	199	!!
	200	!! ** Action : - Update (ua,va) with the now vorticity term trend
[503]	201	!! - save the trends in (ztrdu,ztrdv) in 2 parts (relative
[3]	202	!! and planetary vorticity trends) ('key_trddyn')
	203	!!
[503]	204	!! References : Sadourny, r., 1975, j. atmos. sciences, 32, 680-689.
[3]	205	!!----------------------------------------------------------------------
[2715]	206	USE wrk_nemo, ONLY: wrk_in_use, wrk_not_released
	207	USE wrk_nemo, ONLY: zwx => wrk_2d_1 , zwy => wrk_2d_2 , zwz => wrk_2d_3 ! 2D workspace
	208	!
[643]	209	INTEGER , INTENT(in ) :: kt ! ocean time-step index
	210	INTEGER , INTENT(in ) :: kvor ! =ncor (planetary) ; =ntot (total) ;
[1438]	211	! ! =nrvm (relative vorticity or metric)
[643]	212	REAL(wp), INTENT(inout), DIMENSION(jpi,jpj,jpk) :: pua ! total u-trend
	213	REAL(wp), INTENT(inout), DIMENSION(jpi,jpj,jpk) :: pva ! total v-trend
[2715]	214	!
	215	INTEGER :: ji, jj, jk ! dummy loop indices
	216	REAL(wp) :: zx1, zy1, zfact2, zx2, zy2 ! local scalars
[3]	217	!!----------------------------------------------------------------------
	218
[2715]	219	IF( wrk_in_use(2, 1,2,3) ) THEN
	220	CALL ctl_stop('dyn:vor_ene: requested workspace arrays unavailable') ; RETURN
	221	ENDIF
	222
[52]	223	IF( kt == nit000 ) THEN
	224	IF(lwp) WRITE(numout,*)
[455]	225	IF(lwp) WRITE(numout,*) 'dyn:vor_ene : vorticity term: energy conserving scheme'
	226	IF(lwp) WRITE(numout,*) '~~~~~~~~~~~'
[52]	227	ENDIF
[3]	228
[1438]	229	zfact2 = 0.5 * 0.5 ! Local constant initialization
[216]	230
[455]	231	!CDIR PARALLEL DO PRIVATE( zwx, zwy, zwz )
[3]	232	! ! ===============
	233	DO jk = 1, jpkm1 ! Horizontal slab
	234	! ! ===============
[1438]	235	!
[3]	236	! Potential vorticity and horizontal fluxes
	237	! -----------------------------------------
[643]	238	SELECT CASE( kvor ) ! vorticity considered
	239	CASE ( 1 ) ; zwz(:,:) = ff(:,:) ! planetary vorticity (Coriolis)
	240	CASE ( 2 ) ; zwz(:,:) = rotn(:,:,jk) ! relative vorticity
	241	CASE ( 3 ) ! metric term
	242	DO jj = 1, jpjm1
	243	DO ji = 1, fs_jpim1 ! vector opt.
	244	zwz(ji,jj) = ( ( vn(ji+1,jj ,jk) + vn (ji,jj,jk) ) * ( e2v(ji+1,jj ) - e2v(ji,jj) ) &
	245	& - ( un(ji ,jj+1,jk) + un (ji,jj,jk) ) * ( e1u(ji ,jj+1) - e1u(ji,jj) ) ) &
	246	& * 0.5 / ( e1f(ji,jj) * e2f(ji,jj) )
	247	END DO
	248	END DO
	249	CASE ( 4 ) ; zwz(:,:) = ( rotn(:,:,jk) + ff(:,:) ) ! total (relative + planetary vorticity)
	250	CASE ( 5 ) ! total (coriolis + metric)
	251	DO jj = 1, jpjm1
	252	DO ji = 1, fs_jpim1 ! vector opt.
	253	zwz(ji,jj) = ( ff (ji,jj) &
	254	& + ( ( vn(ji+1,jj ,jk) + vn (ji,jj,jk) ) * ( e2v(ji+1,jj ) - e2v(ji,jj) ) &
	255	& - ( un(ji ,jj+1,jk) + un (ji,jj,jk) ) * ( e1u(ji ,jj+1) - e1u(ji,jj) ) ) &
	256	& * 0.5 / ( e1f(ji,jj) * e2f(ji,jj) ) &
	257	& )
	258	END DO
	259	END DO
[455]	260	END SELECT
	261
	262	IF( ln_sco ) THEN
	263	zwz(:,:) = zwz(:,:) / fse3f(:,:,jk)
[3]	264	zwx(:,:) = e2u(:,:) * fse3u(:,:,jk) * un(:,:,jk)
	265	zwy(:,:) = e1v(:,:) * fse3v(:,:,jk) * vn(:,:,jk)
	266	ELSE
	267	zwx(:,:) = e2u(:,:) * un(:,:,jk)
	268	zwy(:,:) = e1v(:,:) * vn(:,:,jk)
	269	ENDIF
	270
	271	! Compute and add the vorticity term trend
	272	! ----------------------------------------
	273	DO jj = 2, jpjm1
	274	DO ji = fs_2, fs_jpim1 ! vector opt.
	275	zy1 = zwy(ji,jj-1) + zwy(ji+1,jj-1)
	276	zy2 = zwy(ji,jj ) + zwy(ji+1,jj )
	277	zx1 = zwx(ji-1,jj) + zwx(ji-1,jj+1)
	278	zx2 = zwx(ji ,jj) + zwx(ji ,jj+1)
[455]	279	pua(ji,jj,jk) = pua(ji,jj,jk) + zfact2 / e1u(ji,jj) * ( zwz(ji ,jj-1) * zy1 + zwz(ji,jj) * zy2 )
	280	pva(ji,jj,jk) = pva(ji,jj,jk) - zfact2 / e2v(ji,jj) * ( zwz(ji-1,jj ) * zx1 + zwz(ji,jj) * zx2 )
[3]	281	END DO
	282	END DO
	283	! ! ===============
	284	END DO ! End of slab
	285	! ! ===============
[2715]	286	IF( wrk_not_released(2, 1,2,3) ) CALL ctl_stop('dyn:vor_ene: failed to release workspace arrays')
	287	!
[455]	288	END SUBROUTINE vor_ene
[216]	289
	290
[455]	291	SUBROUTINE vor_mix( kt )
[3]	292	!!----------------------------------------------------------------------
[455]	293	!! * ROUTINE vor_mix *
[3]	294	!!
	295	!! ** Purpose : Compute the now total vorticity trend and add it to
	296	!! the general trend of the momentum equation.
	297	!!
	298	!! ** Method : Trend evaluated using now fields (centered in time)
	299	!! Mixte formulation : conserves the potential enstrophy of a hori-
	300	!! zontally non-divergent flow for (rotzu x uh), the relative vor-
	301	!! ticity term and the horizontal kinetic energy for (f x uh), the
	302	!! coriolis term. the now trend of the vorticity term is given by:
[455]	303	!! * s-coordinate (ln_sco=T), the e3. are inside the derivatives:
[3]	304	!! voru = 1/e1u mj-1(rotn/e3f) mj-1[ mi(e1v*e3v vn) ]
	305	!! +1/e1u mj-1[ f/e3f mi(e1v*e3v vn) ]
	306	!! vorv = 1/e2v mi-1(rotn/e3f) mi-1[ mj(e2u*e3u un) ]
	307	!! +1/e2v mi-1[ f/e3f mj(e2u*e3u un) ]
	308	!! * z-coordinate (default key), e3t=e3u=e3v, the trend becomes:
	309	!! voru = 1/e1u mj-1(rotn) mj-1[ mi(e1v vn) ]
	310	!! +1/e1u mj-1[ f mi(e1v vn) ]
	311	!! vorv = 1/e2v mi-1(rotn) mi-1[ mj(e2u un) ]
	312	!! +1/e2v mi-1[ f mj(e2u un) ]
	313	!! Add this now trend to the general momentum trend (ua,va):
	314	!! (ua,va) = (ua,va) + ( voru , vorv )
	315	!!
	316	!! ** Action : - Update (ua,va) arrays with the now vorticity term trend
[503]	317	!! - Save the trends in (ztrdu,ztrdv) in 2 parts (relative
[3]	318	!! and planetary vorticity trends) ('key_trddyn')
	319	!!
[503]	320	!! References : Sadourny, r., 1975, j. atmos. sciences, 32, 680-689.
[3]	321	!!----------------------------------------------------------------------
[2715]	322	USE wrk_nemo, ONLY: wrk_in_use, wrk_not_released
	323	USE wrk_nemo, ONLY: zwx => wrk_2d_4 , zwy => wrk_2d_5 , zwz => wrk_2d_6 , zww => wrk_2d_7 ! 2D workspace
	324	!
[503]	325	INTEGER, INTENT(in) :: kt ! ocean timestep index
[2715]	326	!
[1438]	327	INTEGER :: ji, jj, jk ! dummy loop indices
[2715]	328	REAL(wp) :: zfact1, zua, zcua, zx1, zy1 ! local scalars
	329	REAL(wp) :: zfact2, zva, zcva, zx2, zy2 ! - -
[3]	330	!!----------------------------------------------------------------------
	331
[2715]	332	IF( wrk_in_use(2, 4,5,6,7) ) THEN
	333	CALL ctl_stop('dyn:vor_mix: requested workspace arrays unavailable') ; RETURN
	334	ENDIF
	335
[52]	336	IF( kt == nit000 ) THEN
	337	IF(lwp) WRITE(numout,*)
[455]	338	IF(lwp) WRITE(numout,*) 'dyn:vor_mix : vorticity term: mixed energy/enstrophy conserving scheme'
	339	IF(lwp) WRITE(numout,*) '~~~~~~~~~~~'
[52]	340	ENDIF
[3]	341
[1438]	342	zfact1 = 0.5 * 0.25 ! Local constant initialization
[3]	343	zfact2 = 0.5 * 0.5
	344
[455]	345	!CDIR PARALLEL DO PRIVATE( zwx, zwy, zwz, zww )
[3]	346	! ! ===============
	347	DO jk = 1, jpkm1 ! Horizontal slab
	348	! ! ===============
[1438]	349	!
[3]	350	! Relative and planetary potential vorticity and horizontal fluxes
	351	! ----------------------------------------------------------------
[455]	352	IF( ln_sco ) THEN
[643]	353	IF( ln_dynadv_vec ) THEN
	354	zww(:,:) = rotn(:,:,jk) / fse3f(:,:,jk)
	355	ELSE
	356	DO jj = 1, jpjm1
	357	DO ji = 1, fs_jpim1 ! vector opt.
	358	zww(ji,jj) = ( ( vn(ji+1,jj ,jk) + vn (ji,jj,jk) ) * ( e2v(ji+1,jj ) - e2v(ji,jj) ) &
	359	& - ( un(ji ,jj+1,jk) + un (ji,jj,jk) ) * ( e1u(ji ,jj+1) - e1u(ji,jj) ) ) &
	360	& * 0.5 / ( e1f(ji,jj) * e2f (ji,jj) * fse3f(ji,jj,jk) )
	361	END DO
	362	END DO
	363	ENDIF
[3]	364	zwz(:,:) = ff (:,:) / fse3f(:,:,jk)
	365	zwx(:,:) = e2u(:,:) * fse3u(:,:,jk) * un(:,:,jk)
	366	zwy(:,:) = e1v(:,:) * fse3v(:,:,jk) * vn(:,:,jk)
	367	ELSE
[643]	368	IF( ln_dynadv_vec ) THEN
	369	zww(:,:) = rotn(:,:,jk)
	370	ELSE
	371	DO jj = 1, jpjm1
	372	DO ji = 1, fs_jpim1 ! vector opt.
	373	zww(ji,jj) = ( ( vn(ji+1,jj ,jk) + vn (ji,jj,jk) ) * ( e2v(ji+1,jj ) - e2v(ji,jj) ) &
	374	& - ( un(ji ,jj+1,jk) + un (ji,jj,jk) ) * ( e1u(ji ,jj+1) - e1u(ji,jj) ) ) &
	375	& * 0.5 / ( e1f(ji,jj) * e2f (ji,jj) )
	376	END DO
	377	END DO
	378	ENDIF
	379	zwz(:,:) = ff (:,:)
[3]	380	zwx(:,:) = e2u(:,:) * un(:,:,jk)
	381	zwy(:,:) = e1v(:,:) * vn(:,:,jk)
	382	ENDIF
	383
	384	! Compute and add the vorticity term trend
	385	! ----------------------------------------
	386	DO jj = 2, jpjm1
	387	DO ji = fs_2, fs_jpim1 ! vector opt.
	388	zy1 = ( zwy(ji,jj-1) + zwy(ji+1,jj-1) ) / e1u(ji,jj)
	389	zy2 = ( zwy(ji,jj ) + zwy(ji+1,jj ) ) / e1u(ji,jj)
	390	zx1 = ( zwx(ji-1,jj) + zwx(ji-1,jj+1) ) / e2v(ji,jj)
	391	zx2 = ( zwx(ji ,jj) + zwx(ji ,jj+1) ) / e2v(ji,jj)
	392	! enstrophy conserving formulation for relative vorticity term
	393	zua = zfact1 * ( zww(ji ,jj-1) + zww(ji,jj) ) * ( zy1 + zy2 )
	394	zva =-zfact1 * ( zww(ji-1,jj ) + zww(ji,jj) ) * ( zx1 + zx2 )
	395	! energy conserving formulation for planetary vorticity term
	396	zcua = zfact2 * ( zwz(ji ,jj-1) * zy1 + zwz(ji,jj) * zy2 )
	397	zcva =-zfact2 * ( zwz(ji-1,jj ) * zx1 + zwz(ji,jj) * zx2 )
[503]	398	! mixed vorticity trend added to the momentum trends
[3]	399	ua(ji,jj,jk) = ua(ji,jj,jk) + zcua + zua
	400	va(ji,jj,jk) = va(ji,jj,jk) + zcva + zva
	401	END DO
	402	END DO
	403	! ! ===============
	404	END DO ! End of slab
	405	! ! ===============
[2715]	406	IF( wrk_not_released(2, 4,5,6,7) ) CALL ctl_stop('dyn:vor_mix: failed to release workspace arrays')
	407	!
[455]	408	END SUBROUTINE vor_mix
[216]	409
	410
[643]	411	SUBROUTINE vor_ens( kt, kvor, pua, pva )
[3]	412	!!----------------------------------------------------------------------
[455]	413	!! * ROUTINE vor_ens *
[3]	414	!!
	415	!! ** Purpose : Compute the now total vorticity trend and add it to
	416	!! the general trend of the momentum equation.
	417	!!
	418	!! ** Method : Trend evaluated using now fields (centered in time)
	419	!! and the Sadourny (1975) flux FORM formulation : conserves the
	420	!! potential enstrophy of a horizontally non-divergent flow. the
	421	!! trend of the vorticity term is given by:
[455]	422	!! * s-coordinate (ln_sco=T), the e3. are inside the derivative:
[3]	423	!! voru = 1/e1u mj-1[ (rotn+f)/e3f ] mj-1[ mi(e1v*e3v vn) ]
	424	!! vorv = 1/e2v mi-1[ (rotn+f)/e3f ] mi-1[ mj(e2u*e3u un) ]
	425	!! * z-coordinate (default key), e3t=e3u=e3v, the trend becomes:
	426	!! voru = 1/e1u mj-1[ rotn+f ] mj-1[ mi(e1v vn) ]
	427	!! vorv = 1/e2v mi-1[ rotn+f ] mi-1[ mj(e2u un) ]
	428	!! Add this trend to the general momentum trend (ua,va):
	429	!! (ua,va) = (ua,va) + ( voru , vorv )
	430	!!
	431	!! ** Action : - Update (ua,va) arrays with the now vorticity term trend
[503]	432	!! - Save the trends in (ztrdu,ztrdv) in 2 parts (relative
[3]	433	!! and planetary vorticity trends) ('key_trddyn')
	434	!!
[503]	435	!! References : Sadourny, r., 1975, j. atmos. sciences, 32, 680-689.
[3]	436	!!----------------------------------------------------------------------
[2715]	437	USE wrk_nemo, ONLY: wrk_in_use, wrk_not_released
	438	USE wrk_nemo, ONLY: zwx => wrk_2d_4, zwy => wrk_2d_5, zwz => wrk_2d_6 ! 2D workspace
	439	!
[643]	440	INTEGER , INTENT(in ) :: kt ! ocean time-step index
	441	INTEGER , INTENT(in ) :: kvor ! =ncor (planetary) ; =ntot (total) ;
	442	! ! =nrvm (relative vorticity or metric)
	443	REAL(wp), INTENT(inout), DIMENSION(jpi,jpj,jpk) :: pua ! total u-trend
	444	REAL(wp), INTENT(inout), DIMENSION(jpi,jpj,jpk) :: pva ! total v-trend
[2715]	445	!
[503]	446	INTEGER :: ji, jj, jk ! dummy loop indices
	447	REAL(wp) :: zfact1, zuav, zvau ! temporary scalars
[3]	448	!!----------------------------------------------------------------------
	449
[2715]	450	IF( wrk_in_use(2, 4,5,6) ) THEN
	451	CALL ctl_stop('dyn:vor_ens: requested workspace arrays unavailable') ; RETURN
	452	END IF
	453
[52]	454	IF( kt == nit000 ) THEN
	455	IF(lwp) WRITE(numout,*)
[455]	456	IF(lwp) WRITE(numout,*) 'dyn:vor_ens : vorticity term: enstrophy conserving scheme'
	457	IF(lwp) WRITE(numout,*) '~~~~~~~~~~~'
[52]	458	ENDIF
[3]	459
[1438]	460	zfact1 = 0.5 * 0.25 ! Local constant initialization
[3]	461
[455]	462	!CDIR PARALLEL DO PRIVATE( zwx, zwy, zwz )
[3]	463	! ! ===============
	464	DO jk = 1, jpkm1 ! Horizontal slab
	465	! ! ===============
[1438]	466	!
[3]	467	! Potential vorticity and horizontal fluxes
	468	! -----------------------------------------
[643]	469	SELECT CASE( kvor ) ! vorticity considered
	470	CASE ( 1 ) ; zwz(:,:) = ff(:,:) ! planetary vorticity (Coriolis)
	471	CASE ( 2 ) ; zwz(:,:) = rotn(:,:,jk) ! relative vorticity
	472	CASE ( 3 ) ! metric term
	473	DO jj = 1, jpjm1
	474	DO ji = 1, fs_jpim1 ! vector opt.
	475	zwz(ji,jj) = ( ( vn(ji+1,jj ,jk) + vn (ji,jj,jk) ) * ( e2v(ji+1,jj ) - e2v(ji,jj) ) &
	476	& - ( un(ji ,jj+1,jk) + un (ji,jj,jk) ) * ( e1u(ji ,jj+1) - e1u(ji,jj) ) ) &
	477	& * 0.5 / ( e1f(ji,jj) * e2f(ji,jj) )
	478	END DO
	479	END DO
	480	CASE ( 4 ) ; zwz(:,:) = ( rotn(:,:,jk) + ff(:,:) ) ! total (relative + planetary vorticity)
	481	CASE ( 5 ) ! total (coriolis + metric)
	482	DO jj = 1, jpjm1
	483	DO ji = 1, fs_jpim1 ! vector opt.
	484	zwz(ji,jj) = ( ff (ji,jj) &
	485	& + ( ( vn(ji+1,jj ,jk) + vn (ji,jj,jk) ) * ( e2v(ji+1,jj ) - e2v(ji,jj) ) &
	486	& - ( un(ji ,jj+1,jk) + un (ji,jj,jk) ) * ( e1u(ji ,jj+1) - e1u(ji,jj) ) ) &
[1438]	487	& * 0.5 / ( e1f(ji,jj) * e2f(ji,jj) ) &
[643]	488	& )
	489	END DO
	490	END DO
[455]	491	END SELECT
[1438]	492	!
[455]	493	IF( ln_sco ) THEN
[3]	494	DO jj = 1, jpj ! caution: don't use (:,:) for this loop
	495	DO ji = 1, jpi ! it causes optimization problems on NEC in auto-tasking
[455]	496	zwz(ji,jj) = zwz(ji,jj) / fse3f(ji,jj,jk)
	497	zwx(ji,jj) = e2u(ji,jj) * fse3u(ji,jj,jk) * un(ji,jj,jk)
	498	zwy(ji,jj) = e1v(ji,jj) * fse3v(ji,jj,jk) * vn(ji,jj,jk)
[3]	499	END DO
	500	END DO
	501	ELSE
	502	DO jj = 1, jpj ! caution: don't use (:,:) for this loop
	503	DO ji = 1, jpi ! it causes optimization problems on NEC in auto-tasking
[455]	504	zwx(ji,jj) = e2u(ji,jj) * un(ji,jj,jk)
	505	zwy(ji,jj) = e1v(ji,jj) * vn(ji,jj,jk)
[3]	506	END DO
	507	END DO
	508	ENDIF
[1438]	509	!
[3]	510	! Compute and add the vorticity term trend
	511	! ----------------------------------------
	512	DO jj = 2, jpjm1
	513	DO ji = fs_2, fs_jpim1 ! vector opt.
[455]	514	zuav = zfact1 / e1u(ji,jj) * ( zwy(ji ,jj-1) + zwy(ji+1,jj-1) &
[503]	515	& + zwy(ji ,jj ) + zwy(ji+1,jj ) )
[455]	516	zvau =-zfact1 / e2v(ji,jj) * ( zwx(ji-1,jj ) + zwx(ji-1,jj+1) &
[503]	517	& + zwx(ji ,jj ) + zwx(ji ,jj+1) )
[455]	518	pua(ji,jj,jk) = pua(ji,jj,jk) + zuav * ( zwz(ji ,jj-1) + zwz(ji,jj) )
	519	pva(ji,jj,jk) = pva(ji,jj,jk) + zvau * ( zwz(ji-1,jj ) + zwz(ji,jj) )
[3]	520	END DO
	521	END DO
	522	! ! ===============
	523	END DO ! End of slab
	524	! ! ===============
[2715]	525	IF( wrk_not_released(2, 4,5,6) ) CALL ctl_stop('dyn:vor_ens: failed to release workspace arrays')
	526	!
[455]	527	END SUBROUTINE vor_ens
[216]	528
	529
[643]	530	SUBROUTINE vor_een( kt, kvor, pua, pva )
[108]	531	!!----------------------------------------------------------------------
[455]	532	!! * ROUTINE vor_een *
[108]	533	!!
	534	!! ** Purpose : Compute the now total vorticity trend and add it to
	535	!! the general trend of the momentum equation.
	536	!!
	537	!! ** Method : Trend evaluated using now fields (centered in time)
[1438]	538	!! and the Arakawa and Lamb (1980) flux form formulation : conserves
[108]	539	!! both the horizontal kinetic energy and the potential enstrophy
[1438]	540	!! when horizontal divergence is zero (see the NEMO documentation)
	541	!! Add this trend to the general momentum trend (ua,va).
[108]	542	!!
	543	!! ** Action : - Update (ua,va) with the now vorticity term trend
[503]	544	!! - save the trends in (ztrdu,ztrdv) in 2 parts (relative
[108]	545	!! and planetary vorticity trends) ('key_trddyn')
	546	!!
[503]	547	!! References : Arakawa and Lamb 1980, Mon. Wea. Rev., 109, 18-36
	548	!!----------------------------------------------------------------------
[2715]	549	USE wrk_nemo, ONLY: wrk_in_use, wrk_not_released
	550	USE wrk_nemo, ONLY: zwx => wrk_2d_1 , zwy => wrk_2d_2 , zwz => wrk_2d_3 ! 2D workspace
	551	USE wrk_nemo, ONLY: ztnw => wrk_2d_4 , ztne => wrk_2d_5
	552	USE wrk_nemo, ONLY: ztsw => wrk_2d_6 , ztse => wrk_2d_7
	553	#if defined key_vvl
	554	USE wrk_nemo, ONLY: ze3f => wrk_3d_1 ! 3D workspace (lk_vvl=T)
	555	#endif
	556	!
[643]	557	INTEGER , INTENT(in ) :: kt ! ocean time-step index
	558	INTEGER , INTENT(in ) :: kvor ! =ncor (planetary) ; =ntot (total) ;
[1438]	559	! ! =nrvm (relative vorticity or metric)
[643]	560	REAL(wp), INTENT(inout), DIMENSION(jpi,jpj,jpk) :: pua ! total u-trend
	561	REAL(wp), INTENT(inout), DIMENSION(jpi,jpj,jpk) :: pva ! total v-trend
[218]	562	!!
[1438]	563	INTEGER :: ji, jj, jk ! dummy loop indices
[2715]	564	INTEGER :: ierr ! local integer
	565	REAL(wp) :: zfac12, zua, zva ! local scalars
	566	#if ! defined key_vvl
	567	REAL(wp), ALLOCATABLE, DIMENSION(:,:,:), SAVE :: ze3f ! lk_vvl=F, ze3f=1/e3f saved one for all
[1438]	568	#endif
[108]	569	!!----------------------------------------------------------------------
	570
[2715]	571	IF( wrk_in_use(2, 1,2,3,4,5,6,7) .OR. wrk_in_use(3, 1) ) THEN
	572	CALL ctl_stop('dyn:vor_een: requested workspace arrays unavailable') ; RETURN
	573	ENDIF
	574
[108]	575	IF( kt == nit000 ) THEN
	576	IF(lwp) WRITE(numout,*)
[455]	577	IF(lwp) WRITE(numout,*) 'dyn:vor_een : vorticity term: energy and enstrophy conserving scheme'
	578	IF(lwp) WRITE(numout,*) '~~~~~~~~~~~'
[2715]	579	IF( .NOT.lk_vvl ) THEN
	580	ALLOCATE( ze3f(jpi,jpj,jpk) , STAT=ierr )
	581	IF( lk_mpp ) CALL mpp_sum ( ierr )
	582	IF( ierr /= 0 ) CALL ctl_stop( 'STOP', 'dyn:vor_een : unable to allocate arrays' )
	583	ENDIF
[1438]	584	ENDIF
[108]	585
[1438]	586	IF( kt == nit000 .OR. lk_vvl ) THEN ! reciprocal of e3 at F-point (masked averaging of e3t)
[108]	587	DO jk = 1, jpk
	588	DO jj = 1, jpjm1
	589	DO ji = 1, jpim1
	590	ze3f(ji,jj,jk) = ( fse3t(ji,jj+1,jk)tmask(ji,jj+1,jk) + fse3t(ji+1,jj+1,jk)tmask(ji+1,jj+1,jk) &
	591	& + fse3t(ji,jj ,jk)tmask(ji,jj ,jk) + fse3t(ji+1,jj ,jk)tmask(ji+1,jj ,jk) ) * 0.25
[2715]	592	IF( ze3f(ji,jj,jk) /= 0._wp ) ze3f(ji,jj,jk) = 1._wp / ze3f(ji,jj,jk)
[108]	593	END DO
	594	END DO
	595	END DO
	596	CALL lbc_lnk( ze3f, 'F', 1. )
	597	ENDIF
	598
[2715]	599	zfac12 = 1._wp / 12._wp ! Local constant initialization
[216]	600
[108]	601
[455]	602	!CDIR PARALLEL DO PRIVATE( zwx, zwy, zwz, ztnw, ztne, ztsw, ztse )
[108]	603	! ! ===============
	604	DO jk = 1, jpkm1 ! Horizontal slab
	605	! ! ===============
	606
	607	! Potential vorticity and horizontal fluxes
	608	! -----------------------------------------
[643]	609	SELECT CASE( kvor ) ! vorticity considered
[1438]	610	CASE ( 1 ) ! planetary vorticity (Coriolis)
	611	zwz(:,:) = ff(:,:) * ze3f(:,:,jk)
	612	CASE ( 2 ) ! relative vorticity
	613	zwz(:,:) = rotn(:,:,jk) * ze3f(:,:,jk)
[643]	614	CASE ( 3 ) ! metric term
	615	DO jj = 1, jpjm1
	616	DO ji = 1, fs_jpim1 ! vector opt.
	617	zwz(ji,jj) = ( ( vn(ji+1,jj ,jk) + vn (ji,jj,jk) ) * ( e2v(ji+1,jj ) - e2v(ji,jj) ) &
	618	& - ( un(ji ,jj+1,jk) + un (ji,jj,jk) ) * ( e1u(ji ,jj+1) - e1u(ji,jj) ) ) &
	619	& * 0.5 / ( e1f(ji,jj) * e2f(ji,jj) ) * ze3f(ji,jj,jk)
	620	END DO
	621	END DO
[1516]	622	CALL lbc_lnk( zwz, 'F', 1. )
	623	CASE ( 4 ) ! total (relative + planetary vorticity)
[1438]	624	zwz(:,:) = ( rotn(:,:,jk) + ff(:,:) ) * ze3f(:,:,jk)
[643]	625	CASE ( 5 ) ! total (coriolis + metric)
	626	DO jj = 1, jpjm1
	627	DO ji = 1, fs_jpim1 ! vector opt.
	628	zwz(ji,jj) = ( ff (ji,jj) &
	629	& + ( ( vn(ji+1,jj ,jk) + vn (ji,jj,jk) ) * ( e2v(ji+1,jj ) - e2v(ji,jj) ) &
	630	& - ( un(ji ,jj+1,jk) + un (ji,jj,jk) ) * ( e1u(ji ,jj+1) - e1u(ji,jj) ) ) &
[1438]	631	& * 0.5 / ( e1f(ji,jj) * e2f(ji,jj) ) &
[643]	632	& ) * ze3f(ji,jj,jk)
	633	END DO
	634	END DO
[1516]	635	CALL lbc_lnk( zwz, 'F', 1. )
[455]	636	END SELECT
	637
[108]	638	zwx(:,:) = e2u(:,:) * fse3u(:,:,jk) * un(:,:,jk)
	639	zwy(:,:) = e1v(:,:) * fse3v(:,:,jk) * vn(:,:,jk)
	640
	641	! Compute and add the vorticity term trend
	642	! ----------------------------------------
[1438]	643	jj = 2
	644	ztne(1,:) = 0 ; ztnw(1,:) = 0 ; ztse(1,:) = 0 ; ztsw(1,:) = 0
[108]	645	DO ji = 2, jpi
	646	ztne(ji,jj) = zwz(ji-1,jj ) + zwz(ji ,jj ) + zwz(ji ,jj-1)
	647	ztnw(ji,jj) = zwz(ji-1,jj-1) + zwz(ji-1,jj ) + zwz(ji ,jj )
	648	ztse(ji,jj) = zwz(ji ,jj ) + zwz(ji ,jj-1) + zwz(ji-1,jj-1)
	649	ztsw(ji,jj) = zwz(ji ,jj-1) + zwz(ji-1,jj-1) + zwz(ji-1,jj )
	650	END DO
	651	DO jj = 3, jpj
[1694]	652	DO ji = fs_2, jpi ! vector opt. ok because we start at jj = 3
[108]	653	ztne(ji,jj) = zwz(ji-1,jj ) + zwz(ji ,jj ) + zwz(ji ,jj-1)
	654	ztnw(ji,jj) = zwz(ji-1,jj-1) + zwz(ji-1,jj ) + zwz(ji ,jj )
	655	ztse(ji,jj) = zwz(ji ,jj ) + zwz(ji ,jj-1) + zwz(ji-1,jj-1)
	656	ztsw(ji,jj) = zwz(ji ,jj-1) + zwz(ji-1,jj-1) + zwz(ji-1,jj )
	657	END DO
	658	END DO
	659	DO jj = 2, jpjm1
	660	DO ji = fs_2, fs_jpim1 ! vector opt.
	661	zua = + zfac12 / e1u(ji,jj) * ( ztne(ji,jj ) * zwy(ji ,jj ) + ztnw(ji+1,jj) * zwy(ji+1,jj ) &
	662	& + ztse(ji,jj ) * zwy(ji ,jj-1) + ztsw(ji+1,jj) * zwy(ji+1,jj-1) )
	663	zva = - zfac12 / e2v(ji,jj) * ( ztsw(ji,jj+1) * zwx(ji-1,jj+1) + ztse(ji,jj+1) * zwx(ji ,jj+1) &
	664	& + ztnw(ji,jj ) * zwx(ji-1,jj ) + ztne(ji,jj ) * zwx(ji ,jj ) )
[455]	665	pua(ji,jj,jk) = pua(ji,jj,jk) + zua
	666	pva(ji,jj,jk) = pva(ji,jj,jk) + zva
[108]	667	END DO
	668	END DO
	669	! ! ===============
	670	END DO ! End of slab
	671	! ! ===============
[2715]	672	IF( wrk_not_released(2, 1,2,3,4,5,6,7) .OR. &
	673	wrk_not_released(3, 1) ) CALL ctl_stop('dyn:vor_een: failed to release workspace arrays')
	674	!
[455]	675	END SUBROUTINE vor_een
[216]	676
	677
[2528]	678	SUBROUTINE dyn_vor_init
[3]	679	!!---------------------------------------------------------------------
[2528]	680	!! * ROUTINE dyn_vor_init *
[3]	681	!!
	682	!! ** Purpose : Control the consistency between cpp options for
[1438]	683	!! tracer advection schemes
[3]	684	!!----------------------------------------------------------------------
[2715]	685	INTEGER :: ioptio ! local integer
	686	!!
[1601]	687	NAMELIST/namdyn_vor/ ln_dynvor_ens, ln_dynvor_ene, ln_dynvor_mix, ln_dynvor_een
[3]	688	!!----------------------------------------------------------------------
	689
[1601]	690	REWIND ( numnam ) ! Read Namelist namdyn_vor : Vorticity scheme options
	691	READ ( numnam, namdyn_vor )
[3]	692
[503]	693	IF(lwp) THEN ! Namelist print
[3]	694	WRITE(numout,*)
[2528]	695	WRITE(numout,*) 'dyn_vor_init : vorticity term : read namelist and control the consistency'
	696	WRITE(numout,*) '~~~~~~~~~~~~'
[1601]	697	WRITE(numout,*) ' Namelist namdyn_vor : oice of the vorticity term scheme'
[503]	698	WRITE(numout,*) ' energy conserving scheme ln_dynvor_ene = ', ln_dynvor_ene
	699	WRITE(numout,*) ' enstrophy conserving scheme ln_dynvor_ens = ', ln_dynvor_ens
	700	WRITE(numout,*) ' mixed enstrophy/energy conserving scheme ln_dynvor_mix = ', ln_dynvor_mix
	701	WRITE(numout,*) ' enstrophy and energy conserving scheme ln_dynvor_een = ', ln_dynvor_een
[52]	702	ENDIF
	703
[503]	704	ioptio = 0 ! Control of vorticity scheme options
	705	IF( ln_dynvor_ene ) ioptio = ioptio + 1
	706	IF( ln_dynvor_ens ) ioptio = ioptio + 1
	707	IF( ln_dynvor_mix ) ioptio = ioptio + 1
	708	IF( ln_dynvor_een ) ioptio = ioptio + 1
	709	IF( lk_esopa ) ioptio = 1
	710
	711	IF( ioptio /= 1 ) CALL ctl_stop( ' use ONE and ONLY one vorticity scheme' )
	712
[643]	713	! ! Set nvor (type of scheme for vorticity)
[503]	714	IF( ln_dynvor_ene ) nvor = 0
	715	IF( ln_dynvor_ens ) nvor = 1
	716	IF( ln_dynvor_mix ) nvor = 2
	717	IF( ln_dynvor_een ) nvor = 3
	718	IF( lk_esopa ) nvor = -1
	719
[643]	720	! ! Set ncor, nrvm, ntot (type of vorticity)
	721	IF(lwp) WRITE(numout,*)
	722	ncor = 1
	723	IF( ln_dynadv_vec ) THEN
	724	IF(lwp) WRITE(numout,*) ' Vector form advection : vorticity = Coriolis + relative vorticity'
	725	nrvm = 2
	726	ntot = 4
	727	ELSE
	728	IF(lwp) WRITE(numout,*) ' Flux form advection : vorticity = Coriolis + metric term'
	729	nrvm = 3
	730	ntot = 5
	731	ENDIF
	732
[503]	733	IF(lwp) THEN ! Print the choice
	734	WRITE(numout,*)
[643]	735	IF( nvor == 0 ) WRITE(numout,*) ' vorticity scheme : energy conserving scheme'
	736	IF( nvor == 1 ) WRITE(numout,*) ' vorticity scheme : enstrophy conserving scheme'
	737	IF( nvor == 2 ) WRITE(numout,*) ' vorticity scheme : mixed enstrophy/energy conserving scheme'
	738	IF( nvor == 3 ) WRITE(numout,*) ' vorticity scheme : energy and enstrophy conserving scheme'
[503]	739	IF( nvor == -1 ) WRITE(numout,*) ' esopa test: use all lateral physics options'
[3]	740	ENDIF
[503]	741	!
[2528]	742	END SUBROUTINE dyn_vor_init
[3]	743
[503]	744	!!==============================================================================
[3]	745	END MODULE dynvor

Note: See TracBrowser for help on using the repository browser.

Download in other formats: