source: NEMO/branches/2018/dev_r9838_ENHANCE04_RK3/src/NST/agrif_oce_update.F90 @ 10030

Last change on this file since 10030 was 10030, checked in by gm, 3 years ago

#1911 (ENHANCE-04): RK3 branch - step II.3 remove e3uw_$ e3vw_$, except e3.w_0 and use only after e3 in dyn/trazdf

  • Property svn:keywords set to Id
File size: 57.3 KB
Line 
1#define TWO_WAY        /* TWO WAY NESTING */
2#undef DECAL_FEEDBACK  /* SEPARATION of INTERFACES*/
3#undef VOL_REFLUX      /* VOLUME REFLUXING*/
4 
5MODULE agrif_oce_update
6   !!======================================================================
7   !!                   ***  MODULE  agrif_oce_interp  ***
8   !! AGRIF: update package for the ocean dynamics (OPA)
9   !!======================================================================
10   !! History :  2.0  !  2002-06  (L. Debreu)  Original code
11   !!            3.2  !  2009-04  (R. Benshila)
12   !!            3.6  !  2014-09  (R. Benshila)
13   !!----------------------------------------------------------------------
14
15   !!----------------------------------------------------------------------
16   !!   Agrif_Update_Tra   : T-S agrif update
17   !!   Agrif_Update_Dyn   : dynamics agrif update
18   !!   Agrif_Update_ssh   : sea surface height update
19   !!   Agrif_Update_Tke   :
20   !!   Agrif_Update_vvl   :
21   !!   dom_vvl_update_UVF :
22   !!   updateTS           :
23   !!   updateu            :
24   !!   correct_u_bdy      :
25   !!   updatev            :
26   !!   correct_v_bdy      :
27   !!   updateu2d          :
28   !!   updatev2d          :
29   !!   updateSSH          :
30   !!   updateub2b         :
31   !!   reflux_sshu        :
32   !!   updatevb2b         :
33   !!   reflux_sshv        :
34   !!   update_scales      :
35   !!   updateEN           :
36   !!   updateAVT          :
37   !!   updateAVM          :
38   !!   updatee3t          :
39   !!----------------------------------------------------------------------
40
41#if defined key_agrif 
42   !!----------------------------------------------------------------------
43   !!   'key_agrif'                                              AGRIF zoom
44   !!----------------------------------------------------------------------
45   USE par_oce        ! ocean parameter
46   USE oce            ! ocean variables
47   USE dom_oce        ! ocean domain
48   USE zdf_oce        ! vertical physics: ocean variables
49   USE agrif_oce      !
50   !
51   USE in_out_manager ! I/O manager
52   USE lib_mpp        ! MPP library
53   USE domvvl         ! Need interpolation routines
54
55   IMPLICIT NONE
56   PRIVATE
57
58   PUBLIC   Agrif_Update_Tra, Agrif_Update_Dyn, Agrif_Update_vvl, Agrif_Update_ssh
59   PUBLIC   Update_Scales
60
61   !!----------------------------------------------------------------------
62   !! NEMO/NST 4.0 , NEMO Consortium (2018)
63   !! $Id$
64   !! Software governed by the CeCILL licence (./LICENSE)
65   !!----------------------------------------------------------------------
66CONTAINS
67
68   SUBROUTINE Agrif_Update_Tra( )
69      !!----------------------------------------------------------------------
70      !!                   *** ROUTINE Agrif_Update_Tra ***
71      !!----------------------------------------------------------------------
72      !
73      IF (Agrif_Root()) RETURN
74      !
75#if defined TWO_WAY 
76      IF (lwp.AND.lk_agrif_debug) Write(*,*) 'Update tracers  from grid Number',Agrif_Fixed()
77
78      Agrif_UseSpecialValueInUpdate = .TRUE.
79      Agrif_SpecialValueFineGrid    = 0._wp
80      !
81# if ! defined DECAL_FEEDBACK
82      CALL Agrif_Update_Variable(tsn_id, procname=updateTS)
83! near boundary update:
84!      CALL Agrif_Update_Variable(tsn_id,locupdate=(/0,2/), procname=updateTS)
85# else
86      CALL Agrif_Update_Variable(tsn_id, locupdate=(/1,0/),procname=updateTS)
87! near boundary update:
88!      CALL Agrif_Update_Variable(tsn_id,locupdate=(/1,2/), procname=updateTS)
89# endif
90      !
91      Agrif_UseSpecialValueInUpdate = .FALSE.
92      !
93#endif
94      !
95   END SUBROUTINE Agrif_Update_Tra
96
97
98   SUBROUTINE Agrif_Update_Dyn( )
99      !!----------------------------------------------------------------------
100      !!                   *** ROUTINE Agrif_Update_Dyn ***
101      !!----------------------------------------------------------------------
102      !
103      IF (Agrif_Root()) RETURN
104      !
105#if defined TWO_WAY
106      IF (lwp.AND.lk_agrif_debug) Write(*,*) 'Update momentum from grid Number',Agrif_Fixed()
107
108      Agrif_UseSpecialValueInUpdate = .FALSE.
109      Agrif_SpecialValueFineGrid = 0.
110      !     
111# if ! defined DECAL_FEEDBACK
112      CALL Agrif_Update_Variable(un_update_id,procname = updateU)
113      CALL Agrif_Update_Variable(vn_update_id,procname = updateV)
114! near boundary update:
115!      CALL Agrif_Update_Variable(un_update_id,locupdate=(/0,1/),procname = updateU)
116!      CALL Agrif_Update_Variable(vn_update_id,locupdate=(/0,1/),procname = updateV)
117# else
118      CALL Agrif_Update_Variable(un_update_id,locupdate1=(/0,-1/),locupdate2=(/1,-2/),procname = updateU)
119      CALL Agrif_Update_Variable(vn_update_id,locupdate1=(/1,-2/),locupdate2=(/0,-1/),procname = updateV)
120! near boundary update:
121!      CALL Agrif_Update_Variable(un_update_id,locupdate1=(/0,1/),locupdate2=(/1,1/),procname = updateU)
122!      CALL Agrif_Update_Variable(vn_update_id,locupdate1=(/1,1/),locupdate2=(/0,1/),procname = updateV)
123# endif
124
125# if ! defined DECAL_FEEDBACK
126      CALL Agrif_Update_Variable(e1u_id,procname = updateU2d)
127      CALL Agrif_Update_Variable(e2v_id,procname = updateV2d) 
128# else
129      CALL Agrif_Update_Variable(e1u_id,locupdate1=(/0,-1/),locupdate2=(/1,-2/),procname = updateU2d)
130      CALL Agrif_Update_Variable(e2v_id,locupdate1=(/1,-2/),locupdate2=(/0,-1/),procname = updateV2d) 
131# endif
132      !
133# if ! defined DECAL_FEEDBACK
134      ! Account for updated thicknesses at boundary edges
135      IF (.NOT.ln_linssh) THEN
136!         CALL Agrif_Update_Variable(un_update_id,locupdate1=(/0,0/),locupdate2=(/0,0/),procname = correct_u_bdy)
137!         CALL Agrif_Update_Variable(vn_update_id,locupdate1=(/0,0/),locupdate2=(/0,0/),procname = correct_v_bdy)
138      ENDIF
139# endif
140      !
141      IF ( ln_dynspg_ts .AND. ln_bt_fw ) THEN
142         ! Update time integrated transports
143#  if ! defined DECAL_FEEDBACK
144         CALL Agrif_Update_Variable(ub2b_update_id,procname = updateub2b)
145         CALL Agrif_Update_Variable(vb2b_update_id,procname = updatevb2b)
146#  else
147         CALL Agrif_Update_Variable(ub2b_update_id,locupdate1=(/0,-1/),locupdate2=(/1,-2/),procname = updateub2b)
148         CALL Agrif_Update_Variable(vb2b_update_id,locupdate1=(/1,-2/),locupdate2=(/0,-1/),procname = updatevb2b)
149#  endif
150      END IF
151#endif
152      !
153   END SUBROUTINE Agrif_Update_Dyn
154
155
156   SUBROUTINE Agrif_Update_ssh( )
157      !!----------------------------------------------------------------------
158      !!                   *** ROUTINE Agrif_Update_ssh ***
159      !!----------------------------------------------------------------------
160      !
161      IF (Agrif_Root()) RETURN
162      !
163#if defined TWO_WAY
164      !
165      Agrif_UseSpecialValueInUpdate = .TRUE.
166      Agrif_SpecialValueFineGrid = 0.
167# if ! defined DECAL_FEEDBACK
168      CALL Agrif_Update_Variable(sshn_id,procname = updateSSH)
169# else
170      CALL Agrif_Update_Variable(sshn_id,locupdate=(/1,0/),procname = updateSSH)
171# endif
172      !
173      Agrif_UseSpecialValueInUpdate = .FALSE.
174      !
175#  if defined VOL_REFLUX
176      IF ( ln_dynspg_ts.AND.ln_bt_fw ) THEN
177         ! Refluxing on ssh:
178#  if defined DECAL_FEEDBACK
179         CALL Agrif_Update_Variable(ub2b_update_id,locupdate1=(/0, 0/),locupdate2=(/1, 1/),procname = reflux_sshu)
180         CALL Agrif_Update_Variable(vb2b_update_id,locupdate1=(/1, 1/),locupdate2=(/0, 0/),procname = reflux_sshv)
181#  else
182         CALL Agrif_Update_Variable(ub2b_update_id,locupdate1=(/-1,-1/),locupdate2=(/ 0, 0/),procname = reflux_sshu)
183         CALL Agrif_Update_Variable(vb2b_update_id,locupdate1=(/ 0, 0/),locupdate2=(/-1,-1/),procname = reflux_sshv)
184#  endif
185      END IF
186#  endif
187      !
188#endif
189      !
190   END SUBROUTINE Agrif_Update_ssh
191
192
193   SUBROUTINE Agrif_Update_Tke( )
194      !!----------------------------------------------------------------------
195      !!                   *** ROUTINE Agrif_Update_Tke ***
196      !!----------------------------------------------------------------------
197      !
198      IF (Agrif_Root()) RETURN
199      !       
200#  if defined TWO_WAY
201      !
202      Agrif_UseSpecialValueInUpdate = .TRUE.
203      Agrif_SpecialValueFineGrid = 0.
204      !
205      CALL Agrif_Update_Variable( en_id, locupdate=(/0,0/), procname=updateEN  )
206      CALL Agrif_Update_Variable(avt_id, locupdate=(/0,0/), procname=updateAVT )
207      CALL Agrif_Update_Variable(avm_id, locupdate=(/0,0/), procname=updateAVM )
208      !
209      Agrif_UseSpecialValueInUpdate = .FALSE.
210      !
211#  endif
212      !
213   END SUBROUTINE Agrif_Update_Tke
214
215
216   SUBROUTINE Agrif_Update_vvl( )
217      !!----------------------------------------------------------------------
218      !!                   *** ROUTINE Agrif_Update_vvl ***
219      !!----------------------------------------------------------------------
220      !
221      IF ( Agrif_Root() )   RETURN
222      !
223#if defined TWO_WAY 
224      !
225      IF (lwp.AND.lk_agrif_debug) Write(*,*) 'Update e3 from grid Number',Agrif_Fixed(), 'Step', Agrif_Nb_Step()
226      !
227      Agrif_UseSpecialValueInUpdate = .TRUE.
228      Agrif_SpecialValueFineGrid = 0.
229      !
230      ! No interface separation here, update vertical grid at T points
231      ! everywhere over the overlapping regions (one account for refluxing in that case):
232      CALL Agrif_Update_Variable(e3t_id, procname=updatee3t) 
233      !
234      Agrif_UseSpecialValueInUpdate = .FALSE.
235      !
236      CALL Agrif_ChildGrid_To_ParentGrid()
237      CALL dom_vvl_update_UVF
238      CALL Agrif_ParentGrid_To_ChildGrid()
239      !
240#endif
241      !
242   END SUBROUTINE Agrif_Update_vvl
243
244
245   SUBROUTINE dom_vvl_update_UVF
246      !!----------------------------------------------------------------------
247      !!                   *** ROUTINE dom_vvl_update_UVF ***
248      !!----------------------------------------------------------------------
249      INTEGER ::   jk      ! dummy loop index
250      REAL(wp)::   zcoef   ! local scalar
251      REAL(wp), DIMENSION(jpi,jpj) ::   zssht_h, zsshu_h, zsshv_h, zsshf_h
252      !!----------------------------------------------------------------------
253      !
254      IF (lwp.AND.lk_agrif_debug)   Write(*,*) 'Finalize e3 on grid Number', Agrif_Fixed(), 'Step', Agrif_Nb_Step()
255
256      ! Save "old" scale factor (prior update) for subsequent asselin correction of prognostic variables
257      ! -----------------------
258      e3u_a(:,:,:) = e3u_n(:,:,:)
259      e3v_a(:,:,:) = e3v_n(:,:,:)
260      hu_a (:,:)   = hu_n (:,:)
261      hv_a (:,:)   = hv_n (:,:)
262      !
263      !                          !* NOW fields :
264      CALL ssh2e3_now                  ! set:  ht,  hu,  hv, r1_hu, r1_hv
265      !                                !      e3t, e3u, e3v, e3f               (from 1 to jpkm1)
266      !                                !      e3w, gdept_n, gdepw_n, gde3w_n   (from 1 to jpk  )
267
268      !                          !* BEFORE fields :
269      IF (.NOT.( lk_agrif_fstep .AND. l_1st_euler ) ) THEN
270         !
271         CALL ssh2e3_before               ! set:       hu,  hv, r1_hu, r1_hv
272         !                                !      e3t, e3u, e3v                 (from 1 to jpkm1)
273         !                                !      e3w                           (from 1 to jpk  )
274      ENDIF
275      !
276   END SUBROUTINE dom_vvl_update_UVF
277
278# if defined key_vertical
279
280   SUBROUTINE updateTS( tabres, i1, i2, j1, j2, k1, k2, n1, n2, before )
281      !!----------------------------------------------------------------------
282      !!                   *** ROUTINE updateT ***
283      !!----------------------------------------------------------------------
284      INTEGER, INTENT(in) :: i1,i2,j1,j2,k1,k2,n1,n2
285      REAL(wp),DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) :: tabres
286      LOGICAL, INTENT(in) :: before
287      !!
288      INTEGER :: ji,jj,jk,jn
289      REAL(wp), DIMENSION(i1:i2,j1:j2,1:jpk,n1:n2) :: tabres_child
290      REAL(wp) :: h_in(k1:k2)
291      REAL(wp) :: h_out(1:jpk)
292      INTEGER  :: N_in, N_out
293      REAL(wp) :: zrho_xy, h_diff
294      REAL(wp) :: tabin(k1:k2,n1:n2)
295      !!----------------------------------------------------------------------
296      !
297      IF (before) THEN
298         AGRIF_SpecialValue = -999._wp
299         zrho_xy = Agrif_rhox() * Agrif_rhoy() 
300         DO jn = n1, n2-1
301            DO jk = k1, k2
302               DO jj = j1, j2
303                  DO ji = i1, i2
304                     tabres(ji,jj,jk,jn) =  tmask(ji,jj,jk)    * (tsn(ji,jj,jk,jn) * e3t_n(ji,jj,jk) )   &
305                                         + (tmask(ji,jj,jk)-1) * 999._wp
306                  END DO
307               END DO
308            END DO
309         END DO
310         DO jk = k1, k2
311            DO jj = j1, j2
312               DO ji = i1, i2
313                  tabres(ji,jj,jk,n2) =  tmask(ji,jj,jk)    * e3t_n(ji,jj,jk)   &
314                     &                + (tmask(ji,jj,jk)-1) * 999._wp
315               END DO
316            END DO
317         END DO
318      ELSE
319         tabres_child(:,:,:,:) = 0.
320         AGRIF_SpecialValue = 0._wp
321         DO jj = j1 , j2
322            DO ji = i1, i2
323               N_in = 0
324               DO jk = k1, k2 !k2 = jpk of child grid
325                  IF ( tabres(ji,jj,jk,n2) == 0  )   EXIT
326                  N_in = N_in + 1
327                  tabin(jk,:) = tabres(ji,jj,jk,n1:n2-1)/tabres(ji,jj,jk,n2)
328                  h_in (N_in) = tabres(ji,jj,jk,n2)
329               END DO
330               N_out = 0
331               DO jk = 1, jpk ! jpk of parent grid
332                  IF (tmask(ji,jj,jk) < -900)   EXIT ! TODO: Will not work with ISF
333                  N_out = N_out + 1
334                  h_out(N_out) = e3t_n(ji,jj,jk) 
335               END DO
336               IF (N_in > 0) THEN !Remove this?
337                  h_diff = sum(h_out(1:N_out))-sum(h_in(1:N_in))
338                  IF (h_diff < -1.e-4) THEN
339                     print *,'CHECK YOUR bathy T points ...',ji,jj,h_diff,sum(h_in(1:N_in)),sum(h_out(1:N_out))
340                     print *,h_in(1:N_in)
341                     print *,h_out(1:N_out)
342                     STOP
343                  ENDIF
344                  DO jn = n1, n2-1
345                     CALL reconstructandremap(tabin(1:N_in,jn),h_in(1:N_in),tabres_child(ji,jj,1:N_out,jn),h_out(1:N_out),N_in,N_out)
346                  END DO
347               ENDIF
348            END DO
349         END DO
350
351         IF ( .NOT.( lk_agrif_fstep .AND. l_1st_euler ) ) THEN       ! Add asselin part
352            DO jn = n1, n2-1
353               DO jk = 1, jpk
354                  DO jj = j1, j2
355                     DO ji = i1, i2
356                        IF( tabres_child(ji,jj,jk,jn) /= 0. ) THEN
357                           tsb(ji,jj,jk,jn) = tsb(ji,jj,jk,jn)   & 
358                              &             + rn_atfp * ( tabres_child(ji,jj,jk,jn) - tsn(ji,jj,jk,jn) ) * tmask(ji,jj,jk)
359                        ENDIF
360                     END DO
361                  END DO
362               END DO
363            END DO
364         ENDIF
365         DO jn = n1, n2-1
366            DO jk = 1, jpk
367               DO jj = j1, j2
368                  DO ji = i1, i2
369                     IF( tabres_child(ji,jj,jk,jn) /= 0. ) THEN
370                        tsn(ji,jj,jk,jn) = tabres_child(ji,jj,jk,jn) * tmask(ji,jj,jk)
371                     END IF
372                  END DO
373               END DO
374            END DO
375         END DO
376      ENDIF
377      !
378   END SUBROUTINE updateTS
379
380# else
381
382   SUBROUTINE updateTS( tabres, i1, i2, j1, j2, k1, k2, n1, n2, before )
383      !!----------------------------------------------------------------------
384      !!                   *** ROUTINE ROUTINE updateT ***
385      !!----------------------------------------------------------------------
386      INTEGER, INTENT(in) :: i1,i2,j1,j2,k1,k2,n1,n2
387      REAL(wp),DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) :: tabres
388      LOGICAL, INTENT(in) :: before
389      !
390      INTEGER :: ji,jj,jk,jn
391      REAL(wp) :: ztb, ztnu, ztno
392      !!----------------------------------------------------------------------
393      !
394      IF (before) THEN
395         DO jn = 1,jpts
396            DO jk=k1,k2
397               DO jj=j1,j2
398                  DO ji=i1,i2
399!> jc tmp
400                     tabres(ji,jj,jk,jn) = tsn(ji,jj,jk,jn)  * e3t_n(ji,jj,jk) / e3t_0(ji,jj,jk)
401!                     tabres(ji,jj,jk,jn) = tsn(ji,jj,jk,jn)  * e3t_n(ji,jj,jk)
402!< jc tmp
403                  END DO
404               END DO
405            END DO
406         END DO
407      ELSE
408!> jc tmp
409         DO jn = 1,jpts
410            tabres(i1:i2,j1:j2,k1:k2,jn) =  tabres(i1:i2,j1:j2,k1:k2,jn) * e3t_0(i1:i2,j1:j2,k1:k2) &
411                                         & * tmask(i1:i2,j1:j2,k1:k2)
412         END DO
413!< jc tmp
414         IF ( .NOT.( lk_agrif_fstep .AND. l_1st_euler ) ) THEN    ! Add asselin part
415            DO jn = 1,jpts
416               DO jk = k1, k2
417                  DO jj = j1, j2
418                     DO ji = i1, i2
419                        IF( tabres(ji,jj,jk,jn) /= 0._wp ) THEN
420                           ztb  = tsb(ji,jj,jk,jn) * e3t_b(ji,jj,jk) ! fse3t_b prior update should be used
421                           ztnu = tabres(ji,jj,jk,jn)
422                           ztno = tsn(ji,jj,jk,jn) * e3t_a(ji,jj,jk)
423                           tsb(ji,jj,jk,jn) = ( ztb + rn_atfp * ( ztnu - ztno) ) / e3t_b(ji,jj,jk) * tmask(ji,jj,jk)
424                        ENDIF
425                     END DO
426                  END DO
427               END DO
428            END DO
429         ENDIF
430         DO jn = 1,jpts
431            DO jk=k1,k2
432               DO jj=j1,j2
433                  DO ji=i1,i2
434                     IF( tabres(ji,jj,jk,jn) /= 0._wp ) THEN
435                        tsn(ji,jj,jk,jn) = tabres(ji,jj,jk,jn) / e3t_n(ji,jj,jk)
436                     END IF
437                  END DO
438               END DO
439            END DO
440         END DO
441         !
442         IF ( l_1st_euler .AND. Agrif_Nb_Step() == 0 ) THEN
443            tsb(i1:i2,j1:j2,k1:k2,1:jpts)  = tsn(i1:i2,j1:j2,k1:k2,1:jpts)
444         ENDIF
445         !
446      ENDIF
447      !
448   END SUBROUTINE updateTS
449
450# endif
451
452# if defined key_vertical
453
454   SUBROUTINE updateu( tabres, i1, i2, j1, j2, k1, k2, n1, n2, before )
455      !!----------------------------------------------------------------------
456      !!                   *** ROUTINE updateu ***
457      !!----------------------------------------------------------------------
458      INTEGER                                     , INTENT(in   ) :: i1, i2, j1, j2, k1, k2, n1, n2
459      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) :: tabres
460      LOGICAL                                     , INTENT(in   ) :: before
461      !
462      INTEGER ::   ji, jj, jk
463      REAL(wp)::   zrhoy
464! VERTICAL REFINEMENT BEGIN
465      REAL(wp), DIMENSION(i1:i2,j1:j2,1:jpk) :: tabres_child
466      REAL(wp) :: h_in(k1:k2)
467      REAL(wp) :: h_out(1:jpk)
468      INTEGER  :: N_in, N_out
469      REAL(wp) :: h_diff, excess, thick
470      REAL(wp) :: tabin(k1:k2)
471! VERTICAL REFINEMENT END
472      !!----------------------------------------------------------------------
473      !
474      IF( before ) THEN
475         zrhoy = Agrif_Rhoy()
476         AGRIF_SpecialValue = -999._wp
477         DO jk=k1,k2
478            DO jj=j1,j2
479               DO ji=i1,i2
480                  tabres(ji,jj,jk,1) = zrhoy * e2u(ji,jj) * e3u_n(ji,jj,jk) * umask(ji,jj,jk) * un(ji,jj,jk)  &
481                                       + (umask(ji,jj,jk)-1)*999._wp
482                  tabres(ji,jj,jk,2) = zrhoy * umask(ji,jj,jk) * e2u(ji,jj) * e3u_n(ji,jj,jk)  &
483                                       + (umask(ji,jj,jk)-1)*999._wp
484               END DO
485            END DO
486         END DO
487      ELSE
488         tabres_child(:,:,:) = 0.
489         AGRIF_SpecialValue = 0._wp
490         DO jj=j1,j2
491            DO ji=i1,i2
492               N_in = 0
493               h_in(:) = 0._wp
494               tabin(:) = 0._wp
495               DO jk=k1,k2 !k2=jpk of child grid
496                  IF( tabres(ji,jj,jk,2) < -900) EXIT
497                  N_in = N_in + 1
498                  tabin(jk) = tabres(ji,jj,jk,1)/tabres(ji,jj,jk,2)
499                  h_in(N_in) = tabres(ji,jj,jk,2)/e2u(ji,jj)
500               END DO
501               N_out = 0
502               DO jk=1,jpk
503                  IF (umask(ji,jj,jk) == 0) EXIT
504                  N_out = N_out + 1
505                  h_out(N_out) = e3u_n(ji,jj,jk)
506               END DO
507               IF (N_in * N_out > 0) THEN
508                  h_diff = sum(h_out(1:N_out))-sum(h_in(1:N_in))
509                  IF (h_diff < -1.e-4) THEN
510!Even if bathy at T points match it's possible for the U points to be deeper in the child grid.
511!In this case we need to move transport from the child grid cells below bed of parent grid into the bottom cell.
512                     excess = 0._wp
513                     DO jk=N_in,1,-1
514                        thick = MIN(-1*h_diff, h_in(jk))
515                        excess = excess + tabin(jk)*thick*e2u(ji,jj)
516                        tabin(jk) = tabin(jk)*(1. - thick/h_in(jk))
517                        h_diff = h_diff + thick
518                        IF ( h_diff == 0) THEN
519                           N_in = jk
520                           h_in(jk) = h_in(jk) - thick
521                           EXIT
522                        ENDIF
523                     END DO
524                  ENDIF
525                  CALL reconstructandremap(tabin(1:N_in),h_in(1:N_in),tabres_child(ji,jj,1:N_out),h_out(1:N_out),N_in,N_out)
526                  tabres_child(ji,jj,N_out) = tabres_child(ji,jj,N_out) + excess/(e2u(ji,jj)*h_out(N_out))
527               ENDIF
528            END DO
529         END DO
530
531         DO jk = 1, jpk
532            DO jj = j1, j2
533               DO ji = i1, i2
534                  IF ( .NOT.( lk_agrif_fstep .AND. l_1st_euler) ) THEN    ! Add asselin part
535!!gm                  IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN   ! Add asselin part
536                     ub(ji,jj,jk) = ub(ji,jj,jk) + rn_atfp * ( tabres_child(ji,jj,jk) - un(ji,jj,jk) ) * umask(ji,jj,jk)
537                  ENDIF
538                  un(ji,jj,jk) = tabres_child(ji,jj,jk) * umask(ji,jj,jk)
539               END DO
540            END DO
541         END DO
542      ENDIF
543      !
544   END SUBROUTINE updateu
545
546#else
547
548   SUBROUTINE updateu( tabres, i1, i2, j1, j2, k1, k2, n1, n2, before )
549      !!----------------------------------------------------------------------
550      !!                   *** ROUTINE updateu ***
551      !!----------------------------------------------------------------------
552      INTEGER                               , INTENT(in   ) :: i1, i2, j1, j2, k1, k2, n1, n2
553      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) :: tabres
554      LOGICAL                                     , INTENT(in   ) :: before
555      !
556      INTEGER  :: ji, jj, jk
557      REAL(wp) :: zrhoy, zub, zunu, zuno
558      !!----------------------------------------------------------------------
559      !
560      IF( before ) THEN
561         zrhoy = Agrif_Rhoy()
562         DO jk = k1, k2
563            tabres(i1:i2,j1:j2,jk,1) = zrhoy * e2u(i1:i2,j1:j2) * e3u_n(i1:i2,j1:j2,jk) * un(i1:i2,j1:j2,jk)
564         END DO
565      ELSE
566         DO jk=k1,k2
567            DO jj=j1,j2
568               DO ji=i1,i2
569                  tabres(ji,jj,jk,1) = tabres(ji,jj,jk,1) * r1_e2u(ji,jj) 
570                  !
571                  IF ( .NOT.( lk_agrif_fstep .AND. l_1st_euler ) ) THEN    ! Add asselin part
572!!gm                  IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN   ! Add asselin part
573                     zub  = ub(ji,jj,jk) * e3u_b(ji,jj,jk)  ! fse3t_b prior update should be used
574                     zuno = un(ji,jj,jk) * e3u_a(ji,jj,jk)
575                     zunu = tabres(ji,jj,jk,1)
576                     ub(ji,jj,jk) = ( zub + rn_atfp * ( zunu - zuno) ) / e3u_b(ji,jj,jk) * umask(ji,jj,jk)
577                  ENDIF
578                  !
579                  un(ji,jj,jk) = tabres(ji,jj,jk,1) / e3u_n(ji,jj,jk) * umask(ji,jj,jk)
580               END DO
581            END DO
582         END DO
583         !
584         IF ( l_1st_euler .AND. Agrif_Nb_Step() == 0 ) THEN
585!!gm         IF  ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN
586            ub(i1:i2,j1:j2,k1:k2)  = un(i1:i2,j1:j2,k1:k2)
587         ENDIF
588         !
589      ENDIF
590      !
591   END SUBROUTINE updateu
592
593# endif
594
595   SUBROUTINE correct_u_bdy( tabres, i1, i2, j1, j2, k1, k2, n1, n2, before, nb, ndir )
596      !!----------------------------------------------------------------------
597      !!                   *** ROUTINE correct_u_bdy ***
598      !!----------------------------------------------------------------------
599      INTEGER                                     , INTENT(in   ) :: i1, i2, j1, j2, k1, k2, n1, n2
600      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) :: tabres
601      LOGICAL                                     , INTENT(in   ) :: before
602      INTEGER                                     , INTENT(in   ) :: nb, ndir
603      !!
604      LOGICAL :: western_side, eastern_side 
605      INTEGER ::   jj, jk
606      REAL(wp)::   zcor
607      !!----------------------------------------------------------------------
608      !
609      IF( .NOT.before ) THEN
610         !
611         western_side  = (nb == 1).AND.(ndir == 1)
612         eastern_side  = (nb == 1).AND.(ndir == 2)
613         !
614         IF (western_side) THEN
615            DO jj=j1,j2
616               zcor = un_b(i1-1,jj) * hu_a(i1-1,jj) * r1_hu_n(i1-1,jj) - un_b(i1-1,jj)
617               un_b(i1-1,jj) = un_b(i1-1,jj) + zcor
618               DO jk=1,jpkm1
619                  un(i1-1,jj,jk) = un(i1-1,jj,jk) + zcor * umask(i1-1,jj,jk)
620               END DO
621            END DO
622         ENDIF
623         !
624         IF (eastern_side) THEN
625            DO jj=j1,j2
626               zcor = un_b(i2+1,jj) * hu_a(i2+1,jj) * r1_hu_n(i2+1,jj) - un_b(i2+1,jj)
627               un_b(i2+1,jj) = un_b(i2+1,jj) + zcor
628               DO jk=1,jpkm1
629                  un(i2+1,jj,jk) = un(i2+1,jj,jk) + zcor * umask(i2+1,jj,jk)
630               END DO
631            END DO
632         ENDIF
633         !
634      ENDIF
635      !
636   END SUBROUTINE correct_u_bdy
637
638# if  defined key_vertical
639
640   SUBROUTINE updatev( tabres, i1, i2, j1, j2, k1, k2, n1, n2, before )
641      !!----------------------------------------------------------------------
642      !!                   *** ROUTINE updatev ***
643      !!----------------------------------------------------------------------
644      INTEGER                                     , INTENT(in   ) :: i1, i2, j1, j2, k1, k2, n1, n2
645      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) :: tabres
646      LOGICAL                                     , INTENT(in   ) :: before
647      !
648      INTEGER  ::   ji, jj, jk
649      REAL(wp) ::   zrhox
650! VERTICAL REFINEMENT BEGIN
651      REAL(wp), DIMENSION(i1:i2,j1:j2,1:jpk) :: tabres_child
652      REAL(wp) :: h_in(k1:k2)
653      REAL(wp) :: h_out(1:jpk)
654      INTEGER :: N_in, N_out
655      REAL(wp) :: h_diff, excess, thick
656      REAL(wp) :: tabin(k1:k2)
657! VERTICAL REFINEMENT END
658      !!----------------------------------------------------------------------     
659      !
660      IF( before ) THEN
661         zrhox = Agrif_Rhox()
662         AGRIF_SpecialValue = -999._wp
663         DO jk=k1,k2
664            DO jj=j1,j2
665               DO ji=i1,i2
666                  tabres(ji,jj,jk,1) = zrhox * e1v(ji,jj) * e3v_n(ji,jj,jk) * vmask(ji,jj,jk) * vn(ji,jj,jk) &
667                                       + (vmask(ji,jj,jk)-1)*999._wp
668                  tabres(ji,jj,jk,2) = vmask(ji,jj,jk) * zrhox * e1v(ji,jj) * e3v_n(ji,jj,jk) &
669                                       + (vmask(ji,jj,jk)-1)*999._wp
670               END DO
671            END DO
672         END DO
673      ELSE
674         tabres_child(:,:,:) = 0.
675         AGRIF_SpecialValue = 0._wp
676         DO jj=j1,j2
677            DO ji=i1,i2
678               N_in = 0
679               DO jk=k1,k2
680                  IF (tabres(ji,jj,jk,2) < -900) EXIT
681                  N_in = N_in + 1
682                  tabin(jk) = tabres(ji,jj,jk,1)/tabres(ji,jj,jk,2)
683                  h_in(N_in) = tabres(ji,jj,jk,2)/e1v(ji,jj)
684               END DO
685               N_out = 0
686               DO jk=1,jpk
687                  IF (vmask(ji,jj,jk) == 0) EXIT
688                  N_out = N_out + 1
689                  h_out(N_out) = e3v_n(ji,jj,jk)
690               END DO
691               IF (N_in * N_out > 0) THEN
692                  h_diff = sum(h_out(1:N_out))-sum(h_in(1:N_in))
693                  IF (h_diff < -1.e-4) then
694!Even if bathy at T points match it's possible for the U points to be deeper in the child grid.
695!In this case we need to move transport from the child grid cells below bed of parent grid into the bottom cell.
696                     excess = 0._wp
697                     DO jk=N_in,1,-1
698                        thick = MIN(-1*h_diff, h_in(jk))
699                        excess = excess + tabin(jk)*thick*e2u(ji,jj)
700                        tabin(jk) = tabin(jk)*(1. - thick/h_in(jk))
701                        h_diff = h_diff + thick
702                        IF ( h_diff == 0) THEN
703                           N_in = jk
704                           h_in(jk) = h_in(jk) - thick
705                           EXIT
706                        ENDIF
707                     END DO
708                  ENDIF
709                  CALL reconstructandremap(tabin(1:N_in),h_in(1:N_in),tabres_child(ji,jj,1:N_out),h_out(1:N_out),N_in,N_out)
710                  tabres_child(ji,jj,N_out) = tabres_child(ji,jj,N_out) + excess/(e1v(ji,jj)*h_out(N_out))
711               ENDIF
712            END DO
713         END DO
714
715         DO jk=1,jpk
716            DO jj=j1,j2
717               DO ji=i1,i2
718                  !
719                  IF( .NOT.( lk_agrif_fstep .AND. l_1st_euler ) ) THEN  ! Add asselin part
720!!gm                  IF( .NOT.(lk_agrif_fstep.AND.(neuler==0)) ) THEN  ! Add asselin part
721                     vb(ji,jj,jk) = vb(ji,jj,jk) + rn_atfp * ( tabres_child(ji,jj,jk) - vn(ji,jj,jk) ) * vmask(ji,jj,jk)
722                  ENDIF
723                  !
724                  vn(ji,jj,jk) = tabres_child(ji,jj,jk) * vmask(ji,jj,jk)
725               END DO
726            END DO
727         END DO
728      ENDIF
729      !
730   END SUBROUTINE updatev
731
732# else
733
734   SUBROUTINE updatev( tabres, i1, i2, j1, j2, k1, k2, n1, n2, before)
735      !!----------------------------------------------------------------------
736      !!                   *** ROUTINE updatev ***
737      !!----------------------------------------------------------------------
738      INTEGER                                     , INTENT(in   ) :: i1, i2, j1, j2, k1, k2, n1, n2
739      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) :: tabres
740      LOGICAL                                     , INTENT(in   ) :: before
741      !
742      INTEGER  :: ji, jj, jk
743      REAL(wp) :: zrhox, zvb, zvnu, zvno
744      !!----------------------------------------------------------------------     
745      !
746      IF (before) THEN
747         zrhox = Agrif_Rhox()
748         DO jk=k1,k2
749            DO jj=j1,j2
750               DO ji=i1,i2
751                  tabres(ji,jj,jk,1) = zrhox * e1v(ji,jj) * e3v_n(ji,jj,jk) * vn(ji,jj,jk)
752               END DO
753            END DO
754         END DO
755      ELSE
756         DO jk=k1,k2
757            DO jj=j1,j2
758               DO ji=i1,i2
759                  tabres(ji,jj,jk,1) = tabres(ji,jj,jk,1) * r1_e1v(ji,jj)
760                  !
761                  IF ( .NOT.( lk_agrif_fstep .AND. l_1st_euler ) ) THEN    ! Add asselin part
762!!gm                  IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN   ! Add asselin part
763                     zvb  = vb(ji,jj,jk) * e3v_b(ji,jj,jk) ! fse3t_b prior update should be used
764                     zvno = vn(ji,jj,jk) * e3v_a(ji,jj,jk)
765                     zvnu = tabres(ji,jj,jk,1)
766                     vb(ji,jj,jk) = ( zvb + rn_atfp * ( zvnu - zvno) ) / e3v_b(ji,jj,jk) * vmask(ji,jj,jk)
767                  ENDIF
768                  !
769                  vn(ji,jj,jk) = tabres(ji,jj,jk,1) / e3v_n(ji,jj,jk) * vmask(ji,jj,jk)
770               END DO
771            END DO
772         END DO
773         !
774         IF ( l_1st_euler .AND. Agrif_Nb_Step() == 0 ) THEN
775!!gm         IF  ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN
776            vb(i1:i2,j1:j2,k1:k2)  = vn(i1:i2,j1:j2,k1:k2)
777         ENDIF
778         !
779      ENDIF
780      !
781   END SUBROUTINE updatev
782
783# endif
784
785   SUBROUTINE correct_v_bdy( tabres, i1, i2, j1, j2, k1, k2, n1, n2, before, nb, ndir )
786      !!----------------------------------------------------------------------
787      !!                   *** ROUTINE correct_v_bdy ***
788      !!----------------------------------------------------------------------
789      INTEGER                                     , INTENT(in   ) :: i1, i2, j1, j2, k1, k2, n1, n2
790      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) :: tabres
791      LOGICAL                                     , INTENT(in   ) :: before
792      INTEGER                                     , INTENT(in)    :: nb, ndir
793      !!
794      LOGICAL :: southern_side, northern_side 
795      !
796      INTEGER  :: ji, jk
797      REAL(wp) :: zcor
798      !!----------------------------------------------------------------------
799      !
800      IF( .NOT.before ) THEN
801         !
802         southern_side = (nb == 2).AND.(ndir == 1)
803         northern_side = (nb == 2).AND.(ndir == 2)
804         !
805         IF (southern_side) THEN
806            DO ji=i1,i2
807               zcor = vn_b(ji,j1-1) * hv_a(ji,j1-1) * r1_hv_n(ji,j1-1) - vn_b(ji,j1-1)
808               vn_b(ji,j1-1) = vn_b(ji,j1-1) + zcor
809               DO jk=1,jpkm1
810                  vn(ji,j1-1,jk) = vn(ji,j1-1,jk) + zcor * vmask(ji,j1-1,jk)
811               END DO
812            END DO
813         ENDIF
814         !
815         IF (northern_side) THEN
816            DO ji=i1,i2
817               zcor = vn_b(ji,j2+1) * hv_a(ji,j2+1) * r1_hv_n(ji,j2+1) - vn_b(ji,j2+1)
818               vn_b(ji,j2+1) = vn_b(ji,j2+1) + zcor
819               DO jk=1,jpkm1
820                  vn(ji,j2+1,jk) = vn(ji,j2+1,jk) + zcor * vmask(ji,j2+1,jk)
821               END DO
822            END DO
823         ENDIF
824         !
825      ENDIF
826      !
827   END SUBROUTINE correct_v_bdy
828
829
830   SUBROUTINE updateu2d( tabres, i1, i2, j1, j2, before )
831      !!----------------------------------------------------------------------
832      !!                   *** ROUTINE updateu2d ***
833      !!----------------------------------------------------------------------
834      INTEGER                         , INTENT(in   ) ::   i1, i2, j1, j2
835      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) ::   tabres
836      LOGICAL                         , INTENT(in   ) ::   before
837      !
838      INTEGER  :: ji, jj, jk
839      REAL(wp) :: zrhoy
840      REAL(wp) :: zcorr
841      !!----------------------------------------------------------------------
842      !
843      IF( before ) THEN
844         zrhoy = Agrif_Rhoy()
845         DO jj=j1,j2
846            DO ji=i1,i2
847               tabres(ji,jj) = zrhoy * un_b(ji,jj) * hu_n(ji,jj) * e2u(ji,jj)
848            END DO
849         END DO
850      ELSE
851         DO jj=j1,j2
852            DO ji=i1,i2
853               tabres(ji,jj) =  tabres(ji,jj) * r1_e2u(ji,jj) 
854               !   
855               ! Update "now" 3d velocities:
856               spgu(ji,jj) = 0._wp
857               DO jk=1,jpkm1
858                  spgu(ji,jj) = spgu(ji,jj) + e3u_n(ji,jj,jk) * un(ji,jj,jk)
859               END DO
860               !
861               zcorr = (tabres(ji,jj) - spgu(ji,jj)) * r1_hu_n(ji,jj)
862               DO jk=1,jpkm1             
863                  un(ji,jj,jk) = un(ji,jj,jk) + zcorr * umask(ji,jj,jk)           
864               END DO
865               !
866               ! Update barotropic velocities:
867               IF ( .NOT.ln_dynspg_ts .OR. (ln_dynspg_ts.AND.(.NOT.ln_bt_fw)) ) THEN
868                  IF ( .NOT.( lk_agrif_fstep .AND. l_1st_euler ) ) THEN    ! Add asselin part
869!!gm                  IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN ! Add asselin part
870                     zcorr = (tabres(ji,jj) - un_b(ji,jj) * hu_a(ji,jj)) * r1_hu_b(ji,jj)
871                     ub_b(ji,jj) = ub_b(ji,jj) + rn_atfp * zcorr * umask(ji,jj,1)
872                  END IF
873               ENDIF   
874               un_b(ji,jj) = tabres(ji,jj) * r1_hu_n(ji,jj) * umask(ji,jj,1)
875               !       
876               ! Correct "before" velocities to hold correct bt component:
877               spgu(ji,jj) = 0.e0
878               DO jk=1,jpkm1
879                  spgu(ji,jj) = spgu(ji,jj) + e3u_b(ji,jj,jk) * ub(ji,jj,jk)
880               END DO
881               !
882               zcorr = ub_b(ji,jj) - spgu(ji,jj) * r1_hu_b(ji,jj)
883               DO jk=1,jpkm1             
884                  ub(ji,jj,jk) = ub(ji,jj,jk) + zcorr * umask(ji,jj,jk)           
885               END DO
886               !
887            END DO
888         END DO
889         !
890         IF ( l_1st_euler .AND. Agrif_Nb_Step() == 0 ) THEN
891!!gm         IF  ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN
892            ub_b(i1:i2,j1:j2)  = un_b(i1:i2,j1:j2)
893         ENDIF
894      ENDIF
895      !
896   END SUBROUTINE updateu2d
897
898
899   SUBROUTINE updatev2d( tabres, i1, i2, j1, j2, before )
900      !!----------------------------------------------------------------------
901      !!                   *** ROUTINE updatev2d ***
902      !!----------------------------------------------------------------------
903      INTEGER                         , INTENT(in   ) ::   i1, i2, j1, j2
904      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) ::   tabres
905      LOGICAL                         , INTENT(in   ) ::   before
906      !
907      INTEGER  :: ji, jj, jk
908      REAL(wp) :: zrhox, zcorr
909      !!----------------------------------------------------------------------
910      !
911      IF( before ) THEN
912         zrhox = Agrif_Rhox()
913         DO jj=j1,j2
914            DO ji=i1,i2
915               tabres(ji,jj) = zrhox * vn_b(ji,jj) * hv_n(ji,jj) * e1v(ji,jj) 
916            END DO
917         END DO
918      ELSE
919         DO jj=j1,j2
920            DO ji=i1,i2
921               tabres(ji,jj) =  tabres(ji,jj) * r1_e1v(ji,jj) 
922               !   
923               ! Update "now" 3d velocities:
924               spgv(ji,jj) = 0.e0
925               DO jk=1,jpkm1
926                  spgv(ji,jj) = spgv(ji,jj) + e3v_n(ji,jj,jk) * vn(ji,jj,jk)
927               END DO
928               !
929               zcorr = (tabres(ji,jj) - spgv(ji,jj)) * r1_hv_n(ji,jj)
930               DO jk=1,jpkm1             
931                  vn(ji,jj,jk) = vn(ji,jj,jk) + zcorr * vmask(ji,jj,jk)           
932               END DO
933               !
934               ! Update barotropic velocities:
935               IF ( .NOT.ln_dynspg_ts .OR. ( ln_dynspg_ts .AND. .NOT.ln_bt_fw ) ) THEN
936                  IF ( .NOT.( lk_agrif_fstep .AND. l_1st_euler ) ) THEN    ! Add asselin part
937!!gm                  IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN ! Add asselin part
938                     zcorr = (tabres(ji,jj) - vn_b(ji,jj) * hv_a(ji,jj)) * r1_hv_b(ji,jj)
939                     vb_b(ji,jj) = vb_b(ji,jj) + rn_atfp * zcorr * vmask(ji,jj,1)
940                  END IF
941               ENDIF             
942               vn_b(ji,jj) = tabres(ji,jj) * r1_hv_n(ji,jj) * vmask(ji,jj,1)
943               !       
944               ! Correct "before" velocities to hold correct bt component:
945               spgv(ji,jj) = 0.e0
946               DO jk=1,jpkm1
947                  spgv(ji,jj) = spgv(ji,jj) + e3v_b(ji,jj,jk) * vb(ji,jj,jk)
948               END DO
949               !
950               zcorr = vb_b(ji,jj) - spgv(ji,jj) * r1_hv_b(ji,jj)
951               DO jk=1,jpkm1             
952                  vb(ji,jj,jk) = vb(ji,jj,jk) + zcorr * vmask(ji,jj,jk)           
953               END DO
954               !
955            END DO
956         END DO
957         !
958         IF ( l_1st_euler .AND. Agrif_Nb_Step() == 0 ) THEN
959!!gm         IF  ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN
960            vb_b(i1:i2,j1:j2)  = vn_b(i1:i2,j1:j2)
961         ENDIF
962         !
963      ENDIF
964      !
965   END SUBROUTINE updatev2d
966
967
968   SUBROUTINE updateSSH( tabres, i1, i2, j1, j2, before )
969      !!----------------------------------------------------------------------
970      !!                   *** ROUTINE updateSSH ***
971      !!----------------------------------------------------------------------
972      INTEGER                         , INTENT(in   ) ::   i1, i2, j1, j2
973      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) ::   tabres
974      LOGICAL                         , INTENT(in   ) ::   before
975      !
976      INTEGER :: ji, jj
977      !!----------------------------------------------------------------------
978      !
979      IF( before ) THEN
980         DO jj = j1, j2
981            DO ji = i1, i2
982               tabres(ji,jj) = ssh(ji,jj,Nnn)
983            END DO
984         END DO
985      ELSE
986         IF ( .NOT.( lk_agrif_fstep .AND. l_1st_euler ) ) THEN
987!!gm         IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN
988            DO jj = j1, j2
989               DO ji = i1, i2
990                  ssh(ji,jj,Nbb) = ssh(ji,jj,Nbb) + rn_atfp * ( tabres(ji,jj) - ssh(ji,jj,Nnn) ) * tmask(ji,jj,1)
991               END DO
992            END DO
993         ENDIF
994         !
995         DO jj = j1, j2
996            DO ji = i1, i2
997               ssh(ji,jj,Nnn) = tabres(ji,jj) * tmask(ji,jj,1)
998            END DO
999         END DO
1000         !
1001         IF ( l_1st_euler .AND. Agrif_Nb_Step() == 0 ) THEN
1002!!gm         IF  ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN
1003            ssh(i1:i2,j1:j2,Nbb)  = ssh(i1:i2,j1:j2,Nnn)
1004         ENDIF
1005         !
1006      ENDIF
1007      !
1008   END SUBROUTINE updateSSH
1009
1010
1011   SUBROUTINE updateub2b( tabres, i1, i2, j1, j2, before )
1012      !!----------------------------------------------------------------------
1013      !!                      *** ROUTINE updateub2b ***
1014      !!----------------------------------------------------------------------
1015      INTEGER                            , INTENT(in) ::   i1, i2, j1, j2
1016      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) ::   tabres
1017      LOGICAL                            , INTENT(in) ::   before
1018      !!
1019      INTEGER :: ji, jj
1020      REAL(wp) :: zrhoy, za1, zcor
1021      !!---------------------------------------------
1022      !
1023      IF (before) THEN
1024         zrhoy = Agrif_Rhoy()
1025         DO jj=j1,j2
1026            DO ji=i1,i2
1027               tabres(ji,jj) = ub2_i_b(ji,jj) * e2u(ji,jj)
1028            END DO
1029         END DO
1030         tabres = zrhoy * tabres
1031      ELSE
1032         !
1033         tabres(i1:i2,j1:j2) = tabres(i1:i2,j1:j2) * r1_e2u(i1:i2,j1:j2)
1034         !
1035         za1 = 1._wp / REAL(Agrif_rhot(), wp)
1036         DO jj=j1,j2
1037            DO ji=i1,i2
1038               zcor=tabres(ji,jj) - ub2_b(ji,jj)
1039               ! Update time integrated fluxes also in case of multiply nested grids:
1040               ub2_i_b(ji,jj) = ub2_i_b(ji,jj) + za1 * zcor 
1041               ! Update corrective fluxes:
1042               un_bf(ji,jj)  = un_bf(ji,jj) + zcor
1043               ! Update half step back fluxes:
1044               ub2_b(ji,jj) = tabres(ji,jj)
1045            END DO
1046         END DO
1047      ENDIF
1048      !
1049   END SUBROUTINE updateub2b
1050
1051
1052   SUBROUTINE reflux_sshu( tabres, i1, i2, j1, j2, before, nb, ndir )
1053      !!----------------------------------------------------------------------
1054      !!                   *** ROUTINE reflux_sshu ***
1055      !!----------------------------------------------------------------------
1056      INTEGER                         , INTENT(in   ) ::   i1, i2, j1, j2
1057      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) ::   tabres
1058      LOGICAL                         , INTENT(in   ) ::   before
1059      INTEGER                         , INTENT(in   ) ::   nb, ndir
1060      !
1061      LOGICAL ::   western_side, eastern_side 
1062      INTEGER ::   ji, jj
1063      REAL(wp)::   zrhoy, za1, zcor
1064      !!----------------------------------------------------------------------
1065      !
1066      IF (before) THEN
1067         zrhoy = Agrif_Rhoy()
1068         DO jj=j1,j2
1069            DO ji=i1,i2
1070               tabres(ji,jj) = ub2_i_b(ji,jj) * e2u(ji,jj)
1071            END DO
1072         END DO
1073         tabres = zrhoy * tabres
1074      ELSE
1075         !
1076         tabres(i1:i2,j1:j2) = tabres(i1:i2,j1:j2) * r1_e2u(i1:i2,j1:j2)
1077         !
1078         western_side  = (nb == 1).AND.(ndir == 1)
1079         eastern_side  = (nb == 1).AND.(ndir == 2)
1080         !
1081         IF ( western_side ) THEN
1082            DO jj=j1,j2
1083               zcor = rn_Dt * r1_e1e2t(i1  ,jj) * e2u(i1,jj) * (ub2_b(i1,jj)-tabres(i1,jj)) 
1084               ssh(i1  ,jj,Nnn) = ssh(i1  ,jj,Nnn) + zcor
1085               IF ( .NOT.( lk_agrif_fstep .AND. l_1st_euler ) )   ssh(i1  ,jj,Nbb) = ssh(i1  ,jj,Nbb) + rn_atfp * zcor
1086!!gm               IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) sshb(i1  ,jj) = sshb(i1  ,jj) + rn_atfp * zcor
1087            END DO
1088         ENDIF
1089         IF ( eastern_side ) THEN
1090            DO jj=j1,j2
1091               zcor = - rn_Dt * r1_e1e2t(i2+1,jj) * e2u(i2,jj) * (ub2_b(i2,jj)-tabres(i2,jj))
1092               ssh(i2+1,jj,Nnn) = ssh(i2+1,jj,Nnn) + zcor
1093               IF (.NOT.( lk_agrif_fstep .AND. l_1st_euler ) )   ssh(i2+1,jj,Nbb) = ssh(i2+1,jj,Nbb) + rn_atfp * zcor
1094!!gm               IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) sshb(i2+1,jj) = sshb(i2+1,jj) + rn_atfp * zcor
1095            END DO
1096         ENDIF
1097         !
1098      ENDIF
1099      !
1100   END SUBROUTINE reflux_sshu
1101
1102
1103   SUBROUTINE updatevb2b( tabres, i1, i2, j1, j2, before )
1104      !!----------------------------------------------------------------------
1105      !!                    *** ROUTINE updatevb2b ***
1106      !!----------------------------------------------------------------------
1107      INTEGER                         , INTENT(in   ) ::   i1, i2, j1, j2
1108      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) ::   tabres
1109      LOGICAL                         , INTENT(in   ) ::   before
1110      !
1111      INTEGER :: ji, jj
1112      REAL(wp) :: zrhox, za1, zcor
1113      !!---------------------------------------------------------------------
1114      !
1115      IF( before ) THEN
1116         zrhox = Agrif_Rhox()
1117         DO jj=j1,j2
1118            DO ji=i1,i2
1119               tabres(ji,jj) = vb2_i_b(ji,jj) * e1v(ji,jj) 
1120            END DO
1121         END DO
1122         tabres = zrhox * tabres
1123      ELSE
1124         !
1125         tabres(i1:i2,j1:j2) = tabres(i1:i2,j1:j2) * r1_e1v(i1:i2,j1:j2)
1126         !
1127         za1 = 1._wp / REAL(Agrif_rhot(), wp)
1128         DO jj=j1,j2
1129            DO ji=i1,i2
1130               zcor=tabres(ji,jj) - vb2_b(ji,jj)
1131               ! Update time integrated fluxes also in case of multiply nested grids:
1132               vb2_i_b(ji,jj) = vb2_i_b(ji,jj) + za1 * zcor 
1133               ! Update corrective fluxes:
1134               vn_bf(ji,jj)  = vn_bf(ji,jj) + zcor
1135               ! Update half step back fluxes:
1136               vb2_b(ji,jj) = tabres(ji,jj)
1137            END DO
1138         END DO
1139      ENDIF
1140      !
1141   END SUBROUTINE updatevb2b
1142
1143
1144   SUBROUTINE reflux_sshv( tabres, i1, i2, j1, j2, before, nb, ndir )
1145      !!----------------------------------------------------------------------
1146      !!                   *** ROUTINE reflux_sshv ***
1147      !!----------------------------------------------------------------------
1148      INTEGER                         , INTENT(in   ) ::   i1, i2, j1, j2
1149      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) ::   tabres
1150      LOGICAL                         , INTENT(in   ) ::   before
1151      INTEGER                         , INTENT(in   ) ::   nb, ndir
1152      !!
1153      LOGICAL :: southern_side, northern_side 
1154      INTEGER :: ji, jj
1155      REAL(wp) :: zrhox, za1, zcor
1156      !!----------------------------------------------------------------------
1157      !
1158      IF (before) THEN
1159         zrhox = Agrif_Rhox()
1160         DO jj=j1,j2
1161            DO ji=i1,i2
1162               tabres(ji,jj) = vb2_i_b(ji,jj) * e1v(ji,jj) 
1163            END DO
1164         END DO
1165         tabres = zrhox * tabres
1166      ELSE
1167         !
1168         tabres(i1:i2,j1:j2) = tabres(i1:i2,j1:j2) * r1_e1v(i1:i2,j1:j2)
1169         !
1170         southern_side = (nb == 2).AND.(ndir == 1)
1171         northern_side = (nb == 2).AND.(ndir == 2)
1172         !
1173         IF (southern_side) THEN
1174            DO ji=i1,i2
1175               zcor = rn_Dt * r1_e1e2t(ji,j1  ) * e1v(ji,j1  ) * ( vb2_b(ji,j1)-tabres(ji,j1) )
1176               ssh(ji,j1  ,Nnn) = ssh(ji,j1  ,Nnn) + zcor
1177               IF ( .NOT.( lk_agrif_fstep .AND. l_1st_euler ) )   ssh(ji,j1  ,Nbb) = ssh(ji,j1,Nbb) + rn_atfp * zcor
1178!!gm               IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) sshb(ji,j1  ) = sshb(ji,j1) + rn_atfp * zcor
1179            END DO
1180         ENDIF
1181         IF (northern_side) THEN               
1182            DO ji=i1,i2
1183               zcor = - rn_Dt * r1_e1e2t(ji,j2+1) * e1v(ji,j2  ) * ( vb2_b(ji,j2)-tabres(ji,j2) )
1184               ssh(ji,j2+1,Nnn) = ssh(ji,j2+1,Nnn) + zcor
1185               IF ( .NOT.( lk_agrif_fstep .AND. l_1st_euler ) )   ssh(ji,j2+1,Nbb) = ssh(ji,j2+1,Nbb) + rn_atfp * zcor
1186!!gm               IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) sshb(ji,j2+1) = sshb(ji,j2+1) + rn_atfp * zcor
1187            END DO
1188         ENDIF
1189         !
1190      ENDIF
1191      !
1192   END SUBROUTINE reflux_sshv
1193
1194
1195   SUBROUTINE update_scales( tabres, i1, i2, j1, j2, k1, k2, n1,n2, before )
1196      !!----------------------------------------------------------------------
1197      !!                      *** ROUTINE updateT ***
1198      !
1199      ! ====>>>>>>>>>>    currently not used
1200      !
1201      !!----------------------------------------------------------------------
1202      INTEGER                                    , INTENT(in   ) ::   i1, i2, j1, j2, k1, k2, n1, n2
1203      REAL(wp),DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) ::   tabres
1204      LOGICAL                                    , INTENT(in   ) ::   before
1205      !!
1206      INTEGER :: ji,jj,jk
1207      REAL(wp) :: ztemp
1208      !!----------------------------------------------------------------------
1209
1210      IF (before) THEN
1211         DO jk=k1,k2
1212            DO jj=j1,j2
1213               DO ji=i1,i2
1214                  tabres(ji,jj,jk,1) = e1t(ji,jj)*e2t(ji,jj)*tmask(ji,jj,jk)
1215                  tabres(ji,jj,jk,2) = e1t(ji,jj)*tmask(ji,jj,jk)
1216                  tabres(ji,jj,jk,3) = e2t(ji,jj)*tmask(ji,jj,jk)
1217               END DO
1218            END DO
1219         END DO
1220         tabres(:,:,:,1)=tabres(:,:,:,1)*Agrif_Rhox()*Agrif_Rhoy()
1221         tabres(:,:,:,2)=tabres(:,:,:,2)*Agrif_Rhox()
1222         tabres(:,:,:,3)=tabres(:,:,:,3)*Agrif_Rhoy()
1223      ELSE
1224         DO jk=k1,k2
1225            DO jj=j1,j2
1226               DO ji=i1,i2
1227                  IF( tabres(ji,jj,jk,1) .NE. 0. ) THEN
1228                     print *,'VAL = ',ji,jj,jk,tabres(ji,jj,jk,1),e1t(ji,jj)*e2t(ji,jj)*tmask(ji,jj,jk)
1229                     print *,'VAL2 = ',ji,jj,jk,tabres(ji,jj,jk,2),e1t(ji,jj)*tmask(ji,jj,jk)
1230                     print *,'VAL3 = ',ji,jj,jk,tabres(ji,jj,jk,3),e2t(ji,jj)*tmask(ji,jj,jk)
1231                     ztemp = sqrt(tabres(ji,jj,jk,1)/(tabres(ji,jj,jk,2)*tabres(ji,jj,jk,3)))
1232                     print *,'CORR = ',ztemp-1.
1233                     print *,'NEW VALS = ',tabres(ji,jj,jk,2)*ztemp,tabres(ji,jj,jk,3)*ztemp, &
1234                           tabres(ji,jj,jk,2)*ztemp*tabres(ji,jj,jk,3)*ztemp
1235                     e1t(ji,jj) = tabres(ji,jj,jk,2)*ztemp
1236                     e2t(ji,jj) = tabres(ji,jj,jk,3)*ztemp
1237                  END IF
1238               END DO
1239            END DO
1240         END DO
1241      ENDIF
1242      !
1243   END SUBROUTINE update_scales
1244
1245
1246   SUBROUTINE updateEN( ptab, i1, i2, j1, j2, k1, k2, before )
1247      !!----------------------------------------------------------------------
1248      !!                      *** ROUTINE updateen ***
1249      !!----------------------------------------------------------------------
1250      INTEGER                               , INTENT(in   ) ::   i1, i2, j1, j2, k1, k2
1251      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) ::   ptab
1252      LOGICAL                               , INTENT(in   ) ::   before
1253      !!----------------------------------------------------------------------
1254      !
1255      IF( before ) THEN
1256         ptab (i1:i2,j1:j2,k1:k2) = en(i1:i2,j1:j2,k1:k2)
1257      ELSE
1258         en(i1:i2,j1:j2,k1:k2) = ptab (i1:i2,j1:j2,k1:k2) 
1259      ENDIF
1260      !
1261   END SUBROUTINE updateEN
1262
1263
1264   SUBROUTINE updateAVT( ptab, i1, i2, j1, j2, k1, k2, before )
1265      !!----------------------------------------------------------------------
1266      !!                      *** ROUTINE updateavt ***
1267      !!----------------------------------------------------------------------
1268      INTEGER                               , INTENT(in   ) ::   i1, i2, j1, j2, k1, k2
1269      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) ::   ptab
1270      LOGICAL                               , INTENT(in   ) ::   before
1271      !!----------------------------------------------------------------------
1272      !
1273      IF( before ) THEN   ;   ptab (i1:i2,j1:j2,k1:k2) = avt_k(i1:i2,j1:j2,k1:k2)
1274      ELSE                ;   avt_k(i1:i2,j1:j2,k1:k2) = ptab (i1:i2,j1:j2,k1:k2) 
1275      ENDIF
1276      !
1277   END SUBROUTINE updateAVT
1278
1279
1280   SUBROUTINE updateAVM( ptab, i1, i2, j1, j2, k1, k2, before )
1281      !!----------------------------------------------------------------------
1282      !!                      *** ROUTINE updateavm ***
1283      !!----------------------------------------------------------------------
1284      INTEGER                               , INTENT(in   ) ::   i1, i2, j1, j2, k1, k2
1285      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) ::   ptab
1286      LOGICAL                               , INTENT(in   ) ::   before
1287      !!----------------------------------------------------------------------
1288      !
1289      IF( before ) THEN   ;   ptab (i1:i2,j1:j2,k1:k2) = avm_k(i1:i2,j1:j2,k1:k2)
1290      ELSE                ;   avm_k(i1:i2,j1:j2,k1:k2) = ptab (i1:i2,j1:j2,k1:k2) 
1291      ENDIF
1292      !
1293   END SUBROUTINE updateAVM
1294
1295
1296   SUBROUTINE updatee3t(ptab_dum, i1, i2, j1, j2, k1, k2, before )
1297      !!----------------------------------------------------------------------
1298      !!                   *** ROUTINE updatee3t ***
1299      !!----------------------------------------------------------------------
1300      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2) :: ptab_dum
1301      INTEGER, INTENT(in) :: i1, i2, j1, j2, k1, k2
1302      LOGICAL, INTENT(in) :: before
1303      !
1304      REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: ptab
1305      INTEGER :: ji,jj,jk
1306      REAL(wp) :: zcoef
1307      !!---------------------------------------------
1308      !
1309      IF (.NOT.before) THEN
1310         !
1311         ALLOCATE( ptab(i1:i2,j1:j2,1:jpk) ) 
1312         !
1313         ! Update e3t from ssh (z* case only)
1314         DO jk = 1, jpkm1
1315            DO jj = j1, j2
1316               DO ji = i1, i2
1317                  ptab(ji,jj,jk) = e3t_0(ji,jj,jk) * (1._wp + ssh(ji,jj,Nnn) * r1_ht_0(ji,jj) *tmask(ji,jj,jk) )
1318               END DO
1319            END DO
1320         END DO
1321         !
1322         ! 1) Updates at BEFORE time step:
1323         ! -------------------------------
1324         !
1325         ! Save "old" scale factor (prior update) for subsequent asselin correction
1326         ! of prognostic variables
1327         e3t_a(i1:i2,j1:j2,1:jpkm1) = e3t_n(i1:i2,j1:j2,1:jpkm1)
1328
1329         ! One should also save e3t_b, but lacking of workspace...
1330!         hdivn(i1:i2,j1:j2,1:jpkm1)   = e3t_b(i1:i2,j1:j2,1:jpkm1)
1331
1332         IF ( .NOT.( lk_agrif_fstep .AND. l_1st_euler ) ) THEN
1333!!gm         IF (.NOT.(lk_agrif_fstep.AND.(neuler==0) )) THEN
1334            DO jk = 1, jpkm1
1335               DO jj=j1,j2
1336                  DO ji=i1,i2
1337                     e3t_b(ji,jj,jk) =  e3t_b(ji,jj,jk) + rn_atfp * ( ptab(ji,jj,jk) - e3t_n(ji,jj,jk) )
1338                  END DO
1339               END DO
1340            END DO
1341            !
1342            e3w_b  (i1:i2,j1:j2,1) = e3w_0(i1:i2,j1:j2,1) + e3t_b(i1:i2,j1:j2,1) - e3t_0(i1:i2,j1:j2,1)
1343            gdepw_b(i1:i2,j1:j2,1) = 0.0_wp
1344            gdept_b(i1:i2,j1:j2,1) = 0.5_wp * e3w_b(i1:i2,j1:j2,1)
1345            !
1346            DO jk = 2, jpk
1347               DO jj = j1,j2
1348                  DO ji = i1,i2           
1349                     zcoef = (tmask(ji,jj,jk) - wmask(ji,jj,jk))
1350                     e3w_b(ji,jj,jk)  = e3w_0(ji,jj,jk) + ( 1.0_wp - 0.5_wp * tmask(ji,jj,jk) ) *        & 
1351                     &                                        ( e3t_b(ji,jj,jk-1) - e3t_0(ji,jj,jk-1) )  &
1352                     &                                  +            0.5_wp * tmask(ji,jj,jk)   *        &
1353                     &                                        ( e3t_b(ji,jj,jk  ) - e3t_0(ji,jj,jk  ) )
1354                     gdepw_b(ji,jj,jk) = gdepw_b(ji,jj,jk-1) + e3t_b(ji,jj,jk-1)
1355                     gdept_b(ji,jj,jk) =      zcoef  * ( gdepw_b(ji,jj,jk  ) + 0.5 * e3w_b(ji,jj,jk))  &
1356                         &               + (1-zcoef) * ( gdept_b(ji,jj,jk-1) +       e3w_b(ji,jj,jk)) 
1357                  END DO
1358               END DO
1359            END DO
1360            !
1361         ENDIF       
1362         !
1363         ! 2) Updates at NOW time step:
1364         ! ----------------------------
1365         !
1366         ! Update vertical scale factor at T-points:
1367         e3t_n(i1:i2,j1:j2,1:jpkm1) = ptab(i1:i2,j1:j2,1:jpkm1)
1368         !
1369         ! Update total depth:
1370         ht_n(i1:i2,j1:j2) = 0._wp
1371         DO jk = 1, jpkm1
1372            ht_n(i1:i2,j1:j2) = ht_n(i1:i2,j1:j2) + e3t_n(i1:i2,j1:j2,jk) * tmask(i1:i2,j1:j2,jk)
1373         END DO
1374         !
1375         ! Update vertical scale factor at W-points and depths:
1376         e3w_n (i1:i2,j1:j2,1) = e3w_0(i1:i2,j1:j2,1) + e3t_n(i1:i2,j1:j2,1) - e3t_0(i1:i2,j1:j2,1)
1377         gdept_n(i1:i2,j1:j2,1) = 0.5_wp * e3w_n(i1:i2,j1:j2,1)
1378         gdepw_n(i1:i2,j1:j2,1) = 0.0_wp
1379         gde3w_n(i1:i2,j1:j2,1) = gdept_n(i1:i2,j1:j2,1) - (ht_n(i1:i2,j1:j2)-ht_0(i1:i2,j1:j2)) ! Last term in the rhs is ssh
1380         !
1381         DO jk = 2, jpk
1382            DO jj = j1,j2
1383               DO ji = i1,i2           
1384               zcoef = (tmask(ji,jj,jk) - wmask(ji,jj,jk))
1385               e3w_n(ji,jj,jk)  = e3w_0(ji,jj,jk) + ( 1.0_wp - 0.5_wp * tmask(ji,jj,jk) ) * ( e3t_n(ji,jj,jk-1) - e3t_0(ji,jj,jk-1) )   &
1386               &                                  +            0.5_wp * tmask(ji,jj,jk)   * ( e3t_n(ji,jj,jk  ) - e3t_0(ji,jj,jk  ) )
1387               gdepw_n(ji,jj,jk) = gdepw_n(ji,jj,jk-1) + e3t_n(ji,jj,jk-1)
1388               gdept_n(ji,jj,jk) =      zcoef  * ( gdepw_n(ji,jj,jk  ) + 0.5 * e3w_n(ji,jj,jk))  &
1389                   &               + (1-zcoef) * ( gdept_n(ji,jj,jk-1) +       e3w_n(ji,jj,jk)) 
1390               gde3w_n(ji,jj,jk) = gdept_n(ji,jj,jk) - (ht_n(ji,jj)-ht_0(ji,jj)) ! Last term in the rhs is ssh
1391               END DO
1392            END DO
1393         END DO
1394         !
1395         IF ( l_1st_euler .AND. Agrif_Nb_Step()==0 ) THEN
1396!!gm         IF  ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN
1397            e3t_b (i1:i2,j1:j2,1:jpk)  = e3t_n (i1:i2,j1:j2,1:jpk)
1398            e3w_b (i1:i2,j1:j2,1:jpk)  = e3w_n (i1:i2,j1:j2,1:jpk)
1399            gdepw_b(i1:i2,j1:j2,1:jpk) = gdepw_n(i1:i2,j1:j2,1:jpk)
1400            gdept_b(i1:i2,j1:j2,1:jpk) = gdept_n(i1:i2,j1:j2,1:jpk)
1401         ENDIF
1402         !
1403         DEALLOCATE(ptab)
1404      ENDIF
1405      !
1406   END SUBROUTINE updatee3t
1407
1408#else
1409   !!----------------------------------------------------------------------
1410   !!   Empty module                                          no AGRIF zoom
1411   !!----------------------------------------------------------------------
1412CONTAINS
1413   SUBROUTINE agrif_oce_update_empty
1414      WRITE(*,*)  'agrif_oce_update : You should not have seen this print! error?'
1415   END SUBROUTINE agrif_oce_update_empty
1416#endif
1417
1418   !!======================================================================
1419END MODULE agrif_oce_update
1420
Note: See TracBrowser for help on using the repository browser.