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

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

#1911 (ENHANCE-04): RK3 branch - step II.1 AGRIF compilation OK for time-level dimension on ssh

  • Property svn:keywords set to Id
File size: 57.5 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!      ua(:,:,:) = e3u_b(:,:,:)
261!      va(:,:,:) = e3v_b(:,:,:)
262      hu_a(:,:) = hu_n(:,:)
263      hv_a(:,:) = hv_n(:,:)
264
265      !                          !* NOW fields :
266      CALL ssh2e3_now                  ! set: ht , hu , hv , r1_hu, r1_hv
267      !                                !      e3t, e3w, e3u, e3uw, e3v, e3vw, e3f   (from 1 to jpkm1)
268      !                                !      gdept_n, gdepw_n, gde3w_n
269
270      !                          !* BEFORE fields :
271      IF (.NOT.( lk_agrif_fstep .AND. l_1st_euler ) ) THEN
272         !
273         CALL ssh2e3_before               ! set:      hu , hv , r1_hu, r1_hv
274         !                                !      e3t, e3w, e3u, e3uw, e3v, e3vw        (from 1 to jpkm1)
275         !
276      ENDIF
277      !
278   END SUBROUTINE dom_vvl_update_UVF
279
280# if defined key_vertical
281
282   SUBROUTINE updateTS( tabres, i1, i2, j1, j2, k1, k2, n1, n2, before )
283      !!----------------------------------------------------------------------
284      !!                   *** ROUTINE updateT ***
285      !!----------------------------------------------------------------------
286      INTEGER, INTENT(in) :: i1,i2,j1,j2,k1,k2,n1,n2
287      REAL(wp),DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) :: tabres
288      LOGICAL, INTENT(in) :: before
289      !!
290      INTEGER :: ji,jj,jk,jn
291      REAL(wp), DIMENSION(i1:i2,j1:j2,1:jpk,n1:n2) :: tabres_child
292      REAL(wp) :: h_in(k1:k2)
293      REAL(wp) :: h_out(1:jpk)
294      INTEGER  :: N_in, N_out
295      REAL(wp) :: zrho_xy, h_diff
296      REAL(wp) :: tabin(k1:k2,n1:n2)
297      !!----------------------------------------------------------------------
298      !
299      IF (before) THEN
300         AGRIF_SpecialValue = -999._wp
301         zrho_xy = Agrif_rhox() * Agrif_rhoy() 
302         DO jn = n1, n2-1
303            DO jk = k1, k2
304               DO jj = j1, j2
305                  DO ji = i1, i2
306                     tabres(ji,jj,jk,jn) = (tsn(ji,jj,jk,jn) * e3t_n(ji,jj,jk) ) &
307                                           * tmask(ji,jj,jk) + (tmask(ji,jj,jk)-1)*999._wp
308                  END DO
309               END DO
310            END DO
311         END DO
312         DO jk = k1, k2
313            DO jj = j1, j2
314               DO ji = i1, i2
315                  tabres(ji,jj,jk,n2) = tmask(ji,jj,jk) * e3t_n(ji,jj,jk) &
316                                           + (tmask(ji,jj,jk)-1)*999._wp
317               END DO
318            END DO
319         END DO
320      ELSE
321         tabres_child(:,:,:,:) = 0.
322         AGRIF_SpecialValue = 0._wp
323         DO jj = j1 , j2
324            DO ji = i1, i2
325               N_in = 0
326               DO jk = k1, k2 !k2 = jpk of child grid
327                  IF ( tabres(ji,jj,jk,n2) == 0  )   EXIT
328                  N_in = N_in + 1
329                  tabin(jk,:) = tabres(ji,jj,jk,n1:n2-1)/tabres(ji,jj,jk,n2)
330                  h_in (N_in) = tabres(ji,jj,jk,n2)
331               END DO
332               N_out = 0
333               DO jk = 1, jpk ! jpk of parent grid
334                  IF (tmask(ji,jj,jk) < -900)   EXIT ! TODO: Will not work with ISF
335                  N_out = N_out + 1
336                  h_out(N_out) = e3t_n(ji,jj,jk) 
337               END DO
338               IF (N_in > 0) THEN !Remove this?
339                  h_diff = sum(h_out(1:N_out))-sum(h_in(1:N_in))
340                  IF (h_diff < -1.e-4) THEN
341                     print *,'CHECK YOUR bathy T points ...',ji,jj,h_diff,sum(h_in(1:N_in)),sum(h_out(1:N_out))
342                     print *,h_in(1:N_in)
343                     print *,h_out(1:N_out)
344                     STOP
345                  ENDIF
346                  DO jn = n1, n2-1
347                     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)
348                  END DO
349               ENDIF
350            END DO
351         END DO
352
353         IF ( .NOT.( lk_agrif_fstep .AND. l_1st_euler ) ) THEN       ! Add asselin part
354
355!!gm         IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN
356            DO jn = n1, n2-1
357               DO jk = 1, jpk
358                  DO jj = j1, j2
359                     DO ji = i1, i2
360                        IF( tabres_child(ji,jj,jk,jn) /= 0. ) THEN
361                           tsb(ji,jj,jk,jn) = tsb(ji,jj,jk,jn)   & 
362                              &             + rn_atfp * ( tabres_child(ji,jj,jk,jn) - tsn(ji,jj,jk,jn) ) * tmask(ji,jj,jk)
363                        ENDIF
364                     END DO
365                  END DO
366               END DO
367            END DO
368         ENDIF
369         DO jn = n1, n2-1
370            DO jk = 1, jpk
371               DO jj = j1, j2
372                  DO ji = i1, i2
373                     IF( tabres_child(ji,jj,jk,jn) /= 0. ) THEN
374                        tsn(ji,jj,jk,jn) = tabres_child(ji,jj,jk,jn) * tmask(ji,jj,jk)
375                     END IF
376                  END DO
377               END DO
378            END DO
379         END DO
380      ENDIF
381      !
382   END SUBROUTINE updateTS
383
384# else
385
386   SUBROUTINE updateTS( tabres, i1, i2, j1, j2, k1, k2, n1, n2, before )
387      !!----------------------------------------------------------------------
388      !!                   *** ROUTINE ROUTINE updateT ***
389      !!----------------------------------------------------------------------
390      INTEGER, INTENT(in) :: i1,i2,j1,j2,k1,k2,n1,n2
391      REAL(wp),DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) :: tabres
392      LOGICAL, INTENT(in) :: before
393      !
394      INTEGER :: ji,jj,jk,jn
395      REAL(wp) :: ztb, ztnu, ztno
396      !!----------------------------------------------------------------------
397      !
398      IF (before) THEN
399         DO jn = 1,jpts
400            DO jk=k1,k2
401               DO jj=j1,j2
402                  DO ji=i1,i2
403!> jc tmp
404                     tabres(ji,jj,jk,jn) = tsn(ji,jj,jk,jn)  * e3t_n(ji,jj,jk) / e3t_0(ji,jj,jk)
405!                     tabres(ji,jj,jk,jn) = tsn(ji,jj,jk,jn)  * e3t_n(ji,jj,jk)
406!< jc tmp
407                  END DO
408               END DO
409            END DO
410         END DO
411      ELSE
412!> jc tmp
413         DO jn = 1,jpts
414            tabres(i1:i2,j1:j2,k1:k2,jn) =  tabres(i1:i2,j1:j2,k1:k2,jn) * e3t_0(i1:i2,j1:j2,k1:k2) &
415                                         & * tmask(i1:i2,j1:j2,k1:k2)
416         END DO
417!< jc tmp
418         IF ( .NOT.( lk_agrif_fstep .AND. l_1st_euler ) ) THEN
419!!gm         IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN
420            ! Add asselin part
421            DO jn = 1,jpts
422               DO jk = k1, k2
423                  DO jj = j1, j2
424                     DO ji = i1, i2
425                        IF( tabres(ji,jj,jk,jn) /= 0._wp ) THEN
426                           ztb  = tsb(ji,jj,jk,jn) * e3t_b(ji,jj,jk) ! fse3t_b prior update should be used
427                           ztnu = tabres(ji,jj,jk,jn)
428                           ztno = tsn(ji,jj,jk,jn) * e3t_a(ji,jj,jk)
429                           tsb(ji,jj,jk,jn) = ( ztb + rn_atfp * ( ztnu - ztno) ) / e3t_b(ji,jj,jk) * tmask(ji,jj,jk)
430                        ENDIF
431                     END DO
432                  END DO
433               END DO
434            END DO
435         ENDIF
436         DO jn = 1,jpts
437            DO jk=k1,k2
438               DO jj=j1,j2
439                  DO ji=i1,i2
440                     IF( tabres(ji,jj,jk,jn) /= 0._wp ) THEN
441                        tsn(ji,jj,jk,jn) = tabres(ji,jj,jk,jn) / e3t_n(ji,jj,jk)
442                     END IF
443                  END DO
444               END DO
445            END DO
446         END DO
447         !
448         IF ( l_1st_euler .AND. Agrif_Nb_Step() == 0 ) THEN
449!!gm         IF  ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN
450            tsb(i1:i2,j1:j2,k1:k2,1:jpts)  = tsn(i1:i2,j1:j2,k1:k2,1:jpts)
451         ENDIF
452         !
453      ENDIF
454      !
455   END SUBROUTINE updateTS
456
457# endif
458
459# if defined key_vertical
460
461   SUBROUTINE updateu( tabres, i1, i2, j1, j2, k1, k2, n1, n2, before )
462      !!----------------------------------------------------------------------
463      !!                   *** ROUTINE updateu ***
464      !!----------------------------------------------------------------------
465      INTEGER                                     , INTENT(in   ) :: i1, i2, j1, j2, k1, k2, n1, n2
466      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) :: tabres
467      LOGICAL                                     , INTENT(in   ) :: before
468      !
469      INTEGER ::   ji, jj, jk
470      REAL(wp)::   zrhoy
471! VERTICAL REFINEMENT BEGIN
472      REAL(wp), DIMENSION(i1:i2,j1:j2,1:jpk) :: tabres_child
473      REAL(wp) :: h_in(k1:k2)
474      REAL(wp) :: h_out(1:jpk)
475      INTEGER  :: N_in, N_out
476      REAL(wp) :: h_diff, excess, thick
477      REAL(wp) :: tabin(k1:k2)
478! VERTICAL REFINEMENT END
479      !!----------------------------------------------------------------------
480      !
481      IF( before ) THEN
482         zrhoy = Agrif_Rhoy()
483         AGRIF_SpecialValue = -999._wp
484         DO jk=k1,k2
485            DO jj=j1,j2
486               DO ji=i1,i2
487                  tabres(ji,jj,jk,1) = zrhoy * e2u(ji,jj) * e3u_n(ji,jj,jk) * umask(ji,jj,jk) * un(ji,jj,jk)  &
488                                       + (umask(ji,jj,jk)-1)*999._wp
489                  tabres(ji,jj,jk,2) = zrhoy * umask(ji,jj,jk) * e2u(ji,jj) * e3u_n(ji,jj,jk)  &
490                                       + (umask(ji,jj,jk)-1)*999._wp
491               END DO
492            END DO
493         END DO
494      ELSE
495         tabres_child(:,:,:) = 0.
496         AGRIF_SpecialValue = 0._wp
497         DO jj=j1,j2
498            DO ji=i1,i2
499               N_in = 0
500               h_in(:) = 0._wp
501               tabin(:) = 0._wp
502               DO jk=k1,k2 !k2=jpk of child grid
503                  IF( tabres(ji,jj,jk,2) < -900) EXIT
504                  N_in = N_in + 1
505                  tabin(jk) = tabres(ji,jj,jk,1)/tabres(ji,jj,jk,2)
506                  h_in(N_in) = tabres(ji,jj,jk,2)/e2u(ji,jj)
507               END DO
508               N_out = 0
509               DO jk=1,jpk
510                  IF (umask(ji,jj,jk) == 0) EXIT
511                  N_out = N_out + 1
512                  h_out(N_out) = e3u_n(ji,jj,jk)
513               END DO
514               IF (N_in * N_out > 0) THEN
515                  h_diff = sum(h_out(1:N_out))-sum(h_in(1:N_in))
516                  IF (h_diff < -1.e-4) THEN
517!Even if bathy at T points match it's possible for the U points to be deeper in the child grid.
518!In this case we need to move transport from the child grid cells below bed of parent grid into the bottom cell.
519                     excess = 0._wp
520                     DO jk=N_in,1,-1
521                        thick = MIN(-1*h_diff, h_in(jk))
522                        excess = excess + tabin(jk)*thick*e2u(ji,jj)
523                        tabin(jk) = tabin(jk)*(1. - thick/h_in(jk))
524                        h_diff = h_diff + thick
525                        IF ( h_diff == 0) THEN
526                           N_in = jk
527                           h_in(jk) = h_in(jk) - thick
528                           EXIT
529                        ENDIF
530                     END DO
531                  ENDIF
532                  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)
533                  tabres_child(ji,jj,N_out) = tabres_child(ji,jj,N_out) + excess/(e2u(ji,jj)*h_out(N_out))
534               ENDIF
535            END DO
536         END DO
537
538         DO jk = 1, jpk
539            DO jj = j1, j2
540               DO ji = i1, i2
541                  IF ( .NOT.( lk_agrif_fstep .AND. l_1st_euler) ) THEN    ! Add asselin part
542!!gm                  IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN   ! Add asselin part
543                     ub(ji,jj,jk) = ub(ji,jj,jk) + rn_atfp * ( tabres_child(ji,jj,jk) - un(ji,jj,jk) ) * umask(ji,jj,jk)
544                  ENDIF
545                  un(ji,jj,jk) = tabres_child(ji,jj,jk) * umask(ji,jj,jk)
546               END DO
547            END DO
548         END DO
549      ENDIF
550      !
551   END SUBROUTINE updateu
552
553#else
554
555   SUBROUTINE updateu( tabres, i1, i2, j1, j2, k1, k2, n1, n2, before )
556      !!----------------------------------------------------------------------
557      !!                   *** ROUTINE updateu ***
558      !!----------------------------------------------------------------------
559      INTEGER                               , INTENT(in   ) :: i1, i2, j1, j2, k1, k2, n1, n2
560      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) :: tabres
561      LOGICAL                                     , INTENT(in   ) :: before
562      !
563      INTEGER  :: ji, jj, jk
564      REAL(wp) :: zrhoy, zub, zunu, zuno
565      !!----------------------------------------------------------------------
566      !
567      IF( before ) THEN
568         zrhoy = Agrif_Rhoy()
569         DO jk = k1, k2
570            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)
571         END DO
572      ELSE
573         DO jk=k1,k2
574            DO jj=j1,j2
575               DO ji=i1,i2
576                  tabres(ji,jj,jk,1) = tabres(ji,jj,jk,1) * r1_e2u(ji,jj) 
577                  !
578                  IF ( .NOT.( lk_agrif_fstep .AND. l_1st_euler ) ) THEN    ! Add asselin part
579!!gm                  IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN   ! Add asselin part
580                     zub  = ub(ji,jj,jk) * e3u_b(ji,jj,jk)  ! fse3t_b prior update should be used
581                     zuno = un(ji,jj,jk) * e3u_a(ji,jj,jk)
582                     zunu = tabres(ji,jj,jk,1)
583                     ub(ji,jj,jk) = ( zub + rn_atfp * ( zunu - zuno) ) / e3u_b(ji,jj,jk) * umask(ji,jj,jk)
584                  ENDIF
585                  !
586                  un(ji,jj,jk) = tabres(ji,jj,jk,1) / e3u_n(ji,jj,jk) * umask(ji,jj,jk)
587               END DO
588            END DO
589         END DO
590         !
591         IF ( l_1st_euler .AND. Agrif_Nb_Step() == 0 ) THEN
592!!gm         IF  ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN
593            ub(i1:i2,j1:j2,k1:k2)  = un(i1:i2,j1:j2,k1:k2)
594         ENDIF
595         !
596      ENDIF
597      !
598   END SUBROUTINE updateu
599
600# endif
601
602   SUBROUTINE correct_u_bdy( tabres, i1, i2, j1, j2, k1, k2, n1, n2, before, nb, ndir )
603      !!----------------------------------------------------------------------
604      !!                   *** ROUTINE correct_u_bdy ***
605      !!----------------------------------------------------------------------
606      INTEGER                                     , INTENT(in   ) :: i1, i2, j1, j2, k1, k2, n1, n2
607      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) :: tabres
608      LOGICAL                                     , INTENT(in   ) :: before
609      INTEGER                                     , INTENT(in   ) :: nb, ndir
610      !!
611      LOGICAL :: western_side, eastern_side 
612      INTEGER ::   jj, jk
613      REAL(wp)::   zcor
614      !!----------------------------------------------------------------------
615      !
616      IF( .NOT.before ) THEN
617         !
618         western_side  = (nb == 1).AND.(ndir == 1)
619         eastern_side  = (nb == 1).AND.(ndir == 2)
620         !
621         IF (western_side) THEN
622            DO jj=j1,j2
623               zcor = un_b(i1-1,jj) * hu_a(i1-1,jj) * r1_hu_n(i1-1,jj) - un_b(i1-1,jj)
624               un_b(i1-1,jj) = un_b(i1-1,jj) + zcor
625               DO jk=1,jpkm1
626                  un(i1-1,jj,jk) = un(i1-1,jj,jk) + zcor * umask(i1-1,jj,jk)
627               END DO
628            END DO
629         ENDIF
630         !
631         IF (eastern_side) THEN
632            DO jj=j1,j2
633               zcor = un_b(i2+1,jj) * hu_a(i2+1,jj) * r1_hu_n(i2+1,jj) - un_b(i2+1,jj)
634               un_b(i2+1,jj) = un_b(i2+1,jj) + zcor
635               DO jk=1,jpkm1
636                  un(i2+1,jj,jk) = un(i2+1,jj,jk) + zcor * umask(i2+1,jj,jk)
637               END DO
638            END DO
639         ENDIF
640         !
641      ENDIF
642      !
643   END SUBROUTINE correct_u_bdy
644
645# if  defined key_vertical
646
647   SUBROUTINE updatev( tabres, i1, i2, j1, j2, k1, k2, n1, n2, before )
648      !!----------------------------------------------------------------------
649      !!                   *** ROUTINE updatev ***
650      !!----------------------------------------------------------------------
651      INTEGER                                     , INTENT(in   ) :: i1, i2, j1, j2, k1, k2, n1, n2
652      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) :: tabres
653      LOGICAL                                     , INTENT(in   ) :: before
654      !
655      INTEGER  ::   ji, jj, jk
656      REAL(wp) ::   zrhox
657! VERTICAL REFINEMENT BEGIN
658      REAL(wp), DIMENSION(i1:i2,j1:j2,1:jpk) :: tabres_child
659      REAL(wp) :: h_in(k1:k2)
660      REAL(wp) :: h_out(1:jpk)
661      INTEGER :: N_in, N_out
662      REAL(wp) :: h_diff, excess, thick
663      REAL(wp) :: tabin(k1:k2)
664! VERTICAL REFINEMENT END
665      !!----------------------------------------------------------------------     
666      !
667      IF( before ) THEN
668         zrhox = Agrif_Rhox()
669         AGRIF_SpecialValue = -999._wp
670         DO jk=k1,k2
671            DO jj=j1,j2
672               DO ji=i1,i2
673                  tabres(ji,jj,jk,1) = zrhox * e1v(ji,jj) * e3v_n(ji,jj,jk) * vmask(ji,jj,jk) * vn(ji,jj,jk) &
674                                       + (vmask(ji,jj,jk)-1)*999._wp
675                  tabres(ji,jj,jk,2) = vmask(ji,jj,jk) * zrhox * e1v(ji,jj) * e3v_n(ji,jj,jk) &
676                                       + (vmask(ji,jj,jk)-1)*999._wp
677               END DO
678            END DO
679         END DO
680      ELSE
681         tabres_child(:,:,:) = 0.
682         AGRIF_SpecialValue = 0._wp
683         DO jj=j1,j2
684            DO ji=i1,i2
685               N_in = 0
686               DO jk=k1,k2
687                  IF (tabres(ji,jj,jk,2) < -900) EXIT
688                  N_in = N_in + 1
689                  tabin(jk) = tabres(ji,jj,jk,1)/tabres(ji,jj,jk,2)
690                  h_in(N_in) = tabres(ji,jj,jk,2)/e1v(ji,jj)
691               END DO
692               N_out = 0
693               DO jk=1,jpk
694                  IF (vmask(ji,jj,jk) == 0) EXIT
695                  N_out = N_out + 1
696                  h_out(N_out) = e3v_n(ji,jj,jk)
697               END DO
698               IF (N_in * N_out > 0) THEN
699                  h_diff = sum(h_out(1:N_out))-sum(h_in(1:N_in))
700                  IF (h_diff < -1.e-4) then
701!Even if bathy at T points match it's possible for the U points to be deeper in the child grid.
702!In this case we need to move transport from the child grid cells below bed of parent grid into the bottom cell.
703                     excess = 0._wp
704                     DO jk=N_in,1,-1
705                        thick = MIN(-1*h_diff, h_in(jk))
706                        excess = excess + tabin(jk)*thick*e2u(ji,jj)
707                        tabin(jk) = tabin(jk)*(1. - thick/h_in(jk))
708                        h_diff = h_diff + thick
709                        IF ( h_diff == 0) THEN
710                           N_in = jk
711                           h_in(jk) = h_in(jk) - thick
712                           EXIT
713                        ENDIF
714                     END DO
715                  ENDIF
716                  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)
717                  tabres_child(ji,jj,N_out) = tabres_child(ji,jj,N_out) + excess/(e1v(ji,jj)*h_out(N_out))
718               ENDIF
719            END DO
720         END DO
721
722         DO jk=1,jpk
723            DO jj=j1,j2
724               DO ji=i1,i2
725                  !
726                  IF( .NOT.( lk_agrif_fstep .AND. l_1st_euler ) ) THEN  ! Add asselin part
727!!gm                  IF( .NOT.(lk_agrif_fstep.AND.(neuler==0)) ) THEN  ! Add asselin part
728                     vb(ji,jj,jk) = vb(ji,jj,jk) + rn_atfp * ( tabres_child(ji,jj,jk) - vn(ji,jj,jk) ) * vmask(ji,jj,jk)
729                  ENDIF
730                  !
731                  vn(ji,jj,jk) = tabres_child(ji,jj,jk) * vmask(ji,jj,jk)
732               END DO
733            END DO
734         END DO
735      ENDIF
736      !
737   END SUBROUTINE updatev
738
739# else
740
741   SUBROUTINE updatev( tabres, i1, i2, j1, j2, k1, k2, n1, n2, before)
742      !!----------------------------------------------------------------------
743      !!                   *** ROUTINE updatev ***
744      !!----------------------------------------------------------------------
745      INTEGER                                     , INTENT(in   ) :: i1, i2, j1, j2, k1, k2, n1, n2
746      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) :: tabres
747      LOGICAL                                     , INTENT(in   ) :: before
748      !
749      INTEGER  :: ji, jj, jk
750      REAL(wp) :: zrhox, zvb, zvnu, zvno
751      !!----------------------------------------------------------------------     
752      !
753      IF (before) THEN
754         zrhox = Agrif_Rhox()
755         DO jk=k1,k2
756            DO jj=j1,j2
757               DO ji=i1,i2
758                  tabres(ji,jj,jk,1) = zrhox * e1v(ji,jj) * e3v_n(ji,jj,jk) * vn(ji,jj,jk)
759               END DO
760            END DO
761         END DO
762      ELSE
763         DO jk=k1,k2
764            DO jj=j1,j2
765               DO ji=i1,i2
766                  tabres(ji,jj,jk,1) = tabres(ji,jj,jk,1) * r1_e1v(ji,jj)
767                  !
768                  IF ( .NOT.( lk_agrif_fstep .AND. l_1st_euler ) ) THEN    ! Add asselin part
769!!gm                  IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN   ! Add asselin part
770                     zvb  = vb(ji,jj,jk) * e3v_b(ji,jj,jk) ! fse3t_b prior update should be used
771                     zvno = vn(ji,jj,jk) * e3v_a(ji,jj,jk)
772                     zvnu = tabres(ji,jj,jk,1)
773                     vb(ji,jj,jk) = ( zvb + rn_atfp * ( zvnu - zvno) ) / e3v_b(ji,jj,jk) * vmask(ji,jj,jk)
774                  ENDIF
775                  !
776                  vn(ji,jj,jk) = tabres(ji,jj,jk,1) / e3v_n(ji,jj,jk) * vmask(ji,jj,jk)
777               END DO
778            END DO
779         END DO
780         !
781         IF ( l_1st_euler .AND. Agrif_Nb_Step() == 0 ) THEN
782!!gm         IF  ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN
783            vb(i1:i2,j1:j2,k1:k2)  = vn(i1:i2,j1:j2,k1:k2)
784         ENDIF
785         !
786      ENDIF
787      !
788   END SUBROUTINE updatev
789
790# endif
791
792   SUBROUTINE correct_v_bdy( tabres, i1, i2, j1, j2, k1, k2, n1, n2, before, nb, ndir )
793      !!----------------------------------------------------------------------
794      !!                   *** ROUTINE correct_v_bdy ***
795      !!----------------------------------------------------------------------
796      INTEGER                                     , INTENT(in   ) :: i1, i2, j1, j2, k1, k2, n1, n2
797      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) :: tabres
798      LOGICAL                                     , INTENT(in   ) :: before
799      INTEGER                                     , INTENT(in)    :: nb, ndir
800      !!
801      LOGICAL :: southern_side, northern_side 
802      !
803      INTEGER  :: ji, jk
804      REAL(wp) :: zcor
805      !!----------------------------------------------------------------------
806      !
807      IF( .NOT.before ) THEN
808         !
809         southern_side = (nb == 2).AND.(ndir == 1)
810         northern_side = (nb == 2).AND.(ndir == 2)
811         !
812         IF (southern_side) THEN
813            DO ji=i1,i2
814               zcor = vn_b(ji,j1-1) * hv_a(ji,j1-1) * r1_hv_n(ji,j1-1) - vn_b(ji,j1-1)
815               vn_b(ji,j1-1) = vn_b(ji,j1-1) + zcor
816               DO jk=1,jpkm1
817                  vn(ji,j1-1,jk) = vn(ji,j1-1,jk) + zcor * vmask(ji,j1-1,jk)
818               END DO
819            END DO
820         ENDIF
821         !
822         IF (northern_side) THEN
823            DO ji=i1,i2
824               zcor = vn_b(ji,j2+1) * hv_a(ji,j2+1) * r1_hv_n(ji,j2+1) - vn_b(ji,j2+1)
825               vn_b(ji,j2+1) = vn_b(ji,j2+1) + zcor
826               DO jk=1,jpkm1
827                  vn(ji,j2+1,jk) = vn(ji,j2+1,jk) + zcor * vmask(ji,j2+1,jk)
828               END DO
829            END DO
830         ENDIF
831         !
832      ENDIF
833      !
834   END SUBROUTINE correct_v_bdy
835
836
837   SUBROUTINE updateu2d( tabres, i1, i2, j1, j2, before )
838      !!----------------------------------------------------------------------
839      !!                   *** ROUTINE updateu2d ***
840      !!----------------------------------------------------------------------
841      INTEGER                         , INTENT(in   ) ::   i1, i2, j1, j2
842      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) ::   tabres
843      LOGICAL                         , INTENT(in   ) ::   before
844      !
845      INTEGER  :: ji, jj, jk
846      REAL(wp) :: zrhoy
847      REAL(wp) :: zcorr
848      !!----------------------------------------------------------------------
849      !
850      IF( before ) THEN
851         zrhoy = Agrif_Rhoy()
852         DO jj=j1,j2
853            DO ji=i1,i2
854               tabres(ji,jj) = zrhoy * un_b(ji,jj) * hu_n(ji,jj) * e2u(ji,jj)
855            END DO
856         END DO
857      ELSE
858         DO jj=j1,j2
859            DO ji=i1,i2
860               tabres(ji,jj) =  tabres(ji,jj) * r1_e2u(ji,jj) 
861               !   
862               ! Update "now" 3d velocities:
863               spgu(ji,jj) = 0._wp
864               DO jk=1,jpkm1
865                  spgu(ji,jj) = spgu(ji,jj) + e3u_n(ji,jj,jk) * un(ji,jj,jk)
866               END DO
867               !
868               zcorr = (tabres(ji,jj) - spgu(ji,jj)) * r1_hu_n(ji,jj)
869               DO jk=1,jpkm1             
870                  un(ji,jj,jk) = un(ji,jj,jk) + zcorr * umask(ji,jj,jk)           
871               END DO
872               !
873               ! Update barotropic velocities:
874               IF ( .NOT.ln_dynspg_ts .OR. (ln_dynspg_ts.AND.(.NOT.ln_bt_fw)) ) THEN
875                  IF ( .NOT.( lk_agrif_fstep .AND. l_1st_euler ) ) THEN    ! Add asselin part
876!!gm                  IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN ! Add asselin part
877                     zcorr = (tabres(ji,jj) - un_b(ji,jj) * hu_a(ji,jj)) * r1_hu_b(ji,jj)
878                     ub_b(ji,jj) = ub_b(ji,jj) + rn_atfp * zcorr * umask(ji,jj,1)
879                  END IF
880               ENDIF   
881               un_b(ji,jj) = tabres(ji,jj) * r1_hu_n(ji,jj) * umask(ji,jj,1)
882               !       
883               ! Correct "before" velocities to hold correct bt component:
884               spgu(ji,jj) = 0.e0
885               DO jk=1,jpkm1
886                  spgu(ji,jj) = spgu(ji,jj) + e3u_b(ji,jj,jk) * ub(ji,jj,jk)
887               END DO
888               !
889               zcorr = ub_b(ji,jj) - spgu(ji,jj) * r1_hu_b(ji,jj)
890               DO jk=1,jpkm1             
891                  ub(ji,jj,jk) = ub(ji,jj,jk) + zcorr * umask(ji,jj,jk)           
892               END DO
893               !
894            END DO
895         END DO
896         !
897         IF ( l_1st_euler .AND. Agrif_Nb_Step() == 0 ) THEN
898!!gm         IF  ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN
899            ub_b(i1:i2,j1:j2)  = un_b(i1:i2,j1:j2)
900         ENDIF
901      ENDIF
902      !
903   END SUBROUTINE updateu2d
904
905
906   SUBROUTINE updatev2d( tabres, i1, i2, j1, j2, before )
907      !!----------------------------------------------------------------------
908      !!                   *** ROUTINE updatev2d ***
909      !!----------------------------------------------------------------------
910      INTEGER                         , INTENT(in   ) ::   i1, i2, j1, j2
911      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) ::   tabres
912      LOGICAL                         , INTENT(in   ) ::   before
913      !
914      INTEGER  :: ji, jj, jk
915      REAL(wp) :: zrhox, zcorr
916      !!----------------------------------------------------------------------
917      !
918      IF( before ) THEN
919         zrhox = Agrif_Rhox()
920         DO jj=j1,j2
921            DO ji=i1,i2
922               tabres(ji,jj) = zrhox * vn_b(ji,jj) * hv_n(ji,jj) * e1v(ji,jj) 
923            END DO
924         END DO
925      ELSE
926         DO jj=j1,j2
927            DO ji=i1,i2
928               tabres(ji,jj) =  tabres(ji,jj) * r1_e1v(ji,jj) 
929               !   
930               ! Update "now" 3d velocities:
931               spgv(ji,jj) = 0.e0
932               DO jk=1,jpkm1
933                  spgv(ji,jj) = spgv(ji,jj) + e3v_n(ji,jj,jk) * vn(ji,jj,jk)
934               END DO
935               !
936               zcorr = (tabres(ji,jj) - spgv(ji,jj)) * r1_hv_n(ji,jj)
937               DO jk=1,jpkm1             
938                  vn(ji,jj,jk) = vn(ji,jj,jk) + zcorr * vmask(ji,jj,jk)           
939               END DO
940               !
941               ! Update barotropic velocities:
942               IF ( .NOT.ln_dynspg_ts .OR. ( ln_dynspg_ts .AND. .NOT.ln_bt_fw ) ) THEN
943                  IF ( .NOT.( lk_agrif_fstep .AND. l_1st_euler ) ) THEN    ! Add asselin part
944!!gm                  IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN ! Add asselin part
945                     zcorr = (tabres(ji,jj) - vn_b(ji,jj) * hv_a(ji,jj)) * r1_hv_b(ji,jj)
946                     vb_b(ji,jj) = vb_b(ji,jj) + rn_atfp * zcorr * vmask(ji,jj,1)
947                  END IF
948               ENDIF             
949               vn_b(ji,jj) = tabres(ji,jj) * r1_hv_n(ji,jj) * vmask(ji,jj,1)
950               !       
951               ! Correct "before" velocities to hold correct bt component:
952               spgv(ji,jj) = 0.e0
953               DO jk=1,jpkm1
954                  spgv(ji,jj) = spgv(ji,jj) + e3v_b(ji,jj,jk) * vb(ji,jj,jk)
955               END DO
956               !
957               zcorr = vb_b(ji,jj) - spgv(ji,jj) * r1_hv_b(ji,jj)
958               DO jk=1,jpkm1             
959                  vb(ji,jj,jk) = vb(ji,jj,jk) + zcorr * vmask(ji,jj,jk)           
960               END DO
961               !
962            END DO
963         END DO
964         !
965         IF ( l_1st_euler .AND. Agrif_Nb_Step() == 0 ) THEN
966!!gm         IF  ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN
967            vb_b(i1:i2,j1:j2)  = vn_b(i1:i2,j1:j2)
968         ENDIF
969         !
970      ENDIF
971      !
972   END SUBROUTINE updatev2d
973
974
975   SUBROUTINE updateSSH( tabres, i1, i2, j1, j2, before )
976      !!----------------------------------------------------------------------
977      !!                   *** ROUTINE updateSSH ***
978      !!----------------------------------------------------------------------
979      INTEGER                         , INTENT(in   ) ::   i1, i2, j1, j2
980      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) ::   tabres
981      LOGICAL                         , INTENT(in   ) ::   before
982      !
983      INTEGER :: ji, jj
984      !!----------------------------------------------------------------------
985      !
986      IF( before ) THEN
987         DO jj = j1, j2
988            DO ji = i1, i2
989               tabres(ji,jj) = ssh(ji,jj,Nnn)
990            END DO
991         END DO
992      ELSE
993         IF ( .NOT.( lk_agrif_fstep .AND. l_1st_euler ) ) THEN
994!!gm         IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN
995            DO jj = j1, j2
996               DO ji = i1, i2
997                  ssh(ji,jj,Nbb) = ssh(ji,jj,Nbb) + rn_atfp * ( tabres(ji,jj) - ssh(ji,jj,Nnn) ) * tmask(ji,jj,1)
998               END DO
999            END DO
1000         ENDIF
1001         !
1002         DO jj = j1, j2
1003            DO ji = i1, i2
1004               ssh(ji,jj,Nnn) = tabres(ji,jj) * tmask(ji,jj,1)
1005            END DO
1006         END DO
1007         !
1008         IF ( l_1st_euler .AND. Agrif_Nb_Step() == 0 ) THEN
1009!!gm         IF  ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN
1010            ssh(i1:i2,j1:j2,Nbb)  = ssh(i1:i2,j1:j2,Nnn)
1011         ENDIF
1012         !
1013      ENDIF
1014      !
1015   END SUBROUTINE updateSSH
1016
1017
1018   SUBROUTINE updateub2b( tabres, i1, i2, j1, j2, before )
1019      !!----------------------------------------------------------------------
1020      !!                      *** ROUTINE updateub2b ***
1021      !!----------------------------------------------------------------------
1022      INTEGER                            , INTENT(in) ::   i1, i2, j1, j2
1023      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) ::   tabres
1024      LOGICAL                            , INTENT(in) ::   before
1025      !!
1026      INTEGER :: ji, jj
1027      REAL(wp) :: zrhoy, za1, zcor
1028      !!---------------------------------------------
1029      !
1030      IF (before) THEN
1031         zrhoy = Agrif_Rhoy()
1032         DO jj=j1,j2
1033            DO ji=i1,i2
1034               tabres(ji,jj) = ub2_i_b(ji,jj) * e2u(ji,jj)
1035            END DO
1036         END DO
1037         tabres = zrhoy * tabres
1038      ELSE
1039         !
1040         tabres(i1:i2,j1:j2) = tabres(i1:i2,j1:j2) * r1_e2u(i1:i2,j1:j2)
1041         !
1042         za1 = 1._wp / REAL(Agrif_rhot(), wp)
1043         DO jj=j1,j2
1044            DO ji=i1,i2
1045               zcor=tabres(ji,jj) - ub2_b(ji,jj)
1046               ! Update time integrated fluxes also in case of multiply nested grids:
1047               ub2_i_b(ji,jj) = ub2_i_b(ji,jj) + za1 * zcor 
1048               ! Update corrective fluxes:
1049               un_bf(ji,jj)  = un_bf(ji,jj) + zcor
1050               ! Update half step back fluxes:
1051               ub2_b(ji,jj) = tabres(ji,jj)
1052            END DO
1053         END DO
1054      ENDIF
1055      !
1056   END SUBROUTINE updateub2b
1057
1058
1059   SUBROUTINE reflux_sshu( tabres, i1, i2, j1, j2, before, nb, ndir )
1060      !!----------------------------------------------------------------------
1061      !!                   *** ROUTINE reflux_sshu ***
1062      !!----------------------------------------------------------------------
1063      INTEGER                         , INTENT(in   ) ::   i1, i2, j1, j2
1064      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) ::   tabres
1065      LOGICAL                         , INTENT(in   ) ::   before
1066      INTEGER                         , INTENT(in   ) ::   nb, ndir
1067      !
1068      LOGICAL ::   western_side, eastern_side 
1069      INTEGER ::   ji, jj
1070      REAL(wp)::   zrhoy, za1, zcor
1071      !!----------------------------------------------------------------------
1072      !
1073      IF (before) THEN
1074         zrhoy = Agrif_Rhoy()
1075         DO jj=j1,j2
1076            DO ji=i1,i2
1077               tabres(ji,jj) = ub2_i_b(ji,jj) * e2u(ji,jj)
1078            END DO
1079         END DO
1080         tabres = zrhoy * tabres
1081      ELSE
1082         !
1083         tabres(i1:i2,j1:j2) = tabres(i1:i2,j1:j2) * r1_e2u(i1:i2,j1:j2)
1084         !
1085         western_side  = (nb == 1).AND.(ndir == 1)
1086         eastern_side  = (nb == 1).AND.(ndir == 2)
1087         !
1088         IF ( western_side ) THEN
1089            DO jj=j1,j2
1090               zcor = rn_Dt * r1_e1e2t(i1  ,jj) * e2u(i1,jj) * (ub2_b(i1,jj)-tabres(i1,jj)) 
1091               ssh(i1  ,jj,Nnn) = ssh(i1  ,jj,Nnn) + zcor
1092               IF ( .NOT.( lk_agrif_fstep .AND. l_1st_euler ) )   ssh(i1  ,jj,Nbb) = ssh(i1  ,jj,Nbb) + rn_atfp * zcor
1093!!gm               IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) sshb(i1  ,jj) = sshb(i1  ,jj) + rn_atfp * zcor
1094            END DO
1095         ENDIF
1096         IF ( eastern_side ) THEN
1097            DO jj=j1,j2
1098               zcor = - rn_Dt * r1_e1e2t(i2+1,jj) * e2u(i2,jj) * (ub2_b(i2,jj)-tabres(i2,jj))
1099               ssh(i2+1,jj,Nnn) = ssh(i2+1,jj,Nnn) + zcor
1100               IF (.NOT.( lk_agrif_fstep .AND. l_1st_euler ) )   ssh(i2+1,jj,Nbb) = ssh(i2+1,jj,Nbb) + rn_atfp * zcor
1101!!gm               IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) sshb(i2+1,jj) = sshb(i2+1,jj) + rn_atfp * zcor
1102            END DO
1103         ENDIF
1104         !
1105      ENDIF
1106      !
1107   END SUBROUTINE reflux_sshu
1108
1109
1110   SUBROUTINE updatevb2b( tabres, i1, i2, j1, j2, before )
1111      !!----------------------------------------------------------------------
1112      !!                    *** ROUTINE updatevb2b ***
1113      !!----------------------------------------------------------------------
1114      INTEGER                         , INTENT(in   ) ::   i1, i2, j1, j2
1115      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) ::   tabres
1116      LOGICAL                         , INTENT(in   ) ::   before
1117      !
1118      INTEGER :: ji, jj
1119      REAL(wp) :: zrhox, za1, zcor
1120      !!---------------------------------------------------------------------
1121      !
1122      IF( before ) THEN
1123         zrhox = Agrif_Rhox()
1124         DO jj=j1,j2
1125            DO ji=i1,i2
1126               tabres(ji,jj) = vb2_i_b(ji,jj) * e1v(ji,jj) 
1127            END DO
1128         END DO
1129         tabres = zrhox * tabres
1130      ELSE
1131         !
1132         tabres(i1:i2,j1:j2) = tabres(i1:i2,j1:j2) * r1_e1v(i1:i2,j1:j2)
1133         !
1134         za1 = 1._wp / REAL(Agrif_rhot(), wp)
1135         DO jj=j1,j2
1136            DO ji=i1,i2
1137               zcor=tabres(ji,jj) - vb2_b(ji,jj)
1138               ! Update time integrated fluxes also in case of multiply nested grids:
1139               vb2_i_b(ji,jj) = vb2_i_b(ji,jj) + za1 * zcor 
1140               ! Update corrective fluxes:
1141               vn_bf(ji,jj)  = vn_bf(ji,jj) + zcor
1142               ! Update half step back fluxes:
1143               vb2_b(ji,jj) = tabres(ji,jj)
1144            END DO
1145         END DO
1146      ENDIF
1147      !
1148   END SUBROUTINE updatevb2b
1149
1150
1151   SUBROUTINE reflux_sshv( tabres, i1, i2, j1, j2, before, nb, ndir )
1152      !!----------------------------------------------------------------------
1153      !!                   *** ROUTINE reflux_sshv ***
1154      !!----------------------------------------------------------------------
1155      INTEGER                         , INTENT(in   ) ::   i1, i2, j1, j2
1156      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) ::   tabres
1157      LOGICAL                         , INTENT(in   ) ::   before
1158      INTEGER                         , INTENT(in   ) ::   nb, ndir
1159      !!
1160      LOGICAL :: southern_side, northern_side 
1161      INTEGER :: ji, jj
1162      REAL(wp) :: zrhox, za1, zcor
1163      !!----------------------------------------------------------------------
1164      !
1165      IF (before) THEN
1166         zrhox = Agrif_Rhox()
1167         DO jj=j1,j2
1168            DO ji=i1,i2
1169               tabres(ji,jj) = vb2_i_b(ji,jj) * e1v(ji,jj) 
1170            END DO
1171         END DO
1172         tabres = zrhox * tabres
1173      ELSE
1174         !
1175         tabres(i1:i2,j1:j2) = tabres(i1:i2,j1:j2) * r1_e1v(i1:i2,j1:j2)
1176         !
1177         southern_side = (nb == 2).AND.(ndir == 1)
1178         northern_side = (nb == 2).AND.(ndir == 2)
1179         !
1180         IF (southern_side) THEN
1181            DO ji=i1,i2
1182               zcor = rn_Dt * r1_e1e2t(ji,j1  ) * e1v(ji,j1  ) * ( vb2_b(ji,j1)-tabres(ji,j1) )
1183               ssh(ji,j1  ,Nnn) = ssh(ji,j1  ,Nnn) + zcor
1184               IF ( .NOT.( lk_agrif_fstep .AND. l_1st_euler ) )   ssh(ji,j1  ,Nbb) = ssh(ji,j1,Nbb) + rn_atfp * zcor
1185!!gm               IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) sshb(ji,j1  ) = sshb(ji,j1) + rn_atfp * zcor
1186            END DO
1187         ENDIF
1188         IF (northern_side) THEN               
1189            DO ji=i1,i2
1190               zcor = - rn_Dt * r1_e1e2t(ji,j2+1) * e1v(ji,j2  ) * ( vb2_b(ji,j2)-tabres(ji,j2) )
1191               ssh(ji,j2+1,Nnn) = ssh(ji,j2+1,Nnn) + zcor
1192               IF ( .NOT.( lk_agrif_fstep .AND. l_1st_euler ) )   ssh(ji,j2+1,Nbb) = ssh(ji,j2+1,Nbb) + rn_atfp * zcor
1193!!gm               IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) sshb(ji,j2+1) = sshb(ji,j2+1) + rn_atfp * zcor
1194            END DO
1195         ENDIF
1196         !
1197      ENDIF
1198      !
1199   END SUBROUTINE reflux_sshv
1200
1201
1202   SUBROUTINE update_scales( tabres, i1, i2, j1, j2, k1, k2, n1,n2, before )
1203      !!----------------------------------------------------------------------
1204      !!                      *** ROUTINE updateT ***
1205      !
1206      ! ====>>>>>>>>>>    currently not used
1207      !
1208      !!----------------------------------------------------------------------
1209      INTEGER                                    , INTENT(in   ) ::   i1, i2, j1, j2, k1, k2, n1, n2
1210      REAL(wp),DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) ::   tabres
1211      LOGICAL                                    , INTENT(in   ) ::   before
1212      !!
1213      INTEGER :: ji,jj,jk
1214      REAL(wp) :: ztemp
1215      !!----------------------------------------------------------------------
1216
1217      IF (before) THEN
1218         DO jk=k1,k2
1219            DO jj=j1,j2
1220               DO ji=i1,i2
1221                  tabres(ji,jj,jk,1) = e1t(ji,jj)*e2t(ji,jj)*tmask(ji,jj,jk)
1222                  tabres(ji,jj,jk,2) = e1t(ji,jj)*tmask(ji,jj,jk)
1223                  tabres(ji,jj,jk,3) = e2t(ji,jj)*tmask(ji,jj,jk)
1224               END DO
1225            END DO
1226         END DO
1227         tabres(:,:,:,1)=tabres(:,:,:,1)*Agrif_Rhox()*Agrif_Rhoy()
1228         tabres(:,:,:,2)=tabres(:,:,:,2)*Agrif_Rhox()
1229         tabres(:,:,:,3)=tabres(:,:,:,3)*Agrif_Rhoy()
1230      ELSE
1231         DO jk=k1,k2
1232            DO jj=j1,j2
1233               DO ji=i1,i2
1234                  IF( tabres(ji,jj,jk,1) .NE. 0. ) THEN
1235                     print *,'VAL = ',ji,jj,jk,tabres(ji,jj,jk,1),e1t(ji,jj)*e2t(ji,jj)*tmask(ji,jj,jk)
1236                     print *,'VAL2 = ',ji,jj,jk,tabres(ji,jj,jk,2),e1t(ji,jj)*tmask(ji,jj,jk)
1237                     print *,'VAL3 = ',ji,jj,jk,tabres(ji,jj,jk,3),e2t(ji,jj)*tmask(ji,jj,jk)
1238                     ztemp = sqrt(tabres(ji,jj,jk,1)/(tabres(ji,jj,jk,2)*tabres(ji,jj,jk,3)))
1239                     print *,'CORR = ',ztemp-1.
1240                     print *,'NEW VALS = ',tabres(ji,jj,jk,2)*ztemp,tabres(ji,jj,jk,3)*ztemp, &
1241                           tabres(ji,jj,jk,2)*ztemp*tabres(ji,jj,jk,3)*ztemp
1242                     e1t(ji,jj) = tabres(ji,jj,jk,2)*ztemp
1243                     e2t(ji,jj) = tabres(ji,jj,jk,3)*ztemp
1244                  END IF
1245               END DO
1246            END DO
1247         END DO
1248      ENDIF
1249      !
1250   END SUBROUTINE update_scales
1251
1252
1253   SUBROUTINE updateEN( ptab, i1, i2, j1, j2, k1, k2, before )
1254      !!----------------------------------------------------------------------
1255      !!                      *** ROUTINE updateen ***
1256      !!----------------------------------------------------------------------
1257      INTEGER                               , INTENT(in   ) ::   i1, i2, j1, j2, k1, k2
1258      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) ::   ptab
1259      LOGICAL                               , INTENT(in   ) ::   before
1260      !!----------------------------------------------------------------------
1261      !
1262      IF( before ) THEN
1263         ptab (i1:i2,j1:j2,k1:k2) = en(i1:i2,j1:j2,k1:k2)
1264      ELSE
1265         en(i1:i2,j1:j2,k1:k2) = ptab (i1:i2,j1:j2,k1:k2) 
1266      ENDIF
1267      !
1268   END SUBROUTINE updateEN
1269
1270
1271   SUBROUTINE updateAVT( ptab, i1, i2, j1, j2, k1, k2, before )
1272      !!----------------------------------------------------------------------
1273      !!                      *** ROUTINE updateavt ***
1274      !!----------------------------------------------------------------------
1275      INTEGER                               , INTENT(in   ) ::   i1, i2, j1, j2, k1, k2
1276      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) ::   ptab
1277      LOGICAL                               , INTENT(in   ) ::   before
1278      !!----------------------------------------------------------------------
1279      !
1280      IF( before ) THEN   ;   ptab (i1:i2,j1:j2,k1:k2) = avt_k(i1:i2,j1:j2,k1:k2)
1281      ELSE                ;   avt_k(i1:i2,j1:j2,k1:k2) = ptab (i1:i2,j1:j2,k1:k2) 
1282      ENDIF
1283      !
1284   END SUBROUTINE updateAVT
1285
1286
1287   SUBROUTINE updateAVM( ptab, i1, i2, j1, j2, k1, k2, before )
1288      !!----------------------------------------------------------------------
1289      !!                      *** ROUTINE updateavm ***
1290      !!----------------------------------------------------------------------
1291      INTEGER                               , INTENT(in   ) ::   i1, i2, j1, j2, k1, k2
1292      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) ::   ptab
1293      LOGICAL                               , INTENT(in   ) ::   before
1294      !!----------------------------------------------------------------------
1295      !
1296      IF( before ) THEN   ;   ptab (i1:i2,j1:j2,k1:k2) = avm_k(i1:i2,j1:j2,k1:k2)
1297      ELSE                ;   avm_k(i1:i2,j1:j2,k1:k2) = ptab (i1:i2,j1:j2,k1:k2) 
1298      ENDIF
1299      !
1300   END SUBROUTINE updateAVM
1301
1302
1303   SUBROUTINE updatee3t(ptab_dum, i1, i2, j1, j2, k1, k2, before )
1304      !!----------------------------------------------------------------------
1305      !!                   *** ROUTINE updatee3t ***
1306      !!----------------------------------------------------------------------
1307      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2) :: ptab_dum
1308      INTEGER, INTENT(in) :: i1, i2, j1, j2, k1, k2
1309      LOGICAL, INTENT(in) :: before
1310      !
1311      REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: ptab
1312      INTEGER :: ji,jj,jk
1313      REAL(wp) :: zcoef
1314      !!---------------------------------------------
1315      !
1316      IF (.NOT.before) THEN
1317         !
1318         ALLOCATE( ptab(i1:i2,j1:j2,1:jpk) ) 
1319         !
1320         ! Update e3t from ssh (z* case only)
1321         DO jk = 1, jpkm1
1322            DO jj = j1, j2
1323               DO ji = i1, i2
1324                  ptab(ji,jj,jk) = e3t_0(ji,jj,jk) * (1._wp + ssh(ji,jj,Nnn) * r1_ht_0(ji,jj) *tmask(ji,jj,jk) )
1325               END DO
1326            END DO
1327         END DO
1328         !
1329         ! 1) Updates at BEFORE time step:
1330         ! -------------------------------
1331         !
1332         ! Save "old" scale factor (prior update) for subsequent asselin correction
1333         ! of prognostic variables
1334         e3t_a(i1:i2,j1:j2,1:jpkm1) = e3t_n(i1:i2,j1:j2,1:jpkm1)
1335
1336         ! One should also save e3t_b, but lacking of workspace...
1337!         hdivn(i1:i2,j1:j2,1:jpkm1)   = e3t_b(i1:i2,j1:j2,1:jpkm1)
1338
1339         IF ( .NOT.( lk_agrif_fstep .AND. l_1st_euler ) ) THEN
1340!!gm         IF (.NOT.(lk_agrif_fstep.AND.(neuler==0) )) THEN
1341            DO jk = 1, jpkm1
1342               DO jj=j1,j2
1343                  DO ji=i1,i2
1344                     e3t_b(ji,jj,jk) =  e3t_b(ji,jj,jk) + rn_atfp * ( ptab(ji,jj,jk) - e3t_n(ji,jj,jk) )
1345                  END DO
1346               END DO
1347            END DO
1348            !
1349            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)
1350            gdepw_b(i1:i2,j1:j2,1) = 0.0_wp
1351            gdept_b(i1:i2,j1:j2,1) = 0.5_wp * e3w_b(i1:i2,j1:j2,1)
1352            !
1353            DO jk = 2, jpk
1354               DO jj = j1,j2
1355                  DO ji = i1,i2           
1356                     zcoef = (tmask(ji,jj,jk) - wmask(ji,jj,jk))
1357                     e3w_b(ji,jj,jk)  = e3w_0(ji,jj,jk) + ( 1.0_wp - 0.5_wp * tmask(ji,jj,jk) ) *        & 
1358                     &                                        ( e3t_b(ji,jj,jk-1) - e3t_0(ji,jj,jk-1) )  &
1359                     &                                  +            0.5_wp * tmask(ji,jj,jk)   *        &
1360                     &                                        ( e3t_b(ji,jj,jk  ) - e3t_0(ji,jj,jk  ) )
1361                     gdepw_b(ji,jj,jk) = gdepw_b(ji,jj,jk-1) + e3t_b(ji,jj,jk-1)
1362                     gdept_b(ji,jj,jk) =      zcoef  * ( gdepw_b(ji,jj,jk  ) + 0.5 * e3w_b(ji,jj,jk))  &
1363                         &               + (1-zcoef) * ( gdept_b(ji,jj,jk-1) +       e3w_b(ji,jj,jk)) 
1364                  END DO
1365               END DO
1366            END DO
1367            !
1368         ENDIF       
1369         !
1370         ! 2) Updates at NOW time step:
1371         ! ----------------------------
1372         !
1373         ! Update vertical scale factor at T-points:
1374         e3t_n(i1:i2,j1:j2,1:jpkm1) = ptab(i1:i2,j1:j2,1:jpkm1)
1375         !
1376         ! Update total depth:
1377         ht_n(i1:i2,j1:j2) = 0._wp
1378         DO jk = 1, jpkm1
1379            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)
1380         END DO
1381         !
1382         ! Update vertical scale factor at W-points and depths:
1383         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)
1384         gdept_n(i1:i2,j1:j2,1) = 0.5_wp * e3w_n(i1:i2,j1:j2,1)
1385         gdepw_n(i1:i2,j1:j2,1) = 0.0_wp
1386         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
1387         !
1388         DO jk = 2, jpk
1389            DO jj = j1,j2
1390               DO ji = i1,i2           
1391               zcoef = (tmask(ji,jj,jk) - wmask(ji,jj,jk))
1392               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) )   &
1393               &                                  +            0.5_wp * tmask(ji,jj,jk)   * ( e3t_n(ji,jj,jk  ) - e3t_0(ji,jj,jk  ) )
1394               gdepw_n(ji,jj,jk) = gdepw_n(ji,jj,jk-1) + e3t_n(ji,jj,jk-1)
1395               gdept_n(ji,jj,jk) =      zcoef  * ( gdepw_n(ji,jj,jk  ) + 0.5 * e3w_n(ji,jj,jk))  &
1396                   &               + (1-zcoef) * ( gdept_n(ji,jj,jk-1) +       e3w_n(ji,jj,jk)) 
1397               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
1398               END DO
1399            END DO
1400         END DO
1401         !
1402         IF ( l_1st_euler .AND. Agrif_Nb_Step()==0 ) THEN
1403!!gm         IF  ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN
1404            e3t_b (i1:i2,j1:j2,1:jpk)  = e3t_n (i1:i2,j1:j2,1:jpk)
1405            e3w_b (i1:i2,j1:j2,1:jpk)  = e3w_n (i1:i2,j1:j2,1:jpk)
1406            gdepw_b(i1:i2,j1:j2,1:jpk) = gdepw_n(i1:i2,j1:j2,1:jpk)
1407            gdept_b(i1:i2,j1:j2,1:jpk) = gdept_n(i1:i2,j1:j2,1:jpk)
1408         ENDIF
1409         !
1410         DEALLOCATE(ptab)
1411      ENDIF
1412      !
1413   END SUBROUTINE updatee3t
1414
1415#else
1416   !!----------------------------------------------------------------------
1417   !!   Empty module                                          no AGRIF zoom
1418   !!----------------------------------------------------------------------
1419CONTAINS
1420   SUBROUTINE agrif_oce_update_empty
1421      WRITE(*,*)  'agrif_oce_update : You should not have seen this print! error?'
1422   END SUBROUTINE agrif_oce_update_empty
1423#endif
1424
1425   !!======================================================================
1426END MODULE agrif_oce_update
1427
Note: See TracBrowser for help on using the repository browser.