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.
dynnxt.F90 in NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/DYN – NEMO

source: NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/DYN/dynnxt.F90 @ 10743

Last change on this file since 10743 was 10743, checked in by davestorkey, 5 years ago

First block of changes:

  1. New state variables uu, vv, e3X and corresponding temporary pointers.
  2. Updated code for time-filtering and time-level-swapping uu, vv, e3X.
  3. DYN routines not updated to use new variables at this stage.
  • Property svn:keywords set to Id
File size: 19.9 KB
Line 
1MODULE dynnxt
2   !!=========================================================================
3   !!                       ***  MODULE  dynnxt  ***
4   !! Ocean dynamics: time stepping
5   !!=========================================================================
6   !! History :  OPA  !  1987-02  (P. Andrich, D. L Hostis)  Original code
7   !!                 !  1990-10  (C. Levy, G. Madec)
8   !!            7.0  !  1993-03  (M. Guyon)  symetrical conditions
9   !!            8.0  !  1997-02  (G. Madec & M. Imbard)  opa, release 8.0
10   !!            8.2  !  1997-04  (A. Weaver)  Euler forward step
11   !!             -   !  1997-06  (G. Madec)  lateral boudary cond., lbc routine
12   !!    NEMO    1.0  !  2002-08  (G. Madec)  F90: Free form and module
13   !!             -   !  2002-10  (C. Talandier, A-M. Treguier) Open boundary cond.
14   !!            2.0  !  2005-11  (V. Garnier) Surface pressure gradient organization
15   !!            2.3  !  2007-07  (D. Storkey) Calls to BDY routines.
16   !!            3.2  !  2009-06  (G. Madec, R.Benshila)  re-introduce the vvl option
17   !!            3.3  !  2010-09  (D. Storkey, E.O'Dea) Bug fix for BDY module
18   !!            3.3  !  2011-03  (P. Oddo) Bug fix for time-splitting+(BDY-OBC) and not VVL
19   !!            3.5  !  2013-07  (J. Chanut) Compliant with time splitting changes
20   !!            3.6  !  2014-04  (G. Madec) add the diagnostic of the time filter trends
21   !!            3.7  !  2015-11  (J. Chanut) Free surface simplification
22   !!-------------------------------------------------------------------------
23 
24   !!-------------------------------------------------------------------------
25   !!   dyn_nxt       : obtain the next (after) horizontal velocity
26   !!-------------------------------------------------------------------------
27   USE oce            ! ocean dynamics and tracers
28   USE dom_oce        ! ocean space and time domain
29   USE sbc_oce        ! Surface boundary condition: ocean fields
30   USE sbcrnf         ! river runoffs
31   USE sbcisf         ! ice shelf
32   USE phycst         ! physical constants
33   USE dynadv         ! dynamics: vector invariant versus flux form
34   USE dynspg_ts      ! surface pressure gradient: split-explicit scheme
35   USE domvvl         ! variable volume
36   USE bdy_oce   , ONLY: ln_bdy
37   USE bdydta         ! ocean open boundary conditions
38   USE bdydyn         ! ocean open boundary conditions
39   USE bdyvol         ! ocean open boundary condition (bdy_vol routines)
40   USE trd_oce        ! trends: ocean variables
41   USE trddyn         ! trend manager: dynamics
42   USE trdken         ! trend manager: kinetic energy
43   !
44   USE in_out_manager ! I/O manager
45   USE iom            ! I/O manager library
46   USE lbclnk         ! lateral boundary condition (or mpp link)
47   USE lib_mpp        ! MPP library
48   USE prtctl         ! Print control
49   USE timing         ! Timing
50#if defined key_agrif
51   USE agrif_oce_interp
52#endif
53
54   IMPLICIT NONE
55   PRIVATE
56
57   PUBLIC    dyn_nxt   ! routine called by step.F90
58
59   !!----------------------------------------------------------------------
60   !! NEMO/OCE 4.0 , NEMO Consortium (2018)
61   !! $Id$
62   !! Software governed by the CeCILL license (see ./LICENSE)
63   !!----------------------------------------------------------------------
64CONTAINS
65
66   SUBROUTINE dyn_nxt ( kt, kNm1, kNp1, puu, pvv, pe3t, pe3u, pe3v )
67      !!----------------------------------------------------------------------
68      !!                  ***  ROUTINE dyn_nxt  ***
69      !!                   
70      !! ** Purpose :   Finalize after horizontal velocity. Apply the boundary
71      !!             condition on the after velocity, achieve the time stepping
72      !!             by applying the Asselin filter on now fields and swapping
73      !!             the fields.
74      !!
75      !! ** Method  : * Ensure after velocities transport matches time splitting
76      !!             estimate (ln_dynspg_ts=T)
77      !!
78      !!              * Apply lateral boundary conditions on after velocity
79      !!             at the local domain boundaries through lbc_lnk call,
80      !!             at the one-way open boundaries (ln_bdy=T),
81      !!             at the AGRIF zoom   boundaries (lk_agrif=T)
82      !!
83      !!              * Apply the time filter applied and swap of the dynamics
84      !!             arrays to start the next time step:
85      !!                (ub,vb) = (un,vn) + atfp [ (ub,vb) + (ua,va) - 2 (un,vn) ]
86      !!                (un,vn) = (ua,va).
87      !!             Note that with flux form advection and non linear free surface,
88      !!             the time filter is applied on thickness weighted velocity.
89      !!             As a result, dyn_nxt MUST be called after tra_nxt.
90      !!
91      !! ** Action :   ub,vb   filtered before horizontal velocity of next time-step
92      !!               un,vn   now horizontal velocity of next time-step
93      !!----------------------------------------------------------------------
94      INTEGER, INTENT( in ) ::   kt         ! ocean time-step index
95      INTEGER, INTENT( in ) ::   kNm1, kNp1 ! before and after time level indices
96      REAL(wp), INTENT( inout), DIMENSION(jpi,jpj,jpk) :: puu, pvv ! now velocities to be time filtered
97      REAL(wp), INTENT( inout), DIMENSION(jpi,jpj,jpk) :: pe3t, pe3u, pe3v ! now scale factors to be time filtered
98      !
99      INTEGER  ::   ji, jj, jk   ! dummy loop indices
100      INTEGER  ::   ikt          ! local integers
101      REAL(wp) ::   zue3a, zue3n, zue3b, zcoef    ! local scalars
102      REAL(wp) ::   zve3a, zve3n, zve3b, z1_2dt   !   -      -
103      REAL(wp), ALLOCATABLE, DIMENSION(:,:)   ::   zue, zve
104      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) ::   ze3t_f, ze3u_f, ze3v_f, zua, zva 
105      !!----------------------------------------------------------------------
106      !
107      IF( ln_timing    )   CALL timing_start('dyn_nxt')
108      IF( ln_dynspg_ts )   ALLOCATE( zue(jpi,jpj)     , zve(jpi,jpj)     )
109      IF( l_trddyn     )   ALLOCATE( zua(jpi,jpj,jpk) , zva(jpi,jpj,jpk) )
110      !
111      IF( kt == nit000 ) THEN
112         IF(lwp) WRITE(numout,*)
113         IF(lwp) WRITE(numout,*) 'dyn_nxt : time stepping'
114         IF(lwp) WRITE(numout,*) '~~~~~~~'
115      ENDIF
116
117      IF ( ln_dynspg_ts ) THEN
118         ! Ensure below that barotropic velocities match time splitting estimate
119         ! Compute actual transport and replace it with ts estimate at "after" time step
120         zue(:,:) = e3u(:,:,1,kNp1) * uu(:,:,1,kNp1) * umask(:,:,1)
121         zve(:,:) = e3v(:,:,1,kNp1) * vv(:,:,1,kNp1) * vmask(:,:,1)
122         DO jk = 2, jpkm1
123            zue(:,:) = zue(:,:) + e3u(:,:,jk,kNp1) * uu(:,:,jk,kNp1) * umask(:,:,jk)
124            zve(:,:) = zve(:,:) + e3v(:,:,jk,kNp1) * vv(:,:,jk,kNp1) * vmask(:,:,jk)
125         END DO
126         DO jk = 1, jpkm1
127            uu(:,:,jk,kNp1) = ( uu(:,:,jk,kNp1) - zue(:,:) * r1_hu_a(:,:) + ua_b(:,:) ) * umask(:,:,jk)
128            vv(:,:,jk,kNp1) = ( vv(:,:,jk,kNp1) - zve(:,:) * r1_hv_a(:,:) + va_b(:,:) ) * vmask(:,:,jk)
129         END DO
130         !
131         IF( .NOT.ln_bt_fw ) THEN
132            ! Remove advective velocity from "now velocities"
133            ! prior to asselin filtering     
134            ! In the forward case, this is done below after asselin filtering   
135            ! so that asselin contribution is removed at the same time
136            DO jk = 1, jpkm1
137               puu(:,:,jk) = ( puu(:,:,jk) - un_adv(:,:)*r1_hu_n(:,:) + un_b(:,:) )*umask(:,:,jk)
138               pvv(:,:,jk) = ( pvv(:,:,jk) - vn_adv(:,:)*r1_hv_n(:,:) + vn_b(:,:) )*vmask(:,:,jk)
139            END DO 
140         ENDIF
141      ENDIF
142
143      ! Update after velocity on domain lateral boundaries
144      ! --------------------------------------------------     
145# if defined key_agrif
146      CALL Agrif_dyn( kt )             !* AGRIF zoom boundaries
147# endif
148      !
149      CALL lbc_lnk_multi( 'dynnxt', uu(:,:,:,kNp1), 'U', -1., vv(:,:,:,kNp1), 'V', -1. )     !* local domain boundaries
150      !
151      !                                !* BDY open boundaries
152      IF( ln_bdy .AND. ln_dynspg_exp )   CALL bdy_dyn( kt )
153      IF( ln_bdy .AND. ln_dynspg_ts  )   CALL bdy_dyn( kt, dyn3d_only=.true. )
154
155!!$   Do we need a call to bdy_vol here??
156      !
157      IF( l_trddyn ) THEN             ! prepare the atf trend computation + some diagnostics
158         z1_2dt = 1._wp / (2. * rdt)        ! Euler or leap-frog time step
159         IF( neuler == 0 .AND. kt == nit000 )   z1_2dt = 1._wp / rdt
160         !
161         !                                  ! Kinetic energy and Conversion
162         IF( ln_KE_trd  )   CALL trd_dyn( uu(:,:,:,kNp1), vv(:,:,:,kNp1), jpdyn_ken, kt )
163         !
164         IF( ln_dyn_trd ) THEN              ! 3D output: total momentum trends
165            zua(:,:,:) = ( uu(:,:,:,kNp1) - uu(:,:,:,kNm1) ) * z1_2dt
166            zva(:,:,:) = ( vv(:,:,:,kNp1) - vv(:,:,:,kNm1) ) * z1_2dt
167            CALL iom_put( "utrd_tot", zua )        ! total momentum trends, except the asselin time filter
168            CALL iom_put( "vtrd_tot", zva )
169         ENDIF
170         !
171         zua(:,:,:) = puu(:,:,:)             ! save the now velocity before the asselin filter
172         zva(:,:,:) = pvv(:,:,:)             ! (caution: there will be a shift by 1 timestep in the
173         !                                  !  computation of the asselin filter trends)
174      ENDIF
175
176      ! Time filter and swap of dynamics arrays
177      ! ------------------------------------------
178!! TO BE DELETED
179!!$      IF( neuler == 0 .AND. kt == nit000 ) THEN        !* Euler at first time-step: only swap
180!!$         DO jk = 1, jpkm1
181!!$            un(:,:,jk) = ua(:,:,jk)                         ! un <-- ua
182!!$            vn(:,:,jk) = va(:,:,jk)
183!!$         END DO
184!!$         IF( .NOT.ln_linssh ) THEN                          ! e3._b <-- e3._n
185!!gm BUG ????    I don't understand why it is not : e3._n <-- e3._a 
186!!$            DO jk = 1, jpkm1
187!               e3t_b(:,:,jk) = e3t_n(:,:,jk)
188!               e3u_b(:,:,jk) = e3u_n(:,:,jk)
189!               e3v_b(:,:,jk) = e3v_n(:,:,jk)
190               !
191!!$               e3t_n(:,:,jk) = e3t_a(:,:,jk)
192!!$               e3u_n(:,:,jk) = e3u_a(:,:,jk)
193!!$               e3v_n(:,:,jk) = e3v_a(:,:,jk)
194!!$            END DO
195!!gm BUG end
196!!$         ENDIF
197!! TO BE DELETED
198                                                            !
199         
200      IF( .NOT.( neuler == 0 .AND. kt == nit000 ) ) THEN    !* Leap-Frog : Asselin filter and swap
201         !                                ! =============!
202         IF( ln_linssh ) THEN             ! Fixed volume !
203            !                             ! =============!
204            DO jk = 1, jpkm1                             
205               DO jj = 1, jpj
206                  DO ji = 1, jpi   
207                     puu(ji,jj,jk) = puu(ji,jj,jk) + atfp * ( uu(ji,jj,jk,kNm1) - 2._wp * puu(ji,jj,jk) + uu(ji,jj,jk,kNp1) )
208                     pvv(ji,jj,jk) = pvv(ji,jj,jk) + atfp * ( vv(ji,jj,jk,kNm1) - 2._wp * pvv(ji,jj,jk) + vv(ji,jj,jk,kNp1) )
209                     !
210!! TO BE DELETED
211!!$                     ub(ji,jj,jk) = zuf                      ! ub <-- filtered velocity
212!!$                     vb(ji,jj,jk) = zvf
213!!$                     un(ji,jj,jk) = ua(ji,jj,jk)             ! un <-- ua
214!!$                     vn(ji,jj,jk) = va(ji,jj,jk)
215!! TO BE DELETED
216                  END DO
217               END DO
218            END DO
219            !                             ! ================!
220         ELSE                             ! Variable volume !
221            !                             ! ================!
222            ! Time-filtered scale factor at t-points
223            ! ----------------------------------------------------
224            ALLOCATE( ze3t_f(jpi,jpj,jpk) )
225            DO jk = 1, jpkm1
226               ze3t_f(:,:,jk) = pe3t(:,:,jk) + atfp * ( e3t(:,:,jk,kNm1) - 2._wp * pe3t(:,:,jk) + e3t(:,:,jk,kNp1) )
227            END DO
228            ! Add volume filter correction: compatibility with tracer advection scheme
229            ! => time filter + conservation correction (only at the first level)
230            zcoef = atfp * rdt * r1_rau0
231
232            ze3t_f(:,:,1) = ze3t_f(:,:,1) - zcoef * ( emp_b(:,:) - emp(:,:) ) * tmask(:,:,1)
233
234            IF ( ln_rnf ) THEN
235               IF( ln_rnf_depth ) THEN
236                  DO jk = 1, jpkm1 ! Deal with Rivers separetely, as can be through depth too
237                     DO jj = 1, jpj
238                        DO ji = 1, jpi
239                           IF( jk <=  nk_rnf(ji,jj)  ) THEN
240                               ze3t_f(ji,jj,jk) =   ze3t_f(ji,jj,jk) - zcoef *  ( - rnf_b(ji,jj) + rnf(ji,jj) ) &
241                                      &          * ( pe3t(ji,jj,jk) / h_rnf(ji,jj) ) * tmask(ji,jj,jk)
242                           ENDIF
243                        ENDDO
244                     ENDDO
245                  ENDDO
246               ELSE
247                  ze3t_f(:,:,1) = ze3t_f(:,:,1) - zcoef *  ( -rnf_b(:,:) + rnf(:,:))*tmask(:,:,1)
248               ENDIF
249            END IF
250
251            IF ( ln_isf ) THEN   ! if ice shelf melting
252               DO jk = 1, jpkm1 ! Deal with isf separetely, as can be through depth too
253                  DO jj = 1, jpj
254                     DO ji = 1, jpi
255                        IF( misfkt(ji,jj) <=jk .and. jk < misfkb(ji,jj)  ) THEN
256                           ze3t_f(ji,jj,jk) = ze3t_f(ji,jj,jk) - zcoef * ( fwfisf_b(ji,jj) - fwfisf(ji,jj) ) &
257                                &          * ( pe3t(ji,jj,jk) * r1_hisf_tbl(ji,jj) ) * tmask(ji,jj,jk)
258                        ELSEIF ( jk==misfkb(ji,jj) ) THEN
259                           ze3t_f(ji,jj,jk) = ze3t_f(ji,jj,jk) - zcoef * ( fwfisf_b(ji,jj) - fwfisf(ji,jj) ) &
260                                &          * ( pe3t(ji,jj,jk) * r1_hisf_tbl(ji,jj) ) * ralpha(ji,jj) * tmask(ji,jj,jk)
261                        ENDIF
262                     END DO
263                  END DO
264               END DO
265            END IF
266            !
267            pe3t(:,:,1:jpkm1) = ze3t_f(:,:,1:jpkm1)        ! filtered scale factor at T-points
268            !
269            IF( ln_dynadv_vec ) THEN      ! Asselin filter applied on velocity
270               ! Before filtered scale factor at (u/v)-points
271               CALL dom_vvl_interpol( pe3t(:,:,:), pe3u(:,:,:), 'U' )
272               CALL dom_vvl_interpol( pe3t(:,:,:), pe3v(:,:,:), 'V' )
273               DO jk = 1, jpkm1
274                  DO jj = 1, jpj
275                     DO ji = 1, jpi
276                        puu(ji,jj,jk) = puu(ji,jj,jk) + atfp * ( uu(ji,jj,jk,kNm1) - 2._wp * puu(ji,jj,jk) + uu(ji,jj,jk,kNp1) )
277                        pvv(ji,jj,jk) = pvv(ji,jj,jk) + atfp * ( vv(ji,jj,jk,kNm1) - 2._wp * pvv(ji,jj,jk) + vv(ji,jj,jk,kNp1) )
278                        !
279!! TO BE DELETED
280!!$                        ub(ji,jj,jk) = zuf                      ! ub <-- filtered velocity
281!!$                        vb(ji,jj,jk) = zvf
282!!$                        un(ji,jj,jk) = ua(ji,jj,jk)             ! un <-- ua
283!!$                        vn(ji,jj,jk) = va(ji,jj,jk)
284!! TO BE DELETED
285                     END DO
286                  END DO
287               END DO
288               !
289            ELSE                          ! Asselin filter applied on thickness weighted velocity
290               !
291               ALLOCATE( ze3u_f(jpi,jpj,jpk) , ze3v_f(jpi,jpj,jpk) )
292               ! Now filtered scale factor at (u/v)-points stored in ze3u_f, ze3v_f
293               CALL dom_vvl_interpol( pe3t(:,:,:), ze3u_f, 'U' )
294               CALL dom_vvl_interpol( pe3t(:,:,:), ze3v_f, 'V' )
295               DO jk = 1, jpkm1
296                  DO jj = 1, jpj
297                     DO ji = 1, jpi                 
298                        zue3a = e3u(ji,jj,jk,kNp1) * uu(ji,jj,jk,kNp1)
299                        zve3a = e3v(ji,jj,jk,kNp1) * vv(ji,jj,jk,kNp1)
300                        zue3n = pe3u(ji,jj,jk)     * puu(ji,jj,jk)
301                        zve3n = pe3v(ji,jj,jk)     * pvv(ji,jj,jk)
302                        zue3b = e3u(ji,jj,jk,kNm1) * uu(ji,jj,jk,kNm1)
303                        zve3b = e3v(ji,jj,jk,kNm1) * vv(ji,jj,jk,kNm1)
304                        !
305                        puu(ji,jj,jk) = ( zue3n + atfp * ( zue3b - 2._wp * zue3n  + zue3a ) ) / ze3u_f(ji,jj,jk)
306                        pvv(ji,jj,jk) = ( zve3n + atfp * ( zve3b - 2._wp * zve3n  + zve3a ) ) / ze3v_f(ji,jj,jk)
307                        !
308!! TO BE DELETED
309!!$                        ub(ji,jj,jk) = zuf                     ! ub <-- filtered velocity
310!!$                        vb(ji,jj,jk) = zvf
311!!$                        un(ji,jj,jk) = ua(ji,jj,jk)            ! un <-- ua
312!!$                        vn(ji,jj,jk) = va(ji,jj,jk)
313!! TO BE DELETED
314                     END DO
315                  END DO
316               END DO
317               pe3u(:,:,1:jpkm1) = ze3u_f(:,:,1:jpkm1) 
318               pe3v(:,:,1:jpkm1) = ze3v_f(:,:,1:jpkm1)
319               !
320               DEALLOCATE( ze3u_f , ze3v_f )
321            ENDIF
322            !
323            DEALLOCATE( ze3t_f )
324         ENDIF
325         !
326         IF( ln_dynspg_ts .AND. ln_bt_fw ) THEN
327            ! Revert filtered "now" velocities to time split estimate
328            ! Doing it here also means that asselin filter contribution is removed 
329            zue(:,:) = pe3u(:,:,1) * puu(:,:,1) * umask(:,:,1)
330            zve(:,:) = pe3v(:,:,1) * pvv(:,:,1) * vmask(:,:,1)   
331            DO jk = 2, jpkm1
332               zue(:,:) = zue(:,:) + pe3u(:,:,jk) * puu(:,:,jk) * umask(:,:,jk)
333               zve(:,:) = zve(:,:) + pe3v(:,:,jk) * pvv(:,:,jk) * vmask(:,:,jk)   
334            END DO
335            DO jk = 1, jpkm1
336               puu(:,:,jk) = puu(:,:,jk) - (zue(:,:) * r1_hu_n(:,:) - un_b(:,:)) * umask(:,:,jk)
337               pvv(:,:,jk) = pvv(:,:,jk) - (zve(:,:) * r1_hv_n(:,:) - vn_b(:,:)) * vmask(:,:,jk)
338            END DO
339         ENDIF
340         !
341      ENDIF ! neuler =/0
342      !
343      ! Set "now" and "before" barotropic velocities for next time step:
344      ! JC: Would be more clever to swap variables than to make a full vertical
345      ! integration
346      !
347      ! DS IMMERSE :: This is very confusing as it stands. But should
348      ! hopefully be simpler when we do the time-level swapping for the
349      ! external mode solution.
350      !
351      IF(.NOT.ln_linssh ) THEN
352         hu_b(:,:) = pe3u(:,:,1) * umask(:,:,1)
353         hv_b(:,:) = pe3v(:,:,1) * vmask(:,:,1)
354         DO jk = 2, jpkm1
355            hu_b(:,:) = hu_b(:,:) + pe3u(:,:,jk) * umask(:,:,jk)
356            hv_b(:,:) = hv_b(:,:) + pe3v(:,:,jk) * vmask(:,:,jk)
357         END DO
358         r1_hu_b(:,:) = ssumask(:,:) / ( hu_b(:,:) + 1._wp - ssumask(:,:) )
359         r1_hv_b(:,:) = ssvmask(:,:) / ( hv_b(:,:) + 1._wp - ssvmask(:,:) )
360      ENDIF
361      !
362      un_b(:,:) = e3u(:,:,1,kNp1) * uu(:,:,1,kNp1) * umask(:,:,1)
363      ub_b(:,:) = pe3u(:,:,1)     * puu(:,:,1) * umask(:,:,1)
364      vn_b(:,:) = e3v(:,:,1,kNp1) * vv(:,:,1,kNp1) * vmask(:,:,1)
365      vb_b(:,:) = pe3v(:,:,1)     * pvv(:,:,1) * vmask(:,:,1)
366      DO jk = 2, jpkm1
367         un_b(:,:) = un_b(:,:) + e3u(:,:,jk,kNp1) * uu(:,:,jk,kNp1) * umask(:,:,jk)
368         ub_b(:,:) = ub_b(:,:) + pe3u(:,:,jk)     * puu(:,:,jk) * umask(:,:,jk)
369         vn_b(:,:) = vn_b(:,:) + e3v(:,:,jk,kNp1) * vv(:,:,jk,kNp1) * vmask(:,:,jk)
370         vb_b(:,:) = vb_b(:,:) + pe3v(:,:,jk)     * pvv(:,:,jk) * vmask(:,:,jk)
371      END DO
372      un_b(:,:) = un_b(:,:) * r1_hu_a(:,:)
373      vn_b(:,:) = vn_b(:,:) * r1_hv_a(:,:)
374      ub_b(:,:) = ub_b(:,:) * r1_hu_b(:,:)
375      vb_b(:,:) = vb_b(:,:) * r1_hv_b(:,:)
376      !
377      IF( .NOT.ln_dynspg_ts ) THEN        ! output the barotropic currents
378         CALL iom_put(  "ubar", un_b(:,:) )
379         CALL iom_put(  "vbar", vn_b(:,:) )
380      ENDIF
381      IF( l_trddyn ) THEN                ! 3D output: asselin filter trends on momentum
382         zua(:,:,:) = ( puu(:,:,:) - zua(:,:,:) ) * z1_2dt
383         zva(:,:,:) = ( pvv(:,:,:) - zva(:,:,:) ) * z1_2dt
384         CALL trd_dyn( zua, zva, jpdyn_atf, kt )
385      ENDIF
386      !
387      IF(ln_ctl)   CALL prt_ctl( tab3d_1=uu(:,:,:,kNp1), clinfo1=' nxt  - uu(:,:,:,kNp1): ', mask1=umask,   &
388         &                       tab3d_2=vv(:,:,:,kNp1), clinfo2=' vv(:,:,:,kNp1): '       , mask2=vmask )
389      !
390      IF( ln_dynspg_ts )   DEALLOCATE( zue, zve )
391      IF( l_trddyn     )   DEALLOCATE( zua, zva )
392      IF( ln_timing    )   CALL timing_stop('dyn_nxt')
393      !
394   END SUBROUTINE dyn_nxt
395
396   !!=========================================================================
397END MODULE dynnxt
Note: See TracBrowser for help on using the repository browser.