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 branches/2011/dev_LOCEAN_CMCC_2011/NEMOGCM/NEMO/OPA_SRC/DYN – NEMO

source: branches/2011/dev_LOCEAN_CMCC_2011/NEMOGCM/NEMO/OPA_SRC/DYN/dynvor.F90 @ 3084

Last change on this file since 3084 was 3084, checked in by cetlod, 12 years ago

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