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 NEMO/branches/2020/dev_r13296_HPC-07_mocavero_mpi3/src/OCE/DYN – NEMO

Context Navigation

source: NEMO/branches/2020/dev_r13296_HPC-07_mocavero_mpi3/src/OCE/DYN/dynvor.F90 @ 13630

Last change on this file since 13630 was 13630, checked in by mocavero, 4 years ago
Add neighborhood collectives calls in the NEMO src - ticket #2496
Property svn:keywords set to `Id`
File size: 53.8 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
[9528]	8	!! 5.0 ! 1991-11 (G. Madec) vor_ene, vor_mix: Original code
[2715]	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
[9019]	16	!! 3.3 ! 2010-10 (C. Ethe, G. Madec) reorganisation of initialisation phase
	17	!! 3.7 ! 2014-04 (G. Madec) trend simplification: suppress jpdyn_trd_dat vorticity
	18	!! - ! 2014-06 (G. Madec) suppression of velocity curl from in-core memory
[7646]	19	!! - ! 2016-12 (G. Madec, E. Clementi) add Stokes-Coriolis trends (ln_stcor=T)
[9019]	20	!! 4.0 ! 2017-07 (G. Madec) linear dynamics + trends diag. with Stokes-Coriolis
[9528]	21	!! - ! 2018-03 (G. Madec) add two new schemes (ln_dynvor_enT and ln_dynvor_eet)
	22	!! - ! 2018-04 (G. Madec) add pre-computed gradient for metric term calculation
[503]	23	!!----------------------------------------------------------------------
[3]	24
	25	!!----------------------------------------------------------------------
[9019]	26	!! dyn_vor : Update the momentum trend with the vorticity trend
	27	!! vor_ens : enstrophy conserving scheme (ln_dynvor_ens=T)
	28	!! vor_ene : energy conserving scheme (ln_dynvor_ene=T)
	29	!! vor_een : energy and enstrophy conserving (ln_dynvor_een=T)
	30	!! dyn_vor_init : set and control of the different vorticity option
[3]	31	!!----------------------------------------------------------------------
[503]	32	USE oce ! ocean dynamics and tracers
	33	USE dom_oce ! ocean space and time domain
[3294]	34	USE dommsk ! ocean mask
[9019]	35	USE dynadv ! momentum advection
[4990]	36	USE trd_oce ! trends: ocean variables
	37	USE trddyn ! trend manager: dynamics
[7646]	38	USE sbcwave ! Surface Waves (add Stokes-Coriolis force)
	39	USE sbc_oce , ONLY : ln_stcor ! use Stoke-Coriolis force
[5836]	40	!
[503]	41	USE lbclnk ! ocean lateral boundary conditions (or mpp link)
	42	USE prtctl ! Print control
	43	USE in_out_manager ! I/O manager
[3294]	44	USE lib_mpp ! MPP library
	45	USE timing ! Timing
[3]	46
	47	IMPLICIT NONE
	48	PRIVATE
	49
[2528]	50	PUBLIC dyn_vor ! routine called by step.F90
[5836]	51	PUBLIC dyn_vor_init ! routine called by nemogcm.F90
[3]	52
[4147]	53	! !!* Namelist namdyn_vor: vorticity term
[9528]	54	LOGICAL, PUBLIC :: ln_dynvor_ens !: enstrophy conserving scheme (ENS)
	55	LOGICAL, PUBLIC :: ln_dynvor_ene !: f-point energy conserving scheme (ENE)
	56	LOGICAL, PUBLIC :: ln_dynvor_enT !: t-point energy conserving scheme (ENT)
	57	LOGICAL, PUBLIC :: ln_dynvor_eeT !: t-point energy conserving scheme (EET)
	58	LOGICAL, PUBLIC :: ln_dynvor_een !: energy & enstrophy conserving scheme (EEN)
[5836]	59	INTEGER, PUBLIC :: nn_een_e3f !: e3f=masked averaging of e3t divided by 4 (=0) or by the sum of mask (=1)
[9528]	60	LOGICAL, PUBLIC :: ln_dynvor_mix !: mixed scheme (MIX)
[5836]	61	LOGICAL, PUBLIC :: ln_dynvor_msk !: vorticity multiplied by fmask (=T) or not (=F) (all vorticity schemes)
[3]	62
[9528]	63	INTEGER, PUBLIC :: nvor_scheme !: choice of the type of advection scheme
	64	! ! associated indices:
	65	INTEGER, PUBLIC, PARAMETER :: np_ENS = 0 ! ENS scheme
[5836]	66	INTEGER, PUBLIC, PARAMETER :: np_ENE = 1 ! ENE scheme
[9528]	67	INTEGER, PUBLIC, PARAMETER :: np_ENT = 2 ! ENT scheme (t-point vorticity)
	68	INTEGER, PUBLIC, PARAMETER :: np_EET = 3 ! EET scheme (EEN using e3t)
[5836]	69	INTEGER, PUBLIC, PARAMETER :: np_EEN = 4 ! EEN scheme
[9528]	70	INTEGER, PUBLIC, PARAMETER :: np_MIX = 5 ! MIX scheme
[455]	71
[5836]	72	INTEGER :: ncor, nrvm, ntot ! choice of calculated vorticity
	73	! ! associated indices:
[9528]	74	INTEGER, PUBLIC, PARAMETER :: np_COR = 1 ! Coriolis (planetary)
	75	INTEGER, PUBLIC, PARAMETER :: np_RVO = 2 ! relative vorticity
	76	INTEGER, PUBLIC, PARAMETER :: np_MET = 3 ! metric term
	77	INTEGER, PUBLIC, PARAMETER :: np_CRV = 4 ! relative + planetary (total vorticity)
	78	INTEGER, PUBLIC, PARAMETER :: np_CME = 5 ! Coriolis + metric term
	79
	80	REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: di_e2u_2 ! = di(e2u)/2 used in T-point metric term calculation
	81	REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: dj_e1v_2 ! = dj(e1v)/2 - - - -
[13286]	82	REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: di_e2v_2e1e2f ! = di(e2v)/(2*e1e2f) used in F-point metric term calculation
	83	REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: dj_e1u_2e1e2f ! = dj(e1u)/(2*e1e2f) - - - -
[5836]	84
	85	REAL(wp) :: r1_4 = 0.250_wp ! =1/4
	86	REAL(wp) :: r1_8 = 0.125_wp ! =1/8
	87	REAL(wp) :: r1_12 = 1._wp / 12._wp ! 1/12
	88
[3]	89	!! * Substitutions
[12377]	90	# include "do_loop_substitute.h90"
[13237]	91	# include "domzgr_substitute.h90"
	92
[3]	93	!!----------------------------------------------------------------------
[9598]	94	!! NEMO/OCE 4.0 , NEMO Consortium (2018)
[1152]	95	!! $Id$
[10068]	96	!! Software governed by the CeCILL license (see ./LICENSE)
[3]	97	!!----------------------------------------------------------------------
	98	CONTAINS
	99
[12377]	100	SUBROUTINE dyn_vor( kt, Kmm, puu, pvv, Krhs )
[3]	101	!!----------------------------------------------------------------------
	102	!!
[455]	103	!! ** Purpose : compute the lateral ocean tracer physics.
	104	!!
[12377]	105	!! ** Action : - Update (puu(:,:,:,Krhs),pvv(:,:,:,Krhs)) with the now vorticity term trend
[503]	106	!! - save the trends in (ztrdu,ztrdv) in 2 parts (relative
[4990]	107	!! and planetary vorticity trends) and send them to trd_dyn
	108	!! for futher diagnostics (l_trddyn=T)
[503]	109	!!----------------------------------------------------------------------
[12377]	110	INTEGER , INTENT( in ) :: kt ! ocean time-step index
	111	INTEGER , INTENT( in ) :: Kmm, Krhs ! ocean time level indices
	112	REAL(wp), DIMENSION(jpi,jpj,jpk,jpt), INTENT(inout) :: puu, pvv ! ocean velocity field and RHS of momentum equation
[2715]	113	!
[9019]	114	REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: ztrdu, ztrdv
[455]	115	!!----------------------------------------------------------------------
[2715]	116	!
[9019]	117	IF( ln_timing ) CALL timing_start('dyn_vor')
[3294]	118	!
[9019]	119	IF( l_trddyn ) THEN !== trend diagnostics case : split the added trend in two parts ==!
	120	!
	121	ALLOCATE( ztrdu(jpi,jpj,jpk), ztrdv(jpi,jpj,jpk) )
	122	!
[12377]	123	ztrdu(:,:,:) = puu(:,:,:,Krhs) !* planetary vorticity trend (including Stokes-Coriolis force)
	124	ztrdv(:,:,:) = pvv(:,:,:,Krhs)
[9019]	125	SELECT CASE( nvor_scheme )
[12377]	126	CASE( np_ENS ) ; CALL vor_ens( kt, Kmm, ncor, puu(:,:,:,Kmm) , pvv(:,:,:,Kmm) , puu(:,:,:,Krhs), pvv(:,:,:,Krhs) ) ! enstrophy conserving scheme
	127	IF( ln_stcor ) CALL vor_ens( kt, Kmm, ncor, usd, vsd, puu(:,:,:,Krhs), pvv(:,:,:,Krhs) ) ! add the Stokes-Coriolis trend
	128	CASE( np_ENE, np_MIX ) ; CALL vor_ene( kt, Kmm, ncor, puu(:,:,:,Kmm) , pvv(:,:,:,Kmm) , puu(:,:,:,Krhs), pvv(:,:,:,Krhs) ) ! energy conserving scheme
	129	IF( ln_stcor ) CALL vor_ene( kt, Kmm, ncor, usd, vsd, puu(:,:,:,Krhs), pvv(:,:,:,Krhs) ) ! add the Stokes-Coriolis trend
	130	CASE( np_ENT ) ; CALL vor_enT( kt, Kmm, ncor, puu(:,:,:,Kmm) , pvv(:,:,:,Kmm) , puu(:,:,:,Krhs), pvv(:,:,:,Krhs) ) ! energy conserving scheme (T-pts)
	131	IF( ln_stcor ) CALL vor_enT( kt, Kmm, ncor, usd, vsd, puu(:,:,:,Krhs), pvv(:,:,:,Krhs) ) ! add the Stokes-Coriolis trend
	132	CASE( np_EET ) ; CALL vor_eeT( kt, Kmm, ncor, puu(:,:,:,Kmm) , pvv(:,:,:,Kmm) , puu(:,:,:,Krhs), pvv(:,:,:,Krhs) ) ! energy conserving scheme (een with e3t)
	133	IF( ln_stcor ) CALL vor_eeT( kt, Kmm, ncor, usd, vsd, puu(:,:,:,Krhs), pvv(:,:,:,Krhs) ) ! add the Stokes-Coriolis trend
	134	CASE( np_EEN ) ; CALL vor_een( kt, Kmm, ncor, puu(:,:,:,Kmm) , pvv(:,:,:,Kmm) , puu(:,:,:,Krhs), pvv(:,:,:,Krhs) ) ! energy & enstrophy scheme
	135	IF( ln_stcor ) CALL vor_een( kt, Kmm, ncor, usd, vsd, puu(:,:,:,Krhs), pvv(:,:,:,Krhs) ) ! add the Stokes-Coriolis trend
[9019]	136	END SELECT
[12377]	137	ztrdu(:,:,:) = puu(:,:,:,Krhs) - ztrdu(:,:,:)
	138	ztrdv(:,:,:) = pvv(:,:,:,Krhs) - ztrdv(:,:,:)
	139	CALL trd_dyn( ztrdu, ztrdv, jpdyn_pvo, kt, Kmm )
[9019]	140	!
	141	IF( n_dynadv /= np_LIN_dyn ) THEN !* relative vorticity or metric trend (only in non-linear case)
[12377]	142	ztrdu(:,:,:) = puu(:,:,:,Krhs)
	143	ztrdv(:,:,:) = pvv(:,:,:,Krhs)
[9019]	144	SELECT CASE( nvor_scheme )
[12377]	145	CASE( np_ENT ) ; CALL vor_enT( kt, Kmm, nrvm, puu(:,:,:,Kmm) , pvv(:,:,:,Kmm) , puu(:,:,:,Krhs), pvv(:,:,:,Krhs) ) ! energy conserving scheme (T-pts)
	146	CASE( np_EET ) ; CALL vor_eeT( kt, Kmm, nrvm, puu(:,:,:,Kmm) , pvv(:,:,:,Kmm) , puu(:,:,:,Krhs), pvv(:,:,:,Krhs) ) ! energy conserving scheme (een with e3t)
	147	CASE( np_ENE ) ; CALL vor_ene( kt, Kmm, nrvm, puu(:,:,:,Kmm) , pvv(:,:,:,Kmm) , puu(:,:,:,Krhs), pvv(:,:,:,Krhs) ) ! energy conserving scheme
	148	CASE( np_ENS, np_MIX ) ; CALL vor_ens( kt, Kmm, nrvm, puu(:,:,:,Kmm) , pvv(:,:,:,Kmm) , puu(:,:,:,Krhs), pvv(:,:,:,Krhs) ) ! enstrophy conserving scheme
	149	CASE( np_EEN ) ; CALL vor_een( kt, Kmm, nrvm, puu(:,:,:,Kmm) , pvv(:,:,:,Kmm) , puu(:,:,:,Krhs), pvv(:,:,:,Krhs) ) ! energy & enstrophy scheme
[9019]	150	END SELECT
[12377]	151	ztrdu(:,:,:) = puu(:,:,:,Krhs) - ztrdu(:,:,:)
	152	ztrdv(:,:,:) = pvv(:,:,:,Krhs) - ztrdv(:,:,:)
	153	CALL trd_dyn( ztrdu, ztrdv, jpdyn_rvo, kt, Kmm )
[9019]	154	ENDIF
	155	!
	156	DEALLOCATE( ztrdu, ztrdv )
	157	!
	158	ELSE !== total vorticity trend added to the general trend ==!
	159	!
	160	SELECT CASE ( nvor_scheme ) !== vorticity trend added to the general trend ==!
[9528]	161	CASE( np_ENT ) !* energy conserving scheme (T-pts)
[12377]	162	CALL vor_enT( kt, Kmm, ntot, puu(:,:,:,Kmm) , pvv(:,:,:,Kmm) , puu(:,:,:,Krhs), pvv(:,:,:,Krhs) ) ! total vorticity trend
	163	IF( ln_stcor ) CALL vor_enT( kt, Kmm, ncor, usd, vsd, puu(:,:,:,Krhs), pvv(:,:,:,Krhs) ) ! add the Stokes-Coriolis trend
[9528]	164	CASE( np_EET ) !* energy conserving scheme (een scheme using e3t)
[12377]	165	CALL vor_eeT( kt, Kmm, ntot, puu(:,:,:,Kmm) , pvv(:,:,:,Kmm) , puu(:,:,:,Krhs), pvv(:,:,:,Krhs) ) ! total vorticity trend
	166	IF( ln_stcor ) CALL vor_eeT( kt, Kmm, ncor, usd, vsd, puu(:,:,:,Krhs), pvv(:,:,:,Krhs) ) ! add the Stokes-Coriolis trend
[9019]	167	CASE( np_ENE ) !* energy conserving scheme
[12377]	168	CALL vor_ene( kt, Kmm, ntot, puu(:,:,:,Kmm) , pvv(:,:,:,Kmm) , puu(:,:,:,Krhs), pvv(:,:,:,Krhs) ) ! total vorticity trend
	169	IF( ln_stcor ) CALL vor_ene( kt, Kmm, ncor, usd, vsd, puu(:,:,:,Krhs), pvv(:,:,:,Krhs) ) ! add the Stokes-Coriolis trend
[9019]	170	CASE( np_ENS ) !* enstrophy conserving scheme
[12377]	171	CALL vor_ens( kt, Kmm, ntot, puu(:,:,:,Kmm) , pvv(:,:,:,Kmm) , puu(:,:,:,Krhs), pvv(:,:,:,Krhs) ) ! total vorticity trend
	172	IF( ln_stcor ) CALL vor_ens( kt, Kmm, ncor, usd, vsd, puu(:,:,:,Krhs), pvv(:,:,:,Krhs) ) ! add the Stokes-Coriolis trend
[9019]	173	CASE( np_MIX ) !* mixed ene-ens scheme
[12377]	174	CALL vor_ens( kt, Kmm, nrvm, puu(:,:,:,Kmm) , pvv(:,:,:,Kmm) , puu(:,:,:,Krhs), pvv(:,:,:,Krhs) ) ! relative vorticity or metric trend (ens)
	175	CALL vor_ene( kt, Kmm, ncor, puu(:,:,:,Kmm) , pvv(:,:,:,Kmm) , puu(:,:,:,Krhs), pvv(:,:,:,Krhs) ) ! planetary vorticity trend (ene)
	176	IF( ln_stcor ) CALL vor_ene( kt, Kmm, ncor, usd, vsd, puu(:,:,:,Krhs), pvv(:,:,:,Krhs) ) ! add the Stokes-Coriolis trend
[9019]	177	CASE( np_EEN ) !* energy and enstrophy conserving scheme
[12377]	178	CALL vor_een( kt, Kmm, ntot, puu(:,:,:,Kmm) , pvv(:,:,:,Kmm) , puu(:,:,:,Krhs), pvv(:,:,:,Krhs) ) ! total vorticity trend
	179	IF( ln_stcor ) CALL vor_een( kt, Kmm, ncor, usd, vsd, puu(:,:,:,Krhs), pvv(:,:,:,Krhs) ) ! add the Stokes-Coriolis trend
[9019]	180	END SELECT
[643]	181	!
[9019]	182	ENDIF
[2715]	183	!
[455]	184	! ! print sum trends (used for debugging)
[12377]	185	IF(sn_cfctl%l_prtctl) CALL prt_ctl( tab3d_1=puu(:,:,:,Krhs), clinfo1=' vor - Ua: ', mask1=umask, &
	186	& tab3d_2=pvv(:,:,:,Krhs), clinfo2= ' Va: ', mask2=vmask, clinfo3='dyn' )
[1438]	187	!
[9019]	188	IF( ln_timing ) CALL timing_stop('dyn_vor')
[3294]	189	!
[455]	190	END SUBROUTINE dyn_vor
	191
	192
[12377]	193	SUBROUTINE vor_enT( kt, Kmm, kvor, pu, pv, pu_rhs, pv_rhs )
[9528]	194	!!----------------------------------------------------------------------
	195	!! * ROUTINE vor_enT *
	196	!!
	197	!! ** Purpose : Compute the now total vorticity trend and add it to
	198	!! the general trend of the momentum equation.
	199	!!
	200	!! ** Method : Trend evaluated using now fields (centered in time)
	201	!! and t-point evaluation of vorticity (planetary and relative).
	202	!! conserves the horizontal kinetic energy.
	203	!! The general trend of momentum is increased due to the vorticity
	204	!! term which is given by:
	205	!! voru = 1/bu mj[ ( mi(mj(bfrvor))+btf_t)/e3t mj[vn] ]
	206	!! vorv = 1/bv mi[ ( mi(mj(bfrvor))+btf_t)/e3f mj[un] ]
	207	!! where rvor is the relative vorticity at f-point
	208	!!
[12377]	209	!! ** Action : - Update (pu_rhs,pv_rhs) with the now vorticity term trend
[9528]	210	!!----------------------------------------------------------------------
	211	INTEGER , INTENT(in ) :: kt ! ocean time-step index
[12377]	212	INTEGER , INTENT(in ) :: Kmm ! ocean time level index
[9528]	213	INTEGER , INTENT(in ) :: kvor ! total, planetary, relative, or metric
	214	REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) :: pu, pv ! now velocities
	215	REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) :: pu_rhs, pv_rhs ! total v-trend
	216	!
	217	INTEGER :: ji, jj, jk ! dummy loop indices
	218	REAL(wp) :: zx1, zy1, zx2, zy2 ! local scalars
[13571]	219	REAL(wp), DIMENSION(jpi,jpj) :: zwx, zwy, zwt ! 2D workspace
	220	REAL(wp), DIMENSION(jpi,jpj,jpkm1) :: zwz ! 3D workspace, jpkm1 -> avoid lbc_lnk on jpk that is not defined
[9528]	221	!!----------------------------------------------------------------------
	222	!
	223	IF( kt == nit000 ) THEN
	224	IF(lwp) WRITE(numout,*)
	225	IF(lwp) WRITE(numout,*) 'dyn:vor_enT : vorticity term: t-point energy conserving scheme'
	226	IF(lwp) WRITE(numout,*) '~~~~~~~~~~~'
	227	ENDIF
	228	!
[10425]	229	!
	230	SELECT CASE( kvor ) !== volume weighted vorticity considered ==!
	231	CASE ( np_RVO ) !* relative vorticity
	232	DO jk = 1, jpkm1 ! Horizontal slab
[13295]	233	DO_2D( 1, 0, 1, 0 )
[12377]	234	zwz(ji,jj,jk) = ( e2v(ji+1,jj) * pv(ji+1,jj,jk) - e2v(ji,jj) * pv(ji,jj,jk) &
	235	& - e1u(ji,jj+1) * pu(ji,jj+1,jk) + e1u(ji,jj) * pu(ji,jj,jk) ) * r1_e1e2f(ji,jj)
	236	END_2D
[9528]	237	IF( ln_dynvor_msk ) THEN ! mask/unmask relative vorticity
[13295]	238	DO_2D( 1, 0, 1, 0 )
[12377]	239	zwz(ji,jj,jk) = zwz(ji,jj,jk) * fmask(ji,jj,jk)
	240	END_2D
[9528]	241	ENDIF
[10425]	242	END DO
	243
[13630]	244	#if defined key_mpi3
	245	CALL lbc_lnk_nc_multi( 'dynvor', zwz, 'F', 1.0_wp )
	246	#else
[13226]	247	CALL lbc_lnk( 'dynvor', zwz, 'F', 1.0_wp )
[13630]	248	#endif
[10425]	249
	250	CASE ( np_CRV ) !* Coriolis + relative vorticity
	251	DO jk = 1, jpkm1 ! Horizontal slab
[13571]	252	DO_2D( 1, 0, 1, 0 ) ! relative vorticity
[12377]	253	zwz(ji,jj,jk) = ( e2v(ji+1,jj) * pv(ji+1,jj,jk) - e2v(ji,jj) * pv(ji,jj,jk) &
	254	& - e1u(ji,jj+1) * pu(ji,jj+1,jk) + e1u(ji,jj) * pu(ji,jj,jk) ) * r1_e1e2f(ji,jj)
	255	END_2D
[10425]	256	IF( ln_dynvor_msk ) THEN ! mask/unmask relative vorticity
[13295]	257	DO_2D( 1, 0, 1, 0 )
[12377]	258	zwz(ji,jj,jk) = zwz(ji,jj,jk) * fmask(ji,jj,jk)
	259	END_2D
[10425]	260	ENDIF
	261	END DO
	262
[13630]	263	#if defined key_mpi3
	264	CALL lbc_lnk_nc_multi( 'dynvor', zwz, 'F', 1.0_wp )
	265	#else
[13226]	266	CALL lbc_lnk( 'dynvor', zwz, 'F', 1.0_wp )
[13630]	267	#endif
[10425]	268
	269	END SELECT
	270
	271	! ! ===============
	272	DO jk = 1, jpkm1 ! Horizontal slab
	273	! ! ===============
	274
	275	SELECT CASE( kvor ) !== volume weighted vorticity considered ==!
	276	CASE ( np_COR ) !* Coriolis (planetary vorticity)
[12377]	277	zwt(:,:) = ff_t(:,:) * e1e2t(:,:)*e3t(:,:,jk,Kmm)
[10425]	278	CASE ( np_RVO ) !* relative vorticity
[13295]	279	DO_2D( 0, 1, 0, 1 )
[12377]	280	zwt(ji,jj) = r1_4 * ( zwz(ji-1,jj ,jk) + zwz(ji,jj ,jk) &
[13237]	281	& + zwz(ji-1,jj-1,jk) + zwz(ji,jj-1,jk) ) &
	282	& * e1e2t(ji,jj)*e3t(ji,jj,jk,Kmm)
[12377]	283	END_2D
[9528]	284	CASE ( np_MET ) !* metric term
[13295]	285	DO_2D( 0, 1, 0, 1 )
[12377]	286	zwt(ji,jj) = ( ( pv(ji,jj,jk) + pv(ji,jj-1,jk) ) * di_e2u_2(ji,jj) &
[13237]	287	& - ( pu(ji,jj,jk) + pu(ji-1,jj,jk) ) * dj_e1v_2(ji,jj) ) &
	288	& * e3t(ji,jj,jk,Kmm)
[12377]	289	END_2D
[9528]	290	CASE ( np_CRV ) !* Coriolis + relative vorticity
[13295]	291	DO_2D( 0, 1, 0, 1 )
[12377]	292	zwt(ji,jj) = ( ff_t(ji,jj) + r1_4 * ( zwz(ji-1,jj ,jk) + zwz(ji,jj ,jk) &
[13237]	293	& + zwz(ji-1,jj-1,jk) + zwz(ji,jj-1,jk) ) ) &
	294	& * e1e2t(ji,jj)*e3t(ji,jj,jk,Kmm)
[12377]	295	END_2D
[9528]	296	CASE ( np_CME ) !* Coriolis + metric
[13295]	297	DO_2D( 0, 1, 0, 1 )
[12377]	298	zwt(ji,jj) = ( ff_t(ji,jj) * e1e2t(ji,jj) &
	299	& + ( pv(ji,jj,jk) + pv(ji,jj-1,jk) ) * di_e2u_2(ji,jj) &
[13237]	300	& - ( pu(ji,jj,jk) + pu(ji-1,jj,jk) ) * dj_e1v_2(ji,jj) ) &
	301	& * e3t(ji,jj,jk,Kmm)
[12377]	302	END_2D
[9528]	303	CASE DEFAULT ! error
	304	CALL ctl_stop('STOP','dyn_vor: wrong value for kvor' )
	305	END SELECT
	306	!
	307	! !== compute and add the vorticity term trend =!
[13295]	308	DO_2D( 0, 0, 0, 0 )
[12377]	309	pu_rhs(ji,jj,jk) = pu_rhs(ji,jj,jk) + r1_4 * r1_e1e2u(ji,jj) / e3u(ji,jj,jk,Kmm) &
	310	& * ( zwt(ji+1,jj) * ( pv(ji+1,jj,jk) + pv(ji+1,jj-1,jk) ) &
	311	& + zwt(ji ,jj) * ( pv(ji ,jj,jk) + pv(ji ,jj-1,jk) ) )
	312	!
	313	pv_rhs(ji,jj,jk) = pv_rhs(ji,jj,jk) - r1_4 * r1_e1e2v(ji,jj) / e3v(ji,jj,jk,Kmm) &
	314	& * ( zwt(ji,jj+1) * ( pu(ji,jj+1,jk) + pu(ji-1,jj+1,jk) ) &
	315	& + zwt(ji,jj ) * ( pu(ji,jj ,jk) + pu(ji-1,jj ,jk) ) )
	316	END_2D
[9528]	317	! ! ===============
	318	END DO ! End of slab
	319	! ! ===============
	320	END SUBROUTINE vor_enT
	321
	322
[12377]	323	SUBROUTINE vor_ene( kt, Kmm, kvor, pu, pv, pu_rhs, pv_rhs )
[455]	324	!!----------------------------------------------------------------------
	325	!! * ROUTINE vor_ene *
	326	!!
[3]	327	!! ** Purpose : Compute the now total vorticity trend and add it to
	328	!! the general trend of the momentum equation.
	329	!!
	330	!! ** Method : Trend evaluated using now fields (centered in time)
[5836]	331	!! and the Sadourny (1975) flux form formulation : conserves the
	332	!! horizontal kinetic energy.
	333	!! The general trend of momentum is increased due to the vorticity
	334	!! term which is given by:
[12377]	335	!! voru = 1/e1u mj-1[ (rvor+f)/e3f mi(e1v*e3v pvv(:,:,:,Kmm)) ]
	336	!! vorv = 1/e2v mi-1[ (rvor+f)/e3f mj(e2u*e3u puu(:,:,:,Kmm)) ]
[5836]	337	!! where rvor is the relative vorticity
[3]	338	!!
[12377]	339	!! ** Action : - Update (pu_rhs,pv_rhs) with the now vorticity term trend
[3]	340	!!
[503]	341	!! References : Sadourny, r., 1975, j. atmos. sciences, 32, 680-689.
[3]	342	!!----------------------------------------------------------------------
[9019]	343	INTEGER , INTENT(in ) :: kt ! ocean time-step index
[12377]	344	INTEGER , INTENT(in ) :: Kmm ! ocean time level index
[9019]	345	INTEGER , INTENT(in ) :: kvor ! total, planetary, relative, or metric
[12377]	346	REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) :: pu, pv ! now velocities
	347	REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) :: pu_rhs, pv_rhs ! total v-trend
[2715]	348	!
[5836]	349	INTEGER :: ji, jj, jk ! dummy loop indices
	350	REAL(wp) :: zx1, zy1, zx2, zy2 ! local scalars
[9019]	351	REAL(wp), DIMENSION(jpi,jpj) :: zwx, zwy, zwz ! 2D workspace
[3]	352	!!----------------------------------------------------------------------
[3294]	353	!
[52]	354	IF( kt == nit000 ) THEN
	355	IF(lwp) WRITE(numout,*)
[455]	356	IF(lwp) WRITE(numout,*) 'dyn:vor_ene : vorticity term: energy conserving scheme'
	357	IF(lwp) WRITE(numout,*) '~~~~~~~~~~~'
[52]	358	ENDIF
[5836]	359	!
[3]	360	! ! ===============
	361	DO jk = 1, jpkm1 ! Horizontal slab
	362	! ! ===============
[1438]	363	!
[5836]	364	SELECT CASE( kvor ) !== vorticity considered ==!
	365	CASE ( np_COR ) !* Coriolis (planetary vorticity)
[7753]	366	zwz(:,:) = ff_f(:,:)
[5836]	367	CASE ( np_RVO ) !* relative vorticity
[13295]	368	DO_2D( 1, 0, 1, 0 )
[12377]	369	zwz(ji,jj) = ( e2v(ji+1,jj ) * pv(ji+1,jj ,jk) - e2v(ji,jj) * pv(ji,jj,jk) &
	370	& - e1u(ji ,jj+1) * pu(ji ,jj+1,jk) + e1u(ji,jj) * pu(ji,jj,jk) ) * r1_e1e2f(ji,jj)
	371	END_2D
[5836]	372	CASE ( np_MET ) !* metric term
[13295]	373	DO_2D( 1, 0, 1, 0 )
[12377]	374	zwz(ji,jj) = ( pv(ji+1,jj ,jk) + pv(ji,jj,jk) ) * di_e2v_2e1e2f(ji,jj) &
	375	& - ( pu(ji ,jj+1,jk) + pu(ji,jj,jk) ) * dj_e1u_2e1e2f(ji,jj)
	376	END_2D
[5836]	377	CASE ( np_CRV ) !* Coriolis + relative vorticity
[13295]	378	DO_2D( 1, 0, 1, 0 )
[12377]	379	zwz(ji,jj) = ff_f(ji,jj) + ( e2v(ji+1,jj) * pv(ji+1,jj,jk) - e2v(ji,jj) * pv(ji,jj,jk) &
	380	& - e1u(ji,jj+1) * pu(ji,jj+1,jk) + e1u(ji,jj) * pu(ji,jj,jk) ) * r1_e1e2f(ji,jj)
	381	END_2D
[5836]	382	CASE ( np_CME ) !* Coriolis + metric
[13295]	383	DO_2D( 1, 0, 1, 0 )
[12377]	384	zwz(ji,jj) = ff_f(ji,jj) + ( pv(ji+1,jj ,jk) + pv(ji,jj,jk) ) * di_e2v_2e1e2f(ji,jj) &
	385	& - ( pu(ji ,jj+1,jk) + pu(ji,jj,jk) ) * dj_e1u_2e1e2f(ji,jj)
	386	END_2D
[5836]	387	CASE DEFAULT ! error
	388	CALL ctl_stop('STOP','dyn_vor: wrong value for kvor' )
[455]	389	END SELECT
[5836]	390	!
	391	IF( ln_dynvor_msk ) THEN !== mask/unmask vorticity ==!
[13295]	392	DO_2D( 1, 0, 1, 0 )
[12377]	393	zwz(ji,jj) = zwz(ji,jj) * fmask(ji,jj,jk)
	394	END_2D
[5836]	395	ENDIF
[455]	396
	397	IF( ln_sco ) THEN
[12377]	398	zwz(:,:) = zwz(:,:) / e3f(:,:,jk)
	399	zwx(:,:) = e2u(:,:) * e3u(:,:,jk,Kmm) * pu(:,:,jk)
	400	zwy(:,:) = e1v(:,:) * e3v(:,:,jk,Kmm) * pv(:,:,jk)
[3]	401	ELSE
[12377]	402	zwx(:,:) = e2u(:,:) * pu(:,:,jk)
	403	zwy(:,:) = e1v(:,:) * pv(:,:,jk)
[3]	404	ENDIF
[5836]	405	! !== compute and add the vorticity term trend =!
[13295]	406	DO_2D( 0, 0, 0, 0 )
[12377]	407	zy1 = zwy(ji,jj-1) + zwy(ji+1,jj-1)
	408	zy2 = zwy(ji,jj ) + zwy(ji+1,jj )
	409	zx1 = zwx(ji-1,jj) + zwx(ji-1,jj+1)
	410	zx2 = zwx(ji ,jj) + zwx(ji ,jj+1)
	411	pu_rhs(ji,jj,jk) = pu_rhs(ji,jj,jk) + r1_4 * r1_e1u(ji,jj) * ( zwz(ji ,jj-1) * zy1 + zwz(ji,jj) * zy2 )
	412	pv_rhs(ji,jj,jk) = pv_rhs(ji,jj,jk) - r1_4 * r1_e2v(ji,jj) * ( zwz(ji-1,jj ) * zx1 + zwz(ji,jj) * zx2 )
	413	END_2D
[3]	414	! ! ===============
	415	END DO ! End of slab
	416	! ! ===============
[455]	417	END SUBROUTINE vor_ene
[216]	418
	419
[12377]	420	SUBROUTINE vor_ens( kt, Kmm, kvor, pu, pv, pu_rhs, pv_rhs )
[3]	421	!!----------------------------------------------------------------------
[455]	422	!! * ROUTINE vor_ens *
[3]	423	!!
	424	!! ** Purpose : Compute the now total vorticity trend and add it to
	425	!! the general trend of the momentum equation.
	426	!!
	427	!! ** Method : Trend evaluated using now fields (centered in time)
	428	!! and the Sadourny (1975) flux FORM formulation : conserves the
	429	!! potential enstrophy of a horizontally non-divergent flow. the
	430	!! trend of the vorticity term is given by:
[12377]	431	!! voru = 1/e1u mj-1[ (rvor+f)/e3f ] mj-1[ mi(e1v*e3v pvv(:,:,:,Kmm)) ]
	432	!! vorv = 1/e2v mi-1[ (rvor+f)/e3f ] mi-1[ mj(e2u*e3u puu(:,:,:,Kmm)) ]
	433	!! Add this trend to the general momentum trend:
	434	!! (u(rhs),v(Krhs)) = (u(rhs),v(Krhs)) + ( voru , vorv )
[3]	435	!!
[12377]	436	!! ** Action : - Update (pu_rhs,pv_rhs)) arrays with the now vorticity term trend
[3]	437	!!
[503]	438	!! References : Sadourny, r., 1975, j. atmos. sciences, 32, 680-689.
[3]	439	!!----------------------------------------------------------------------
[9019]	440	INTEGER , INTENT(in ) :: kt ! ocean time-step index
[12377]	441	INTEGER , INTENT(in ) :: Kmm ! ocean time level index
[9019]	442	INTEGER , INTENT(in ) :: kvor ! total, planetary, relative, or metric
[12377]	443	REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) :: pu, pv ! now velocities
	444	REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) :: pu_rhs, pv_rhs ! total v-trend
[2715]	445	!
[5836]	446	INTEGER :: ji, jj, jk ! dummy loop indices
	447	REAL(wp) :: zuav, zvau ! local scalars
[9019]	448	REAL(wp), DIMENSION(jpi,jpj) :: zwx, zwy, zwz, zww ! 2D workspace
[3]	449	!!----------------------------------------------------------------------
[3294]	450	!
[52]	451	IF( kt == nit000 ) THEN
	452	IF(lwp) WRITE(numout,*)
[455]	453	IF(lwp) WRITE(numout,*) 'dyn:vor_ens : vorticity term: enstrophy conserving scheme'
	454	IF(lwp) WRITE(numout,*) '~~~~~~~~~~~'
[52]	455	ENDIF
[3]	456	! ! ===============
	457	DO jk = 1, jpkm1 ! Horizontal slab
	458	! ! ===============
[1438]	459	!
[5836]	460	SELECT CASE( kvor ) !== vorticity considered ==!
	461	CASE ( np_COR ) !* Coriolis (planetary vorticity)
[7646]	462	zwz(:,:) = ff_f(:,:)
[5836]	463	CASE ( np_RVO ) !* relative vorticity
[13295]	464	DO_2D( 1, 0, 1, 0 )
[12377]	465	zwz(ji,jj) = ( e2v(ji+1,jj ) * pv(ji+1,jj ,jk) - e2v(ji,jj) * pv(ji,jj,jk) &
	466	& - e1u(ji ,jj+1) * pu(ji ,jj+1,jk) + e1u(ji,jj) * pu(ji,jj,jk) ) * r1_e1e2f(ji,jj)
	467	END_2D
[5836]	468	CASE ( np_MET ) !* metric term
[13295]	469	DO_2D( 1, 0, 1, 0 )
[12377]	470	zwz(ji,jj) = ( pv(ji+1,jj ,jk) + pv(ji,jj,jk) ) * di_e2v_2e1e2f(ji,jj) &
	471	& - ( pu(ji ,jj+1,jk) + pu(ji,jj,jk) ) * dj_e1u_2e1e2f(ji,jj)
	472	END_2D
[5836]	473	CASE ( np_CRV ) !* Coriolis + relative vorticity
[13295]	474	DO_2D( 1, 0, 1, 0 )
[12377]	475	zwz(ji,jj) = ff_f(ji,jj) + ( e2v(ji+1,jj ) * pv(ji+1,jj ,jk) - e2v(ji,jj) * pv(ji,jj,jk) &
	476	& - e1u(ji ,jj+1) * pu(ji ,jj+1,jk) + e1u(ji,jj) * pu(ji,jj,jk) ) * r1_e1e2f(ji,jj)
	477	END_2D
[5836]	478	CASE ( np_CME ) !* Coriolis + metric
[13295]	479	DO_2D( 1, 0, 1, 0 )
[12377]	480	zwz(ji,jj) = ff_f(ji,jj) + ( pv(ji+1,jj ,jk) + pv(ji,jj,jk) ) * di_e2v_2e1e2f(ji,jj) &
	481	& - ( pu(ji ,jj+1,jk) + pu(ji,jj,jk) ) * dj_e1u_2e1e2f(ji,jj)
	482	END_2D
[5836]	483	CASE DEFAULT ! error
	484	CALL ctl_stop('STOP','dyn_vor: wrong value for kvor' )
[455]	485	END SELECT
[1438]	486	!
[5836]	487	IF( ln_dynvor_msk ) THEN !== mask/unmask vorticity ==!
[13295]	488	DO_2D( 1, 0, 1, 0 )
[12377]	489	zwz(ji,jj) = zwz(ji,jj) * fmask(ji,jj,jk)
	490	END_2D
[5836]	491	ENDIF
	492	!
	493	IF( ln_sco ) THEN !== horizontal fluxes ==!
[12377]	494	zwz(:,:) = zwz(:,:) / e3f(:,:,jk)
	495	zwx(:,:) = e2u(:,:) * e3u(:,:,jk,Kmm) * pu(:,:,jk)
	496	zwy(:,:) = e1v(:,:) * e3v(:,:,jk,Kmm) * pv(:,:,jk)
[3]	497	ELSE
[12377]	498	zwx(:,:) = e2u(:,:) * pu(:,:,jk)
	499	zwy(:,:) = e1v(:,:) * pv(:,:,jk)
[3]	500	ENDIF
[5836]	501	! !== compute and add the vorticity term trend =!
[13295]	502	DO_2D( 0, 0, 0, 0 )
[12377]	503	zuav = r1_8 * r1_e1u(ji,jj) * ( zwy(ji ,jj-1) + zwy(ji+1,jj-1) &
	504	& + zwy(ji ,jj ) + zwy(ji+1,jj ) )
	505	zvau =-r1_8 * r1_e2v(ji,jj) * ( zwx(ji-1,jj ) + zwx(ji-1,jj+1) &
	506	& + zwx(ji ,jj ) + zwx(ji ,jj+1) )
	507	pu_rhs(ji,jj,jk) = pu_rhs(ji,jj,jk) + zuav * ( zwz(ji ,jj-1) + zwz(ji,jj) )
	508	pv_rhs(ji,jj,jk) = pv_rhs(ji,jj,jk) + zvau * ( zwz(ji-1,jj ) + zwz(ji,jj) )
	509	END_2D
[3]	510	! ! ===============
	511	END DO ! End of slab
	512	! ! ===============
[455]	513	END SUBROUTINE vor_ens
[216]	514
	515
[12377]	516	SUBROUTINE vor_een( kt, Kmm, kvor, pu, pv, pu_rhs, pv_rhs )
[108]	517	!!----------------------------------------------------------------------
[455]	518	!! * ROUTINE vor_een *
[108]	519	!!
	520	!! ** Purpose : Compute the now total vorticity trend and add it to
	521	!! the general trend of the momentum equation.
	522	!!
	523	!! ** Method : Trend evaluated using now fields (centered in time)
[1438]	524	!! and the Arakawa and Lamb (1980) flux form formulation : conserves
[108]	525	!! both the horizontal kinetic energy and the potential enstrophy
[1438]	526	!! when horizontal divergence is zero (see the NEMO documentation)
[12377]	527	!! Add this trend to the general momentum trend (pu_rhs,pv_rhs).
[108]	528	!!
[12377]	529	!! ** Action : - Update (pu_rhs,pv_rhs) with the now vorticity term trend
[108]	530	!!
[503]	531	!! References : Arakawa and Lamb 1980, Mon. Wea. Rev., 109, 18-36
	532	!!----------------------------------------------------------------------
[9019]	533	INTEGER , INTENT(in ) :: kt ! ocean time-step index
[12377]	534	INTEGER , INTENT(in ) :: Kmm ! ocean time level index
[9019]	535	INTEGER , INTENT(in ) :: kvor ! total, planetary, relative, or metric
[12377]	536	REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) :: pu, pv ! now velocities
	537	REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) :: pu_rhs, pv_rhs ! total v-trend
[5836]	538	!
	539	INTEGER :: ji, jj, jk ! dummy loop indices
	540	INTEGER :: ierr ! local integer
	541	REAL(wp) :: zua, zva ! local scalars
[9528]	542	REAL(wp) :: zmsk, ze3f ! local scalars
[13571]	543	REAL(wp), DIMENSION(jpi,jpj) :: zwx , zwy , z1_e3f
	544	REAL(wp), DIMENSION(jpi,jpj) :: ztnw, ztne, ztsw, ztse
	545	REAL(wp), DIMENSION(jpi,jpj,jpkm1) :: zwz ! 3D workspace, jpkm1 -> jpkm1 -> avoid lbc_lnk on jpk that is not defined
[108]	546	!!----------------------------------------------------------------------
[3294]	547	!
[108]	548	IF( kt == nit000 ) THEN
	549	IF(lwp) WRITE(numout,*)
[455]	550	IF(lwp) WRITE(numout,*) 'dyn:vor_een : vorticity term: energy and enstrophy conserving scheme'
	551	IF(lwp) WRITE(numout,*) '~~~~~~~~~~~'
[1438]	552	ENDIF
[5836]	553	!
	554	! ! ===============
	555	DO jk = 1, jpkm1 ! Horizontal slab
	556	! ! ===============
	557	!
	558	SELECT CASE( nn_een_e3f ) ! == reciprocal of e3 at F-point
	559	CASE ( 0 ) ! original formulation (masked averaging of e3t divided by 4)
[13295]	560	DO_2D( 1, 0, 1, 0 )
[13237]	561	ze3f = ( e3t(ji ,jj+1,jk,Kmm)*tmask(ji ,jj+1,jk) &
	562	& + e3t(ji+1,jj+1,jk,Kmm)*tmask(ji+1,jj+1,jk) &
	563	& + e3t(ji ,jj ,jk,Kmm)*tmask(ji ,jj ,jk) &
	564	& + e3t(ji+1,jj ,jk,Kmm)*tmask(ji+1,jj ,jk) )
[12377]	565	IF( ze3f /= 0._wp ) THEN ; z1_e3f(ji,jj) = 4._wp / ze3f
	566	ELSE ; z1_e3f(ji,jj) = 0._wp
	567	ENDIF
	568	END_2D
[5836]	569	CASE ( 1 ) ! new formulation (masked averaging of e3t divided by the sum of mask)
[13295]	570	DO_2D( 1, 0, 1, 0 )
[13237]	571	ze3f = ( e3t(ji ,jj+1,jk,Kmm)*tmask(ji ,jj+1,jk) &
	572	& + e3t(ji+1,jj+1,jk,Kmm)*tmask(ji+1,jj+1,jk) &
	573	& + e3t(ji ,jj ,jk,Kmm)*tmask(ji ,jj ,jk) &
	574	& + e3t(ji+1,jj ,jk,Kmm)*tmask(ji+1,jj ,jk) )
[12377]	575	zmsk = ( tmask(ji,jj+1,jk) + tmask(ji+1,jj+1,jk) &
	576	& + tmask(ji,jj ,jk) + tmask(ji+1,jj ,jk) )
	577	IF( ze3f /= 0._wp ) THEN ; z1_e3f(ji,jj) = zmsk / ze3f
	578	ELSE ; z1_e3f(ji,jj) = 0._wp
	579	ENDIF
	580	END_2D
[5836]	581	END SELECT
	582	!
	583	SELECT CASE( kvor ) !== vorticity considered ==!
	584	CASE ( np_COR ) !* Coriolis (planetary vorticity)
[13295]	585	DO_2D( 1, 0, 1, 0 )
[12377]	586	zwz(ji,jj,jk) = ff_f(ji,jj) * z1_e3f(ji,jj)
	587	END_2D
[5836]	588	CASE ( np_RVO ) !* relative vorticity
[13295]	589	DO_2D( 1, 0, 1, 0 )
[12377]	590	zwz(ji,jj,jk) = ( e2v(ji+1,jj ) * pv(ji+1,jj,jk) - e2v(ji,jj) * pv(ji,jj,jk) &
	591	& - e1u(ji ,jj+1) * pu(ji,jj+1,jk) + e1u(ji,jj) * pu(ji,jj,jk) ) * r1_e1e2f(ji,jj)*z1_e3f(ji,jj)
	592	END_2D
[5836]	593	CASE ( np_MET ) !* metric term
[13295]	594	DO_2D( 1, 0, 1, 0 )
[12377]	595	zwz(ji,jj,jk) = ( ( pv(ji+1,jj,jk) + pv(ji,jj,jk) ) * di_e2v_2e1e2f(ji,jj) &
	596	& - ( pu(ji,jj+1,jk) + pu(ji,jj,jk) ) * dj_e1u_2e1e2f(ji,jj) ) * z1_e3f(ji,jj)
	597	END_2D
[5836]	598	CASE ( np_CRV ) !* Coriolis + relative vorticity
[13295]	599	DO_2D( 1, 0, 1, 0 )
[12377]	600	zwz(ji,jj,jk) = ( ff_f(ji,jj) + ( e2v(ji+1,jj ) * pv(ji+1,jj,jk) - e2v(ji,jj) * pv(ji,jj,jk) &
	601	& - e1u(ji ,jj+1) * pu(ji,jj+1,jk) + e1u(ji,jj) * pu(ji,jj,jk) ) &
	602	& * r1_e1e2f(ji,jj) ) * z1_e3f(ji,jj)
	603	END_2D
[5836]	604	CASE ( np_CME ) !* Coriolis + metric
[13295]	605	DO_2D( 1, 0, 1, 0 )
[12377]	606	zwz(ji,jj,jk) = ( ff_f(ji,jj) + ( pv(ji+1,jj ,jk) + pv(ji,jj,jk) ) * di_e2v_2e1e2f(ji,jj) &
	607	& - ( pu(ji ,jj+1,jk) + pu(ji,jj,jk) ) * dj_e1u_2e1e2f(ji,jj) ) * z1_e3f(ji,jj)
	608	END_2D
[5836]	609	CASE DEFAULT ! error
	610	CALL ctl_stop('STOP','dyn_vor: wrong value for kvor' )
[455]	611	END SELECT
[5836]	612	!
	613	IF( ln_dynvor_msk ) THEN !== mask/unmask vorticity ==!
[13295]	614	DO_2D( 1, 0, 1, 0 )
[12377]	615	zwz(ji,jj,jk) = zwz(ji,jj,jk) * fmask(ji,jj,jk)
	616	END_2D
[5836]	617	ENDIF
[10425]	618	END DO ! End of slab
[5836]	619	!
[13226]	620	CALL lbc_lnk( 'dynvor', zwz, 'F', 1.0_wp )
[10425]	621
	622	DO jk = 1, jpkm1 ! Horizontal slab
[5907]	623	!
[5836]	624	! !== horizontal fluxes ==!
[12377]	625	zwx(:,:) = e2u(:,:) * e3u(:,:,jk,Kmm) * pu(:,:,jk)
	626	zwy(:,:) = e1v(:,:) * e3v(:,:,jk,Kmm) * pv(:,:,jk)
[108]	627
[5836]	628	! !== compute and add the vorticity term trend =!
[1438]	629	jj = 2
	630	ztne(1,:) = 0 ; ztnw(1,:) = 0 ; ztse(1,:) = 0 ; ztsw(1,:) = 0
[5836]	631	DO ji = 2, jpi ! split in 2 parts due to vector opt.
[10425]	632	ztne(ji,jj) = zwz(ji-1,jj ,jk) + zwz(ji ,jj ,jk) + zwz(ji ,jj-1,jk)
	633	ztnw(ji,jj) = zwz(ji-1,jj-1,jk) + zwz(ji-1,jj ,jk) + zwz(ji ,jj ,jk)
	634	ztse(ji,jj) = zwz(ji ,jj ,jk) + zwz(ji ,jj-1,jk) + zwz(ji-1,jj-1,jk)
	635	ztsw(ji,jj) = zwz(ji ,jj-1,jk) + zwz(ji-1,jj-1,jk) + zwz(ji-1,jj ,jk)
[108]	636	END DO
	637	DO jj = 3, jpj
[12377]	638	DO ji = 2, jpi ! vector opt. ok because we start at jj = 3
[10425]	639	ztne(ji,jj) = zwz(ji-1,jj ,jk) + zwz(ji ,jj ,jk) + zwz(ji ,jj-1,jk)
	640	ztnw(ji,jj) = zwz(ji-1,jj-1,jk) + zwz(ji-1,jj ,jk) + zwz(ji ,jj ,jk)
	641	ztse(ji,jj) = zwz(ji ,jj ,jk) + zwz(ji ,jj-1,jk) + zwz(ji-1,jj-1,jk)
	642	ztsw(ji,jj) = zwz(ji ,jj-1,jk) + zwz(ji-1,jj-1,jk) + zwz(ji-1,jj ,jk)
[108]	643	END DO
	644	END DO
[13295]	645	DO_2D( 0, 0, 0, 0 )
[12377]	646	zua = + r1_12 * r1_e1u(ji,jj) * ( ztne(ji,jj ) * zwy(ji ,jj ) + ztnw(ji+1,jj) * zwy(ji+1,jj ) &
	647	& + ztse(ji,jj ) * zwy(ji ,jj-1) + ztsw(ji+1,jj) * zwy(ji+1,jj-1) )
	648	zva = - r1_12 * r1_e2v(ji,jj) * ( ztsw(ji,jj+1) * zwx(ji-1,jj+1) + ztse(ji,jj+1) * zwx(ji ,jj+1) &
	649	& + ztnw(ji,jj ) * zwx(ji-1,jj ) + ztne(ji,jj ) * zwx(ji ,jj ) )
	650	pu_rhs(ji,jj,jk) = pu_rhs(ji,jj,jk) + zua
	651	pv_rhs(ji,jj,jk) = pv_rhs(ji,jj,jk) + zva
	652	END_2D
[108]	653	! ! ===============
	654	END DO ! End of slab
	655	! ! ===============
[455]	656	END SUBROUTINE vor_een
[216]	657
	658
[9528]	659
[12377]	660	SUBROUTINE vor_eeT( kt, Kmm, kvor, pu, pv, pu_rhs, pv_rhs )
[9528]	661	!!----------------------------------------------------------------------
	662	!! * ROUTINE vor_eeT *
	663	!!
	664	!! ** Purpose : Compute the now total vorticity trend and add it to
	665	!! the general trend of the momentum equation.
	666	!!
	667	!! ** Method : Trend evaluated using now fields (centered in time)
	668	!! and the Arakawa and Lamb (1980) vector form formulation using
	669	!! a modified version of Arakawa and Lamb (1980) scheme (see vor_een).
	670	!! The change consists in
[12377]	671	!! Add this trend to the general momentum trend (pu_rhs,pv_rhs).
[9528]	672	!!
[12377]	673	!! ** Action : - Update (pu_rhs,pv_rhs) with the now vorticity term trend
[9528]	674	!!
	675	!! References : Arakawa and Lamb 1980, Mon. Wea. Rev., 109, 18-36
	676	!!----------------------------------------------------------------------
	677	INTEGER , INTENT(in ) :: kt ! ocean time-step index
[12377]	678	INTEGER , INTENT(in ) :: Kmm ! ocean time level index
[9528]	679	INTEGER , INTENT(in ) :: kvor ! total, planetary, relative, or metric
[12377]	680	REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) :: pu, pv ! now velocities
	681	REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) :: pu_rhs, pv_rhs ! total v-trend
[9528]	682	!
	683	INTEGER :: ji, jj, jk ! dummy loop indices
	684	INTEGER :: ierr ! local integer
	685	REAL(wp) :: zua, zva ! local scalars
	686	REAL(wp) :: zmsk, z1_e3t ! local scalars
[13571]	687	REAL(wp), DIMENSION(jpi,jpj) :: zwx , zwy
	688	REAL(wp), DIMENSION(jpi,jpj) :: ztnw, ztne, ztsw, ztse
	689	REAL(wp), DIMENSION(jpi,jpj,jpkm1) :: zwz ! 3D workspace, avoid lbc_lnk on jpk that is not defined
[9528]	690	!!----------------------------------------------------------------------
	691	!
	692	IF( kt == nit000 ) THEN
	693	IF(lwp) WRITE(numout,*)
	694	IF(lwp) WRITE(numout,*) 'dyn:vor_een : vorticity term: energy and enstrophy conserving scheme'
	695	IF(lwp) WRITE(numout,*) '~~~~~~~~~~~'
	696	ENDIF
	697	!
	698	! ! ===============
	699	DO jk = 1, jpkm1 ! Horizontal slab
	700	! ! ===============
	701	!
	702	!
	703	SELECT CASE( kvor ) !== vorticity considered ==!
	704	CASE ( np_COR ) !* Coriolis (planetary vorticity)
[13295]	705	DO_2D( 1, 0, 1, 0 )
[12377]	706	zwz(ji,jj,jk) = ff_f(ji,jj)
	707	END_2D
[9528]	708	CASE ( np_RVO ) !* relative vorticity
[13295]	709	DO_2D( 1, 0, 1, 0 )
[12377]	710	zwz(ji,jj,jk) = ( e2v(ji+1,jj ) * pv(ji+1,jj ,jk) - e2v(ji,jj) * pv(ji,jj,jk) &
	711	& - e1u(ji ,jj+1) * pu(ji ,jj+1,jk) + e1u(ji,jj) * pu(ji,jj,jk) ) &
	712	& * r1_e1e2f(ji,jj)
	713	END_2D
[9528]	714	CASE ( np_MET ) !* metric term
[13295]	715	DO_2D( 1, 0, 1, 0 )
[12377]	716	zwz(ji,jj,jk) = ( pv(ji+1,jj ,jk) + pv(ji,jj,jk) ) * di_e2v_2e1e2f(ji,jj) &
	717	& - ( pu(ji ,jj+1,jk) + pu(ji,jj,jk) ) * dj_e1u_2e1e2f(ji,jj)
	718	END_2D
[9528]	719	CASE ( np_CRV ) !* Coriolis + relative vorticity
[13295]	720	DO_2D( 1, 0, 1, 0 )
[12377]	721	zwz(ji,jj,jk) = ( ff_f(ji,jj) + ( e2v(ji+1,jj ) * pv(ji+1,jj ,jk) - e2v(ji,jj) * pv(ji,jj,jk) &
	722	& - e1u(ji ,jj+1) * pu(ji ,jj+1,jk) + e1u(ji,jj) * pu(ji,jj,jk) ) &
	723	& * r1_e1e2f(ji,jj) )
	724	END_2D
[9528]	725	CASE ( np_CME ) !* Coriolis + metric
[13295]	726	DO_2D( 1, 0, 1, 0 )
[12377]	727	zwz(ji,jj,jk) = ff_f(ji,jj) + ( pv(ji+1,jj ,jk) + pv(ji,jj,jk) ) * di_e2v_2e1e2f(ji,jj) &
	728	& - ( pu(ji ,jj+1,jk) + pu(ji,jj,jk) ) * dj_e1u_2e1e2f(ji,jj)
	729	END_2D
[9528]	730	CASE DEFAULT ! error
	731	CALL ctl_stop('STOP','dyn_vor: wrong value for kvor' )
	732	END SELECT
	733	!
	734	IF( ln_dynvor_msk ) THEN !== mask/unmask vorticity ==!
[13295]	735	DO_2D( 1, 0, 1, 0 )
[12377]	736	zwz(ji,jj,jk) = zwz(ji,jj,jk) * fmask(ji,jj,jk)
	737	END_2D
[9528]	738	ENDIF
[10425]	739	END DO
	740	!
[13226]	741	CALL lbc_lnk( 'dynvor', zwz, 'F', 1.0_wp )
[10425]	742	!
	743	DO jk = 1, jpkm1 ! Horizontal slab
	744
	745	! !== horizontal fluxes ==!
[12377]	746	zwx(:,:) = e2u(:,:) * e3u(:,:,jk,Kmm) * pu(:,:,jk)
	747	zwy(:,:) = e1v(:,:) * e3v(:,:,jk,Kmm) * pv(:,:,jk)
[9528]	748
	749	! !== compute and add the vorticity term trend =!
	750	jj = 2
	751	ztne(1,:) = 0 ; ztnw(1,:) = 0 ; ztse(1,:) = 0 ; ztsw(1,:) = 0
	752	DO ji = 2, jpi ! split in 2 parts due to vector opt.
[12377]	753	z1_e3t = 1._wp / e3t(ji,jj,jk,Kmm)
[10425]	754	ztne(ji,jj) = ( zwz(ji-1,jj ,jk) + zwz(ji ,jj ,jk) + zwz(ji ,jj-1,jk) ) * z1_e3t
	755	ztnw(ji,jj) = ( zwz(ji-1,jj-1,jk) + zwz(ji-1,jj ,jk) + zwz(ji ,jj ,jk) ) * z1_e3t
	756	ztse(ji,jj) = ( zwz(ji ,jj ,jk) + zwz(ji ,jj-1,jk) + zwz(ji-1,jj-1,jk) ) * z1_e3t
	757	ztsw(ji,jj) = ( zwz(ji ,jj-1,jk) + zwz(ji-1,jj-1,jk) + zwz(ji-1,jj ,jk) ) * z1_e3t
[9528]	758	END DO
	759	DO jj = 3, jpj
[12377]	760	DO ji = 2, jpi ! vector opt. ok because we start at jj = 3
	761	z1_e3t = 1._wp / e3t(ji,jj,jk,Kmm)
[10425]	762	ztne(ji,jj) = ( zwz(ji-1,jj ,jk) + zwz(ji ,jj ,jk) + zwz(ji ,jj-1,jk) ) * z1_e3t
	763	ztnw(ji,jj) = ( zwz(ji-1,jj-1,jk) + zwz(ji-1,jj ,jk) + zwz(ji ,jj ,jk) ) * z1_e3t
	764	ztse(ji,jj) = ( zwz(ji ,jj ,jk) + zwz(ji ,jj-1,jk) + zwz(ji-1,jj-1,jk) ) * z1_e3t
	765	ztsw(ji,jj) = ( zwz(ji ,jj-1,jk) + zwz(ji-1,jj-1,jk) + zwz(ji-1,jj ,jk) ) * z1_e3t
[9528]	766	END DO
	767	END DO
[13295]	768	DO_2D( 0, 0, 0, 0 )
[12377]	769	zua = + r1_12 * r1_e1u(ji,jj) * ( ztne(ji,jj ) * zwy(ji ,jj ) + ztnw(ji+1,jj) * zwy(ji+1,jj ) &
	770	& + ztse(ji,jj ) * zwy(ji ,jj-1) + ztsw(ji+1,jj) * zwy(ji+1,jj-1) )
	771	zva = - r1_12 * r1_e2v(ji,jj) * ( ztsw(ji,jj+1) * zwx(ji-1,jj+1) + ztse(ji,jj+1) * zwx(ji ,jj+1) &
	772	& + ztnw(ji,jj ) * zwx(ji-1,jj ) + ztne(ji,jj ) * zwx(ji ,jj ) )
	773	pu_rhs(ji,jj,jk) = pu_rhs(ji,jj,jk) + zua
	774	pv_rhs(ji,jj,jk) = pv_rhs(ji,jj,jk) + zva
	775	END_2D
[9528]	776	! ! ===============
	777	END DO ! End of slab
	778	! ! ===============
	779	END SUBROUTINE vor_eeT
	780
	781
[2528]	782	SUBROUTINE dyn_vor_init
[3]	783	!!---------------------------------------------------------------------
[2528]	784	!! * ROUTINE dyn_vor_init *
[3]	785	!!
	786	!! ** Purpose : Control the consistency between cpp options for
[1438]	787	!! tracer advection schemes
[3]	788	!!----------------------------------------------------------------------
[9528]	789	INTEGER :: ji, jj, jk ! dummy loop indices
	790	INTEGER :: ioptio, ios ! local integer
[2715]	791	!!
[9528]	792	NAMELIST/namdyn_vor/ ln_dynvor_ens, ln_dynvor_ene, ln_dynvor_enT, ln_dynvor_eeT, &
	793	& ln_dynvor_een, nn_een_e3f , ln_dynvor_mix, ln_dynvor_msk
[3]	794	!!----------------------------------------------------------------------
[9528]	795	!
	796	IF(lwp) THEN
	797	WRITE(numout,*)
	798	WRITE(numout,*) 'dyn_vor_init : vorticity term : read namelist and control the consistency'
	799	WRITE(numout,*) '~~~~~~~~~~~~'
	800	ENDIF
	801	!
[4147]	802	READ ( numnam_ref, namdyn_vor, IOSTAT = ios, ERR = 901)
[11536]	803	901 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namdyn_vor in reference namelist' )
[4147]	804	READ ( numnam_cfg, namdyn_vor, IOSTAT = ios, ERR = 902 )
[11536]	805	902 IF( ios > 0 ) CALL ctl_nam ( ios , 'namdyn_vor in configuration namelist' )
[4624]	806	IF(lwm) WRITE ( numond, namdyn_vor )
[9528]	807	!
[503]	808	IF(lwp) THEN ! Namelist print
[7646]	809	WRITE(numout,*) ' Namelist namdyn_vor : choice of the vorticity term scheme'
	810	WRITE(numout,*) ' enstrophy conserving scheme ln_dynvor_ens = ', ln_dynvor_ens
[9528]	811	WRITE(numout,*) ' f-point energy conserving scheme ln_dynvor_ene = ', ln_dynvor_ene
	812	WRITE(numout,*) ' t-point energy conserving scheme ln_dynvor_enT = ', ln_dynvor_enT
	813	WRITE(numout,*) ' energy conserving scheme (een using e3t) ln_dynvor_eeT = ', ln_dynvor_eeT
[7646]	814	WRITE(numout,*) ' enstrophy and energy conserving scheme ln_dynvor_een = ', ln_dynvor_een
	815	WRITE(numout,*) ' e3f = averaging /4 (=0) or /sum(tmask) (=1) nn_een_e3f = ', nn_een_e3f
[9528]	816	WRITE(numout,*) ' mixed enstrophy/energy conserving scheme ln_dynvor_mix = ', ln_dynvor_mix
[7646]	817	WRITE(numout,*) ' masked (=T) or unmasked(=F) vorticity ln_dynvor_msk = ', ln_dynvor_msk
[52]	818	ENDIF
	819
[9528]	820	IF( ln_dynvor_msk ) CALL ctl_stop( 'dyn_vor_init: masked vorticity is not currently not available')
	821
[5836]	822	!!gm this should be removed when choosing a unique strategy for fmask at the coast
[3294]	823	! If energy, enstrophy or mixed advection of momentum in vector form change the value for masks
	824	! at angles with three ocean points and one land point
[5836]	825	IF(lwp) WRITE(numout,*)
[7646]	826	IF(lwp) WRITE(numout,*) ' change fmask value in the angles (T) ln_vorlat = ', ln_vorlat
[3294]	827	IF( ln_vorlat .AND. ( ln_dynvor_ene .OR. ln_dynvor_ens .OR. ln_dynvor_mix ) ) THEN
[13295]	828	DO_3D( 1, 0, 1, 0, 1, jpk )
[12377]	829	IF( tmask(ji,jj+1,jk) + tmask(ji+1,jj+1,jk) &
[12793]	830	& + tmask(ji,jj ,jk) + tmask(ji+1,jj ,jk) == 3._wp ) fmask(ji,jj,jk) = 1._wp
[12377]	831	END_3D
[9528]	832	!
[13630]	833	#if defined key_mpi3
	834	CALL lbc_lnk_nc_multi( 'dynvor', fmask, 'F', 1._wp ) ! Lateral boundary conditions on fmask
	835	#else
[10425]	836	CALL lbc_lnk( 'dynvor', fmask, 'F', 1._wp ) ! Lateral boundary conditions on fmask
[13630]	837	#endif
[9528]	838	!
[3294]	839	ENDIF
[5836]	840	!!gm end
[3294]	841
[5836]	842	ioptio = 0 ! type of scheme for vorticity (set nvor_scheme)
[9528]	843	IF( ln_dynvor_ens ) THEN ; ioptio = ioptio + 1 ; nvor_scheme = np_ENS ; ENDIF
	844	IF( ln_dynvor_ene ) THEN ; ioptio = ioptio + 1 ; nvor_scheme = np_ENE ; ENDIF
	845	IF( ln_dynvor_enT ) THEN ; ioptio = ioptio + 1 ; nvor_scheme = np_ENT ; ENDIF
	846	IF( ln_dynvor_eeT ) THEN ; ioptio = ioptio + 1 ; nvor_scheme = np_EET ; ENDIF
	847	IF( ln_dynvor_een ) THEN ; ioptio = ioptio + 1 ; nvor_scheme = np_EEN ; ENDIF
	848	IF( ln_dynvor_mix ) THEN ; ioptio = ioptio + 1 ; nvor_scheme = np_MIX ; ENDIF
[5836]	849	!
[6140]	850	IF( ioptio /= 1 ) CALL ctl_stop( ' use ONE and ONLY one vorticity scheme' )
[5836]	851	!
	852	IF(lwp) WRITE(numout,*) ! type of calculated vorticity (set ncor, nrvm, ntot)
[9019]	853	ncor = np_COR ! planetary vorticity
	854	SELECT CASE( n_dynadv )
	855	CASE( np_LIN_dyn )
[9190]	856	IF(lwp) WRITE(numout,*) ' ==>>> linear dynamics : total vorticity = Coriolis'
[9019]	857	nrvm = np_COR ! planetary vorticity
	858	ntot = np_COR ! - -
	859	CASE( np_VEC_c2 )
[9190]	860	IF(lwp) WRITE(numout,*) ' ==>>> vector form dynamics : total vorticity = Coriolis + relative vorticity'
[5836]	861	nrvm = np_RVO ! relative vorticity
[9019]	862	ntot = np_CRV ! relative + planetary vorticity
	863	CASE( np_FLX_c2 , np_FLX_ubs )
[9190]	864	IF(lwp) WRITE(numout,*) ' ==>>> flux form dynamics : total vorticity = Coriolis + metric term'
[5836]	865	nrvm = np_MET ! metric term
	866	ntot = np_CME ! Coriolis + metric term
[9528]	867	!
	868	SELECT CASE( nvor_scheme ) ! pre-computed gradients for the metric term:
	869	CASE( np_ENT ) !* T-point metric term : pre-compute di(e2u)/2 and dj(e1v)/2
	870	ALLOCATE( di_e2u_2(jpi,jpj), dj_e1v_2(jpi,jpj) )
[13295]	871	DO_2D( 0, 0, 0, 0 )
[12377]	872	di_e2u_2(ji,jj) = ( e2u(ji,jj) - e2u(ji-1,jj ) ) * 0.5_wp
	873	dj_e1v_2(ji,jj) = ( e1v(ji,jj) - e1v(ji ,jj-1) ) * 0.5_wp
	874	END_2D
[13630]	875	#if defined key_mpi3
	876	CALL lbc_lnk_nc_multi( 'dynvor', di_e2u_2, 'T', -1.0_wp , dj_e1v_2, 'T', -1.0_wp ) ! Lateral boundary conditions
	877	#else
[13226]	878	CALL lbc_lnk_multi( 'dynvor', di_e2u_2, 'T', -1.0_wp , dj_e1v_2, 'T', -1.0_wp ) ! Lateral boundary conditions
[13630]	879	#endif
[9528]	880	!
	881	CASE DEFAULT !* F-point metric term : pre-compute di(e2u)/(2e1e2f) and dj(e1v)/(2e1e2f)
	882	ALLOCATE( di_e2v_2e1e2f(jpi,jpj), dj_e1u_2e1e2f(jpi,jpj) )
[13295]	883	DO_2D( 1, 0, 1, 0 )
[12377]	884	di_e2v_2e1e2f(ji,jj) = ( e2v(ji+1,jj ) - e2v(ji,jj) ) * 0.5 * r1_e1e2f(ji,jj)
	885	dj_e1u_2e1e2f(ji,jj) = ( e1u(ji ,jj+1) - e1u(ji,jj) ) * 0.5 * r1_e1e2f(ji,jj)
	886	END_2D
[13630]	887	#if defined key_mpi3
	888	CALL lbc_lnk_nc_multi( 'dynvor', di_e2v_2e1e2f, 'F', -1.0_wp , dj_e1u_2e1e2f, 'F', -1.0_wp ) ! Lateral boundary conditions
	889	#else
[13226]	890	CALL lbc_lnk_multi( 'dynvor', di_e2v_2e1e2f, 'F', -1.0_wp , dj_e1u_2e1e2f, 'F', -1.0_wp ) ! Lateral boundary conditions
[13630]	891	#endif
[9528]	892	END SELECT
	893	!
[9019]	894	END SELECT
[643]	895
[503]	896	IF(lwp) THEN ! Print the choice
	897	WRITE(numout,*)
[9019]	898	SELECT CASE( nvor_scheme )
[9528]	899	CASE( np_ENS ) ; WRITE(numout,*) ' ==>>> enstrophy conserving scheme (ENS)'
	900	CASE( np_ENE ) ; WRITE(numout,*) ' ==>>> energy conserving scheme (Coriolis at F-points) (ENE)'
	901	CASE( np_ENT ) ; WRITE(numout,*) ' ==>>> energy conserving scheme (Coriolis at T-points) (ENT)'
	902	CASE( np_EET ) ; WRITE(numout,*) ' ==>>> energy conserving scheme (EEN scheme using e3t) (EET)'
	903	CASE( np_EEN ) ; WRITE(numout,*) ' ==>>> energy and enstrophy conserving scheme (EEN)'
	904	CASE( np_MIX ) ; WRITE(numout,*) ' ==>>> mixed enstrophy/energy conserving scheme (MIX)'
[9019]	905	END SELECT
[3]	906	ENDIF
[503]	907	!
[2528]	908	END SUBROUTINE dyn_vor_init
[3]	909
[503]	910	!!==============================================================================
[3]	911	END MODULE dynvor

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

Download in other formats: