source: branches/2017/dev_r7963_nemo_v3_6_AGRIF-3_AGRIFVVL/NEMOGCM/NEMO/NST_SRC/agrif_opa_update.F90 @ 8010

Last change on this file since 8010 was 8010, checked in by jchanut, 3 years ago

AGRIF vvl add on

  • Property svn:keywords set to Id
File size: 35.9 KB
Line 
1#define TWO_WAY        /* TWO WAY NESTING */
2#undef DECAL_FEEDBACK /* SEPARATION of INTERFACES*/
3 
4MODULE agrif_opa_update
5#if defined key_agrif  && ! defined key_offline
6   USE par_oce
7   USE oce
8   USE dom_oce
9   USE agrif_oce
10   USE in_out_manager  ! I/O manager
11   USE lib_mpp
12   USE wrk_nemo 
13   USE dynspg_oce
14   USE zdf_oce        ! vertical physics: ocean variables
15   USE domvvl         ! Need interpolation routines
16
17   IMPLICIT NONE
18   PRIVATE
19
20   PUBLIC Agrif_Update_Tra, Agrif_Update_Dyn, Update_Scales, Agrif_Update_vvl
21
22# if defined key_zdftke
23   PUBLIC Agrif_Update_Tke
24# endif
25   !!----------------------------------------------------------------------
26   !! NEMO/NST 3.6 , NEMO Consortium (2010)
27   !! $Id$
28   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
29   !!----------------------------------------------------------------------
30
31CONTAINS
32
33   SUBROUTINE Agrif_Update_Tra( )
34      !!---------------------------------------------
35      !!   *** ROUTINE Agrif_Update_Tra ***
36      !!---------------------------------------------
37      !
38      IF (Agrif_Root()) RETURN
39      !
40#if defined TWO_WAY 
41      IF (lwp.AND.lk_agrif_debug) Write(*,*) 'Update tracers  from grid Number',Agrif_Fixed(), 'nbcline', nbcline
42
43      Agrif_UseSpecialValueInUpdate = .TRUE.
44      Agrif_SpecialValueFineGrid = 0.
45      !
46      IF (MOD(nbcline,nbclineupdate) == 0) THEN
47# if ! defined DECAL_FEEDBACK
48         CALL Agrif_Update_Variable(tsn_id, procname=updateTS)
49# else
50         CALL Agrif_Update_Variable(tsn_id, locupdate=(/1,0/),procname=updateTS)
51# endif
52      ELSE
53# if ! defined DECAL_FEEDBACK
54         CALL Agrif_Update_Variable(tsn_id,locupdate=(/0,2/), procname=updateTS)
55# else
56         CALL Agrif_Update_Variable(tsn_id,locupdate=(/1,2/), procname=updateTS)
57# endif
58      ENDIF
59      !
60      Agrif_UseSpecialValueInUpdate = .FALSE.
61      !
62#endif
63      !
64   END SUBROUTINE Agrif_Update_Tra
65
66   SUBROUTINE Agrif_Update_Dyn( )
67      !!---------------------------------------------
68      !!   *** ROUTINE Agrif_Update_Dyn ***
69      !!---------------------------------------------
70      !
71      IF (Agrif_Root()) RETURN
72      !
73#if defined TWO_WAY
74      IF (lwp.AND.lk_agrif_debug) Write(*,*) 'Update momentum from grid Number',Agrif_Fixed(), 'nbcline', nbcline
75
76      Agrif_UseSpecialValueInUpdate = .FALSE.
77      Agrif_SpecialValueFineGrid = 0.
78      !     
79      IF (mod(nbcline,nbclineupdate) == 0) THEN
80# if ! defined DECAL_FEEDBACK
81         CALL Agrif_Update_Variable(un_update_id,procname = updateU)
82         CALL Agrif_Update_Variable(vn_update_id,procname = updateV)
83# else
84         CALL Agrif_Update_Variable(un_update_id,locupdate1=(/0,-1/),locupdate2=(/1,-2/),procname = updateU)
85         CALL Agrif_Update_Variable(vn_update_id,locupdate1=(/1,-2/),locupdate2=(/0,-1/),procname = updateV)
86# endif
87      ELSE
88# if ! defined DECAL_FEEDBACK
89         CALL Agrif_Update_Variable(un_update_id,locupdate=(/0,1/),procname = updateU)
90         CALL Agrif_Update_Variable(vn_update_id,locupdate=(/0,1/),procname = updateV)         
91# else
92         CALL Agrif_Update_Variable(un_update_id,locupdate1=(/0,1/),locupdate2=(/1,1/),procname = updateU)
93         CALL Agrif_Update_Variable(vn_update_id,locupdate1=(/1,1/),locupdate2=(/0,1/),procname = updateV)
94# endif
95      ENDIF
96
97# if ! defined DECAL_FEEDBACK
98      CALL Agrif_Update_Variable(e1u_id,procname = updateU2d)
99      CALL Agrif_Update_Variable(e2v_id,procname = updateV2d) 
100# else
101      CALL Agrif_Update_Variable(e1u_id,locupdate1=(/0,-1/),locupdate2=(/1,-2/),procname = updateU2d)
102      CALL Agrif_Update_Variable(e2v_id,locupdate1=(/1,-2/),locupdate2=(/0,-1/),procname = updateV2d) 
103# endif
104
105# if defined key_dynspg_ts
106      IF (ln_bt_fw) THEN
107         ! Update time integrated transports
108         IF (mod(nbcline,nbclineupdate) == 0) THEN
109#  if ! defined DECAL_FEEDBACK
110            CALL Agrif_Update_Variable(ub2b_update_id,procname = updateub2b)
111            CALL Agrif_Update_Variable(vb2b_update_id,procname = updatevb2b)
112#  else
113            CALL Agrif_Update_Variable(ub2b_update_id,locupdate1=(/0,-1/),locupdate2=(/1,-2/),procname = updateub2b)
114            CALL Agrif_Update_Variable(vb2b_update_id,locupdate1=(/1,-2/),locupdate2=(/0,-1/),procname = updatevb2b)
115#  endif
116         ELSE
117#  if ! defined DECAL_FEEDBACK
118            CALL Agrif_Update_Variable(ub2b_update_id,locupdate=(/0,1/),procname = updateub2b)
119            CALL Agrif_Update_Variable(vb2b_update_id,locupdate=(/0,1/),procname = updatevb2b)
120#  else
121            CALL Agrif_Update_Variable(ub2b_update_id,locupdate1=(/0,1/),locupdate2=(/1,1/),procname = updateub2b)
122            CALL Agrif_Update_Variable(vb2b_update_id,locupdate1=(/1,1/),locupdate2=(/0,1/),procname = updatevb2b)
123#  endif
124         ENDIF
125      END IF
126# endif
127      !
128      nbcline = nbcline + 1
129      !
130      Agrif_UseSpecialValueInUpdate = .TRUE.
131      Agrif_SpecialValueFineGrid = 0.
132# if ! defined DECAL_FEEDBACK
133      CALL Agrif_Update_Variable(sshn_id,procname = updateSSH)
134# else
135      CALL Agrif_Update_Variable(sshn_id,locupdate=(/1,0/),procname = updateSSH)
136# endif
137      Agrif_UseSpecialValueInUpdate = .FALSE.
138      !
139#endif
140      !
141   END SUBROUTINE Agrif_Update_Dyn
142
143# if defined key_zdftke
144   SUBROUTINE Agrif_Update_Tke( kt )
145      !!---------------------------------------------
146      !!   *** ROUTINE Agrif_Update_Tke ***
147      !!---------------------------------------------
148      !!
149      INTEGER, INTENT(in) :: kt
150      !       
151      IF( ( Agrif_NbStepint() .NE. 0 ) .AND. (kt /= 0) ) RETURN
152#  if defined TWO_WAY
153
154      Agrif_UseSpecialValueInUpdate = .TRUE.
155      Agrif_SpecialValueFineGrid = 0.
156
157      CALL Agrif_Update_Variable( en_id, locupdate=(/0,0/), procname=updateEN  )
158      CALL Agrif_Update_Variable(avt_id, locupdate=(/0,0/), procname=updateAVT )
159      CALL Agrif_Update_Variable(avm_id, locupdate=(/0,0/), procname=updateAVM )
160
161      Agrif_UseSpecialValueInUpdate = .FALSE.
162
163#  endif
164     
165   END SUBROUTINE Agrif_Update_Tke
166# endif /* key_zdftke */
167
168   SUBROUTINE Agrif_Update_vvl( )
169      !!---------------------------------------------
170      !!   *** ROUTINE Agrif_Update_vvl ***
171      !!---------------------------------------------
172      !
173      IF (Agrif_Root()) RETURN
174      !
175#if defined TWO_WAY 
176      !
177      IF (lwp.AND.lk_agrif_debug) Write(*,*) 'Update e3 from grid Number',Agrif_Fixed(), 'Step', Agrif_Nb_Step()
178      !
179      Agrif_UseSpecialValueInUpdate = .FALSE.
180      Agrif_SpecialValueFineGrid = 0.
181      !
182# if ! defined DECAL_FEEDBACK
183      CALL Agrif_Update_Variable(e3t_id, procname=updatee3t)
184# else
185      CALL Agrif_Update_Variable(e3t_id, locupdate=(/1,0/), procname=updatee3t)
186# endif 
187      !
188      Agrif_UseSpecialValueInUpdate = .FALSE.
189      !
190      CALL Agrif_ChildGrid_To_ParentGrid()
191      CALL dom_vvl_update_UVF
192      CALL Agrif_ParentGrid_To_ChildGrid()
193      !
194#endif
195      !
196   END SUBROUTINE Agrif_Update_vvl
197
198   SUBROUTINE dom_vvl_update_UVF
199      !!---------------------------------------------
200      !!       *** ROUTINE dom_vvl_update_UVF ***
201      !!---------------------------------------------
202#  include "domzgr_substitute.h90"
203      !!
204      INTEGER :: jk
205      REAL(wp):: zcoef
206      !!---------------------------------------------
207
208      IF (lwp.AND.lk_agrif_debug) Write(*,*) 'Finalize e3 on grid Number', &
209                  & Agrif_Fixed(), 'Step', Agrif_Nb_Step()
210
211      ! Save "old" scale factor (prior update) for subsequent asselin correction
212      ! of prognostic variables (needed to update initial state only)
213      fse3u_a(:,:,:) = fse3u_n(:,:,:)
214      fse3v_a(:,:,:) = fse3v_n(:,:,:)
215!      ua(:,:,:) = fse3u_b(:,:,:)
216!      va(:,:,:) = fse3v_b(:,:,:)
217      hu_a(:,:) = hu(:,:)
218      hv_a(:,:) = hv(:,:)
219
220      ! Vertical scale factor interpolations
221      ! ------------------------------------
222      !
223      CALL dom_vvl_interpol( fse3t_n(:,:,:), fse3u_n(:,:,:),  'U'  )
224      CALL dom_vvl_interpol( fse3t_n(:,:,:), fse3v_n(:,:,:),  'V'  )
225      CALL dom_vvl_interpol( fse3u_n(:,:,:), fse3f_n(:,:,:),  'F'  )
226
227      CALL dom_vvl_interpol( fse3u_n(:,:,:), fse3uw_n(:,:,:), 'UW' )
228      CALL dom_vvl_interpol( fse3v_n(:,:,:), fse3vw_n(:,:,:), 'VW' )
229
230      ! Update total depths:
231      hu(:,:) = 0._wp                        ! Ocean depth at U-points
232      hv(:,:) = 0._wp                        ! Ocean depth at V-points
233      DO jk = 1, jpkm1
234         hu(:,:) = hu(:,:) + fse3u_n(:,:,jk) * umask(:,:,jk)
235         hv(:,:) = hv(:,:) + fse3v_n(:,:,jk) * vmask(:,:,jk)
236      END DO
237      !                                      ! Inverse of the local depth
238      hur(:,:) = 1._wp / ( hu(:,:) + 1._wp - umask_i(:,:) ) * umask_i(:,:)
239      hvr(:,:) = 1._wp / ( hv(:,:) + 1._wp - vmask_i(:,:) ) * vmask_i(:,:)
240
241#if ! defined key_dynspg_ts
242      IF  (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN
243#else
244      IF ((.NOT.(lk_agrif_fstep.AND.(neuler==0))).AND.(.NOT.ln_bt_fw)) THEN
245#endif
246         CALL dom_vvl_interpol( fse3t_b(:,:,:), fse3u_b(:,:,:),  'U'  )
247         CALL dom_vvl_interpol( fse3t_b(:,:,:), fse3v_b(:,:,:),  'V'  )
248
249         CALL dom_vvl_interpol( fse3u_b(:,:,:), fse3uw_b(:,:,:), 'UW' )
250         CALL dom_vvl_interpol( fse3v_b(:,:,:), fse3vw_b(:,:,:), 'VW' )
251
252         hu_b(:,:) = 0._wp                        ! Ocean depth at U-points
253         hv_b(:,:) = 0._wp                        ! Ocean depth at V-points
254         DO jk = 1, jpkm1
255            hu_b(:,:) = hu_b(:,:) + fse3u_b(:,:,jk) * umask(:,:,jk)
256            hv_b(:,:) = hv_b(:,:) + fse3v_b(:,:,jk) * vmask(:,:,jk)
257         END DO
258         !                                      ! Inverse of the local depth
259         hur_b(:,:) = 1._wp / ( hu_b(:,:) + 1._wp - umask_i(:,:) ) * umask_i(:,:)
260         hvr_b(:,:) = 1._wp / ( hv_b(:,:) + 1._wp - vmask_i(:,:) ) * vmask_i(:,:)
261      ENDIF
262
263      !
264   END SUBROUTINE dom_vvl_update_UVF
265
266   SUBROUTINE updateTS( tabres, i1, i2, j1, j2, k1, k2, n1, n2, before )
267      !!---------------------------------------------
268      !!           *** ROUTINE updateT ***
269      !!---------------------------------------------
270      INTEGER, INTENT(in) :: i1,i2,j1,j2,k1,k2,n1,n2
271      REAL(wp),DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) :: tabres
272      LOGICAL, INTENT(in) :: before
273      !!
274      INTEGER :: ji,jj,jk,jn
275      REAL(wp) :: ztb, ztnu, ztno
276      !!---------------------------------------------
277      !
278      !
279      IF (before) THEN
280         DO jn = n1,n2
281            DO jk=k1,k2
282               DO jj=j1,j2
283                  DO ji=i1,i2
284!> jc tmp
285                     tabres(ji,jj,jk,jn) = tsn(ji,jj,jk,jn)  * fse3t_n(ji,jj,jk) / e3t_0(ji,jj,jk)
286!                     tabres(ji,jj,jk,jn) = tsn(ji,jj,jk,jn)  * fse3t_n(ji,jj,jk)
287!< jc tmp
288                  END DO
289               END DO
290            END DO
291         END DO
292      ELSE
293!> jc tmp
294         DO jn = n1,n2
295            tabres(i1:i2,j1:j2,k1:k2,jn) =  tabres(i1:i2,j1:j2,k1:k2,jn) * e3t_0(i1:i2,j1:j2,k1:k2)
296         ENDDO
297!< jc tmp
298
299         IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN
300            ! Add asselin part
301            DO jn = n1,n2
302               DO jk=k1,k2
303                  DO jj=j1,j2
304                     DO ji=i1,i2
305                        IF( tabres(ji,jj,jk,jn) .NE. 0. ) THEN
306                           ztb  = tsb(ji,jj,jk,jn) * fse3t_b(ji,jj,jk) ! fse3t_b prior update should be used
307                           ztnu = tabres(ji,jj,jk,jn)
308                           ztno = tsn(ji,jj,jk,jn) * fse3t_a(ji,jj,jk)
309                           tsb(ji,jj,jk,jn) = ( ztb + atfp * ( ztnu - ztno) )  & 
310                                     &        * tmask(ji,jj,jk) / fse3t_b(ji,jj,jk)
311                        ENDIF
312                     ENDDO
313                  ENDDO
314               ENDDO
315            ENDDO
316         ENDIF
317
318         DO jn = n1,n2
319            DO jk=k1,k2
320               DO jj=j1,j2
321                  DO ji=i1,i2
322                     IF( tabres(ji,jj,jk,jn) .NE. 0. ) THEN
323                        tsn(ji,jj,jk,jn) = tabres(ji,jj,jk,jn) / fse3t_n(ji,jj,jk)
324                     END IF
325                  END DO
326               END DO
327            END DO
328         END DO
329         !
330         IF  ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN
331            tsb(i1:i2,j1:j2,k1:k2,n1:n2)  = tsn(i1:i2,j1:j2,k1:k2,n1:n2)
332         ENDIF
333         !
334      ENDIF
335      !
336   END SUBROUTINE updateTS
337
338   SUBROUTINE updateu( tabres, i1, i2, j1, j2, k1, k2, before, nb, ndir )
339      !!---------------------------------------------
340      !!           *** ROUTINE updateu ***
341      !!---------------------------------------------
342#  include "domzgr_substitute.h90"
343      !!
344      INTEGER, INTENT(in) :: i1, i2, j1, j2, k1, k2
345      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: tabres
346      LOGICAL, INTENT(in) :: before
347      INTEGER, INTENT(in) :: nb , ndir
348      !!
349      LOGICAL western_side, eastern_side
350      INTEGER :: ji, jj, jk
351      REAL(wp) :: zrhoy
352      REAL(wp) :: zub, zunu, zuno
353      !!---------------------------------------------
354      !
355      IF (before) THEN
356         zrhoy = Agrif_Rhoy()
357         DO jk=k1,k2
358            DO jj=j1,j2
359               DO ji=i1,i2
360                  tabres(ji,jj,jk) = e2u(ji,jj) * un(ji,jj,jk)
361                  tabres(ji,jj,jk) = tabres(ji,jj,jk) * fse3u_n(ji,jj,jk)
362               END DO
363            END DO
364         END DO
365         tabres = zrhoy * tabres
366      ELSE
367         western_side = (nb == 1).AND.(ndir == 1)
368         eastern_side = (nb == 1).AND.(ndir == 2)
369         DO jk=k1,k2
370            DO jj=j1,j2
371               DO ji=i1,i2
372                  tabres(ji,jj,jk) = tabres(ji,jj,jk) / e2u(ji,jj) 
373                  !
374                  IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN ! Add asselin part
375                     zub  = ub(ji,jj,jk) * fse3u_b(ji,jj,jk)
376                     zuno = un(ji,jj,jk) * fse3u_a(ji,jj,jk)
377                     zunu = tabres(ji,jj,jk)
378                     ub(ji,jj,jk) = ( zub + atfp * ( zunu - zuno) ) &     
379                                    & * umask(ji,jj,jk) / fse3u_b(ji,jj,jk)
380                  ENDIF
381                  !
382                  un(ji,jj,jk) = tabres(ji,jj,jk) * umask(ji,jj,jk) / fse3u_n(ji,jj,jk)
383               END DO
384            END DO
385         END DO
386         !
387!         IF (western_side) THEN
388!            DO jk=k1,k2
389!               DO jj=j1,j2
390!                 un(i1-1,jj,jk) = un(i1-1,jj,jk) * fse3u_a(i1-1,jj,jk) / fse3u_n(i1-1,jj,jk)
391!               END DO
392!            ENDDO
393!         ENDIF
394!         IF (eastern_side) THEN
395!            DO jk=k1,k2
396!               DO jj=j1,j2
397!                 un(i2+1,jj,jk) = un(i2+1,jj,jk) * fse3u_a(i2+1,jj,jk) / fse3u_n(i2+1,jj,jk)
398!               END DO
399!            ENDDO
400!         ENDIF
401         !
402         IF  ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN
403            ub(i1:i2,j1:j2,k1:k2)  = un(i1:i2,j1:j2,k1:k2)
404         ENDIF
405         !
406      ENDIF
407      !
408   END SUBROUTINE updateu
409
410   SUBROUTINE updatev( tabres, i1, i2, j1, j2, k1, k2, before, nb, ndir )
411      !!---------------------------------------------
412      !!           *** ROUTINE updatev ***
413      !!---------------------------------------------
414#  include "domzgr_substitute.h90"
415      !!
416      INTEGER :: i1,i2,j1,j2,k1,k2
417      INTEGER :: ji,jj,jk
418      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2) :: tabres
419      LOGICAL :: before
420      INTEGER, INTENT(in) :: nb , ndir
421      !!
422      LOGICAL :: northern_side, southern_side
423      REAL(wp) :: zrhox
424      REAL(wp) :: zvb, zvnu, zvno
425      !!---------------------------------------------     
426      !
427      IF (before) THEN
428         zrhox = Agrif_Rhox()
429         DO jk=k1,k2
430            DO jj=j1,j2
431               DO ji=i1,i2
432                  tabres(ji,jj,jk) = e1v(ji,jj) * vn(ji,jj,jk)
433                  tabres(ji,jj,jk) = tabres(ji,jj,jk) * fse3v_n(ji,jj,jk)
434               END DO
435            END DO
436         END DO
437         tabres = zrhox * tabres
438      ELSE
439         southern_side = (nb == 2).AND.(ndir == 1)
440         northern_side = (nb == 2).AND.(ndir == 2)
441         DO jk=k1,k2
442            DO jj=j1,j2
443               DO ji=i1,i2
444                  tabres(ji,jj,jk) = tabres(ji,jj,jk) / e1v(ji,jj)
445                  !
446                  IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN ! Add asselin part
447                     zvb  = vb(ji,jj,jk) * fse3v_b(ji,jj,jk)
448                     zvno = vn(ji,jj,jk) * fse3v_a(ji,jj,jk)
449                     zvnu = tabres(ji,jj,jk)
450                     vb(ji,jj,jk) = ( zvb + atfp * ( zvnu - zvno) ) &     
451                                    & * vmask(ji,jj,jk) / fse3v_b(ji,jj,jk)
452                  ENDIF
453                  !
454                  vn(ji,jj,jk) = tabres(ji,jj,jk) * vmask(ji,jj,jk) / fse3v_n(ji,jj,jk)
455               END DO
456            END DO
457         END DO
458         !
459!         IF (southern_side) THEN
460!            DO jk=k1,k2
461!               DO ji=i1,i2
462!                 vn(ji,j1-1,jk) = vn(ji,j1-1,jk) * fse3v_a(ji,j1-1,jk) / fse3v_n(ji,j1-1,jk)
463!               END DO
464!            ENDDO
465!         ENDIF
466!         IF (northern_side) THEN
467!            DO jk=k1,k2
468!               DO ji=i1,i2
469!                 vn(ji,j2+1,jk) = vn(ji,j2+1,jk) * fse3v_a(ji,j2+1,jk) / fse3v_n(ji,j2+1,jk)
470!               END DO
471!            ENDDO
472!         ENDIF
473         !
474         IF  ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN
475            vb(i1:i2,j1:j2,k1:k2)  = vn(i1:i2,j1:j2,k1:k2)
476         ENDIF
477         !
478      ENDIF
479      !
480   END SUBROUTINE updatev
481
482   SUBROUTINE updateu2d( tabres, i1, i2, j1, j2, before, nb, ndir )
483      !!---------------------------------------------
484      !!          *** ROUTINE updateu2d ***
485      !!---------------------------------------------
486#  include "domzgr_substitute.h90"
487      !!
488      INTEGER, INTENT(in) :: i1, i2, j1, j2
489      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres
490      LOGICAL, INTENT(in) :: before
491      INTEGER, INTENT(in) :: nb , ndir
492      !!
493      LOGICAL western_side, eastern_side 
494      INTEGER :: ji, jj, jk
495      REAL(wp) :: zrhoy
496      REAL(wp) :: zcorr
497      !!---------------------------------------------
498      !
499      IF (before) THEN
500         zrhoy = Agrif_Rhoy()
501         DO jj=j1,j2
502            DO ji=i1,i2
503               tabres(ji,jj) = un_b(ji,jj) * hu(ji,jj) * e2u(ji,jj)
504            END DO
505         END DO
506         tabres = zrhoy * tabres
507      ELSE
508         western_side = (nb == 1).AND.(ndir == 1)
509         eastern_side = (nb == 1).AND.(ndir == 2)
510         DO jj=j1,j2
511            DO ji=i1,i2
512               tabres(ji,jj) =  tabres(ji,jj) / e2u(ji,jj) 
513               !   
514               ! Update "now" 3d velocities:
515               spgu(ji,jj) = 0.e0
516               DO jk=1,jpkm1
517                  spgu(ji,jj) = spgu(ji,jj) + fse3u_n(ji,jj,jk) * un(ji,jj,jk)
518               END DO
519               !
520               zcorr = (tabres(ji,jj) - spgu(ji,jj)) * hur(ji,jj)
521               DO jk=1,jpkm1             
522                  un(ji,jj,jk) = un(ji,jj,jk) + zcorr * umask(ji,jj,jk)           
523               END DO
524               !
525               ! Update barotropic velocities:
526#if ! defined key_dynspg_ts
527               IF  (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN
528#else
529               IF ((.NOT.(lk_agrif_fstep.AND.(neuler==0))).AND.(.NOT.ln_bt_fw)) THEN
530#endif
531                  zcorr = (tabres(ji,jj) - un_b(ji,jj) * hu_a(ji,jj)) * hur_b(ji,jj)
532                  ub_b(ji,jj) = ub_b(ji,jj) + atfp * zcorr * umask(ji,jj,1)
533               END IF             
534               un_b(ji,jj) = tabres(ji,jj) * hur(ji,jj) * umask(ji,jj,1)
535               !       
536               ! Correct "before" velocities to hold correct bt component:
537               spgu(ji,jj) = 0.e0
538               DO jk=1,jpkm1
539                  spgu(ji,jj) = spgu(ji,jj) + fse3u_b(ji,jj,jk) * ub(ji,jj,jk)
540               END DO
541               !
542               zcorr = ub_b(ji,jj) - spgu(ji,jj) * hur_b(ji,jj)
543               DO jk=1,jpkm1             
544                  ub(ji,jj,jk) = ub(ji,jj,jk) + zcorr * umask(ji,jj,jk)           
545               END DO
546               !
547            END DO
548         END DO
549!         IF (western_side) THEN
550!            DO jj=j1,j2
551!              un_b(i1-1,jj) = un_b(i1-1,jj) * hu_a(i1-1,jj) * hur(i1-1,jj)
552!            END DO
553!         ENDIF
554!         IF (eastern_side) THEN
555!            DO jj=j1,j2
556!              un_b(i2+1,jj) = un_b(i2+1,jj) * hu_a(i2+1,jj) * hur(i2+1,jj)
557!            ENDDO
558!         ENDIF
559         IF  ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN
560            ub_b(i1:i2,j1:j2)  = un_b(i1:i2,j1:j2)
561         ENDIF
562      ENDIF
563      !
564   END SUBROUTINE updateu2d
565
566   SUBROUTINE updatev2d( tabres, i1, i2, j1, j2, before, nb, ndir  )
567      !!---------------------------------------------
568      !!          *** ROUTINE updatev2d ***
569      !!---------------------------------------------
570      INTEGER, INTENT(in) :: i1, i2, j1, j2
571      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres
572      LOGICAL, INTENT(in) :: before
573      INTEGER, INTENT(in) :: nb , ndir
574      !!
575      LOGICAL :: northern_side, southern_side
576      INTEGER :: ji, jj, jk
577      REAL(wp) :: zrhox
578      REAL(wp) :: zcorr
579      !!---------------------------------------------
580      !
581      IF (before) THEN
582         zrhox = Agrif_Rhox()
583         DO jj=j1,j2
584            DO ji=i1,i2
585               tabres(ji,jj) = vn_b(ji,jj) * hv(ji,jj) * e1v(ji,jj) 
586            END DO
587         END DO
588         tabres = zrhox * tabres
589      ELSE
590         southern_side = (nb == 2).AND.(ndir == 1)
591         northern_side = (nb == 2).AND.(ndir == 2)
592         DO jj=j1,j2
593            DO ji=i1,i2
594               tabres(ji,jj) =  tabres(ji,jj) / e1v(ji,jj) 
595               !   
596               ! Update "now" 3d velocities:
597               spgv(ji,jj) = 0.e0
598               DO jk=1,jpkm1
599                  spgv(ji,jj) = spgv(ji,jj) + fse3v_n(ji,jj,jk) * vn(ji,jj,jk)
600               END DO
601               !
602               zcorr = (tabres(ji,jj) - spgv(ji,jj)) * hvr(ji,jj)
603               DO jk=1,jpkm1             
604                  vn(ji,jj,jk) = vn(ji,jj,jk) + zcorr * vmask(ji,jj,jk)           
605               END DO
606               !
607               ! Update barotropic velocities:
608#if ! defined key_dynspg_ts
609               IF  (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN
610#else
611               IF ((.NOT.(lk_agrif_fstep.AND.(neuler==0))).AND.(.NOT.ln_bt_fw)) THEN
612#endif
613                  zcorr = (tabres(ji,jj) - vn_b(ji,jj)*hv_a(ji,jj)) * hvr_b(ji,jj)
614                  vb_b(ji,jj) = vb_b(ji,jj) + atfp * zcorr * vmask(ji,jj,1)
615               END IF             
616               vn_b(ji,jj) = tabres(ji,jj) * hvr(ji,jj) * vmask(ji,jj,1)
617               !       
618               ! Correct "before" velocities to hold correct bt component:
619               spgv(ji,jj) = 0.e0
620               DO jk=1,jpkm1
621                  spgv(ji,jj) = spgv(ji,jj) + fse3v_b(ji,jj,jk) * vb(ji,jj,jk)
622               END DO
623               !
624               zcorr = vb_b(ji,jj) - spgv(ji,jj) * hvr_b(ji,jj)
625               DO jk=1,jpkm1             
626                  vb(ji,jj,jk) = vb(ji,jj,jk) + zcorr * vmask(ji,jj,jk)           
627               END DO
628               !
629            END DO
630         END DO
631         !
632!         IF (southern_side) THEN
633!            DO ji=i1,i2
634!              vn_b(ji,j1-1) = vn_b(ji,j1-1) * hv_a(ji,j1-1) * hvr(ji,j1-1)
635!            END DO
636!         ENDIF
637!         IF (northern_side) THEN
638!            DO ji=i1,i2
639!              vn_b(ji,j2+1) = vn_b(ji,j2+1) * hv_a(ji,j2+1) * hvr(ji,j2+1)
640!            END DO
641!         ENDIF
642         !
643         IF  ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN
644            vb_b(i1:i2,j1:j2)  = vn_b(i1:i2,j1:j2)
645         ENDIF
646         !
647      ENDIF
648      !
649   END SUBROUTINE updatev2d
650
651
652   SUBROUTINE updateSSH( tabres, i1, i2, j1, j2, before )
653      !!---------------------------------------------
654      !!          *** ROUTINE updateSSH ***
655      !!---------------------------------------------
656      INTEGER, INTENT(in) :: i1, i2, j1, j2
657      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres
658      LOGICAL, INTENT(in) :: before
659      !!
660      INTEGER :: ji, jj
661      !!---------------------------------------------
662      !
663      IF (before) THEN
664         DO jj=j1,j2
665            DO ji=i1,i2
666               tabres(ji,jj) = sshn(ji,jj)
667            END DO
668         END DO
669      ELSE
670#if ! defined key_dynspg_ts
671         IF  (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN
672#else
673         IF ((.NOT.(lk_agrif_fstep.AND.(neuler==0))).AND.(.NOT.ln_bt_fw)) THEN
674#endif
675            DO jj=j1,j2
676               DO ji=i1,i2
677                  sshb(ji,jj) =   sshb(ji,jj) &
678                        & + atfp * ( tabres(ji,jj) - sshn(ji,jj) ) * tmask(ji,jj,1)
679               END DO
680            END DO
681         ENDIF
682         !
683         DO jj=j1,j2
684            DO ji=i1,i2
685               sshn(ji,jj) = tabres(ji,jj) * tmask(ji,jj,1)
686            END DO
687         END DO
688         !
689         IF  ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN
690            sshb(i1:i2,j1:j2)  = sshn(i1:i2,j1:j2)
691         ENDIF
692         !
693      ENDIF
694      !
695   END SUBROUTINE updateSSH
696
697   SUBROUTINE updateub2b( tabres, i1, i2, j1, j2, before )
698      !!---------------------------------------------
699      !!          *** ROUTINE updateub2b ***
700      !!---------------------------------------------
701      INTEGER, INTENT(in) :: i1, i2, j1, j2
702      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres
703      LOGICAL, INTENT(in) :: before
704      !!
705      INTEGER :: ji, jj
706      REAL(wp) :: zrhoy, za1
707      !!---------------------------------------------
708      !
709      IF (before) THEN
710         zrhoy = Agrif_Rhoy()
711         DO jj=j1,j2
712            DO ji=i1,i2
713               tabres(ji,jj) = ub2_i_b(ji,jj) * e2u(ji,jj)
714            END DO
715         END DO
716         tabres = zrhoy * tabres
717      ELSE
718         za1 = 1._wp / REAL(Agrif_rhot(), wp)
719         tabres(i1:i2,j1:j2) = tabres(i1:i2,j1:j2) / e2u(i1:i2,j1:j2)
720         DO jj=j1,j2
721            DO ji=i1,i2 
722               ub2_i_b(ji,jj) = ub2_i_b(ji,jj) & 
723                & + za1 * (tabres(ji,jj) - ub2_b(ji,jj))
724               ub2_b(ji,jj) = tabres(ji,jj)
725            END DO
726         END DO
727      ENDIF
728      !
729   END SUBROUTINE updateub2b
730
731   SUBROUTINE updatevb2b( tabres, i1, i2, j1, j2, before )
732      !!---------------------------------------------
733      !!          *** ROUTINE updatevb2b ***
734      !!---------------------------------------------
735      INTEGER, INTENT(in) :: i1, i2, j1, j2
736      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres
737      LOGICAL, INTENT(in) :: before
738      !!
739      INTEGER :: ji, jj
740      REAL(wp) :: zrhox, za1
741      !!---------------------------------------------
742      !
743      IF (before) THEN
744         zrhox = Agrif_Rhox()
745         DO jj=j1,j2
746            DO ji=i1,i2
747               tabres(ji,jj) = vb2_i_b(ji,jj) * e1v(ji,jj) 
748            END DO
749         END DO
750         tabres = zrhox * tabres
751      ELSE
752         za1 = 1._wp / REAL(Agrif_rhot(), wp)
753         tabres(i1:i2,j1:j2) = tabres(i1:i2,j1:j2) / e1v(i1:i2,j1:j2)
754         DO jj=j1,j2
755            DO ji=i1,i2
756               vb2_i_b(ji,jj) = vb2_i_b(ji,jj) & 
757                & + za1 * (tabres(ji,jj) - vb2_b(ji,jj))
758               vb2_b(ji,jj) = tabres(ji,jj)
759            END DO
760         END DO
761      ENDIF
762      !
763   END SUBROUTINE updatevb2b
764
765
766   SUBROUTINE update_scales( tabres, i1, i2, j1, j2, k1, k2, n1,n2, before )
767      ! currently not used
768      !!---------------------------------------------
769      !!           *** ROUTINE updateT ***
770      !!---------------------------------------------
771#  include "domzgr_substitute.h90"
772
773      INTEGER, INTENT(in) :: i1,i2,j1,j2,k1,k2,n1,n2
774      REAL(wp),DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) :: tabres
775      LOGICAL, iNTENT(in) :: before
776
777      INTEGER :: ji,jj,jk
778      REAL(wp) :: ztemp
779
780      IF (before) THEN
781         DO jk=k1,k2
782            DO jj=j1,j2
783               DO ji=i1,i2
784                  tabres(ji,jj,jk,1) = e1t(ji,jj)*e2t(ji,jj)*tmask(ji,jj,jk)
785                  tabres(ji,jj,jk,2) = e1t(ji,jj)*tmask(ji,jj,jk)
786                  tabres(ji,jj,jk,3) = e2t(ji,jj)*tmask(ji,jj,jk)
787               END DO
788            END DO
789         END DO
790         tabres(:,:,:,1)=tabres(:,:,:,1)*Agrif_Rhox()*Agrif_Rhoy()
791         tabres(:,:,:,2)=tabres(:,:,:,2)*Agrif_Rhox()
792         tabres(:,:,:,3)=tabres(:,:,:,3)*Agrif_Rhoy()
793      ELSE
794         DO jk=k1,k2
795            DO jj=j1,j2
796               DO ji=i1,i2
797                  IF( tabres(ji,jj,jk,1) .NE. 0. ) THEN
798                     print *,'VAL = ',ji,jj,jk,tabres(ji,jj,jk,1),e1t(ji,jj)*e2t(ji,jj)*tmask(ji,jj,jk)
799                     print *,'VAL2 = ',ji,jj,jk,tabres(ji,jj,jk,2),e1t(ji,jj)*tmask(ji,jj,jk)
800                     print *,'VAL3 = ',ji,jj,jk,tabres(ji,jj,jk,3),e2t(ji,jj)*tmask(ji,jj,jk)
801                     ztemp = sqrt(tabres(ji,jj,jk,1)/(tabres(ji,jj,jk,2)*tabres(ji,jj,jk,3)))
802                     print *,'CORR = ',ztemp-1.
803                     print *,'NEW VALS = ',tabres(ji,jj,jk,2)*ztemp,tabres(ji,jj,jk,3)*ztemp, &
804                           tabres(ji,jj,jk,2)*ztemp*tabres(ji,jj,jk,3)*ztemp
805                     e1t(ji,jj) = tabres(ji,jj,jk,2)*ztemp
806                     e2t(ji,jj) = tabres(ji,jj,jk,3)*ztemp
807                  END IF
808               END DO
809            END DO
810         END DO
811      ENDIF
812      !
813   END SUBROUTINE update_scales
814
815# if defined key_zdftke
816   SUBROUTINE updateEN( ptab, i1, i2, j1, j2, k1, k2, before )
817      !!---------------------------------------------
818      !!           *** ROUTINE updateen ***
819      !!---------------------------------------------
820      INTEGER, INTENT(in) :: i1, i2, j1, j2, k1, k2
821      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: ptab
822      LOGICAL, INTENT(in) :: before
823      !!---------------------------------------------
824      !
825      IF (before) THEN
826         ptab (i1:i2,j1:j2,k1:k2) = en(i1:i2,j1:j2,k1:k2)
827      ELSE
828         en(i1:i2,j1:j2,k1:k2) = ptab (i1:i2,j1:j2,k1:k2) 
829      ENDIF
830      !
831   END SUBROUTINE updateEN
832
833
834   SUBROUTINE updateAVT( ptab, i1, i2, j1, j2, k1, k2, before )
835      !!---------------------------------------------
836      !!           *** ROUTINE updateavt ***
837      !!---------------------------------------------
838      INTEGER, INTENT(in) :: i1, i2, j1, j2, k1, k2
839      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: ptab
840      LOGICAL, INTENT(in) :: before
841      !!---------------------------------------------
842      !
843      IF (before) THEN
844         ptab (i1:i2,j1:j2,k1:k2) = avt_k(i1:i2,j1:j2,k1:k2)
845      ELSE
846         avt_k(i1:i2,j1:j2,k1:k2) = ptab (i1:i2,j1:j2,k1:k2) 
847      ENDIF
848      !
849   END SUBROUTINE updateAVT
850
851
852   SUBROUTINE updateAVM( ptab, i1, i2, j1, j2, k1, k2, before )
853      !!---------------------------------------------
854      !!           *** ROUTINE updateavm ***
855      !!---------------------------------------------
856      INTEGER, INTENT(in) :: i1, i2, j1, j2, k1, k2
857      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: ptab
858      LOGICAL, INTENT(in) :: before
859      !!---------------------------------------------
860      !
861      IF (before) THEN
862         ptab (i1:i2,j1:j2,k1:k2) = avm_k(i1:i2,j1:j2,k1:k2)
863      ELSE
864         avm_k(i1:i2,j1:j2,k1:k2) = ptab (i1:i2,j1:j2,k1:k2) 
865      ENDIF
866      !
867   END SUBROUTINE updateAVM
868
869# endif /* key_zdftke */ 
870
871   SUBROUTINE updatee3t( ptab, i1, i2, j1, j2, k1, k2, before )
872      !!---------------------------------------------
873      !!           *** ROUTINE updatee3t ***
874      !!---------------------------------------------
875      INTEGER, INTENT(in) :: i1, i2, j1, j2, k1, k2
876      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: ptab
877      LOGICAL, INTENT(in) :: before
878      INTEGER :: ji,jj,jk
879      REAL(wp) :: zcoef
880      !!---------------------------------------------
881      !
882      IF (before) THEN
883!> jc tmp:
884!         ptab(i1:i2,j1:j2,k1:k2) = fse3t_n(i1:i2,j1:j2,k1:k2)
885         ptab(i1:i2,j1:j2,k1:k2) = fse3t_n(i1:i2,j1:j2,k1:k2) / e3t_0(i1:i2,j1:j2,k1:k2)
886!< jc tmp:
887      ELSE
888         !
889         ! 1) Updates at before time step:
890         ! -------------------------------
891         !
892!> jc tmp:
893         ptab(i1:i2,j1:j2,k1:k2) = ptab(i1:i2,j1:j2,k1:k2) * e3t_0(i1:i2,j1:j2,k1:k2)
894!< jc tmp:
895
896         ! Save "old" scale factor (prior update) for subsequent asselin correction
897         ! of prognostic variables (needed to update initial state only)
898         fse3t_a(i1:i2,j1:j2,k1:k2) = fse3t_n(i1:i2,j1:j2,k1:k2)
899!         hdivb(i1:i2,j1:j2,k1:k2)   = fse3t_b(i1:i2,j1:j2,k1:k2)
900
901#if ! defined key_dynspg_ts
902         IF  (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN
903#else
904         IF ((.NOT.(lk_agrif_fstep.AND.(neuler==0))).AND.(.NOT.ln_bt_fw)) THEN
905#endif
906            DO jk = 1, jpkm1
907               DO jj=j1,j2
908                  DO ji=i1,i2
909                     fse3t_b(ji,jj,jk) =   fse3t_b(ji,jj,jk) &
910                           & + atfp * ( ptab(ji,jj,jk) - fse3t_n(ji,jj,jk) )
911                  END DO
912               END DO
913            END DO
914            !
915            fse3w_b (i1:i2,j1:j2,1) = e3w_0(i1:i2,j1:j2,1) + fse3t_b(i1:i2,j1:j2,1) - e3t_0(i1:i2,j1:j2,1)
916            fsdepw_b(i1:i2,j1:j2,1) = 0.0_wp
917            fsdept_b(i1:i2,j1:j2,1) = 0.5_wp * fse3w_b(i1:i2,j1:j2,1)
918            !
919            DO jk = 2, jpk
920               DO jj = j1,j2
921                  DO ji = i1,i2           
922                     zcoef = (tmask(ji,jj,jk) - wmask(ji,jj,jk))
923                     fse3w_b(ji,jj,jk)  = e3w_0(ji,jj,jk) + ( 1.0_wp - 0.5_wp * tmask(ji,jj,jk) ) *        & 
924                     &                                        ( fse3t_b(ji,jj,jk-1) - e3t_0(ji,jj,jk-1) )  &
925                     &                                    +            0.5_wp * tmask(ji,jj,jk)   *        &
926                     &                                        ( fse3t_b(ji,jj,jk  ) - e3t_0(ji,jj,jk  ) )
927                     fsdepw_b(ji,jj,jk) = fsdepw_b(ji,jj,jk-1) + fse3t_b(ji,jj,jk-1)
928                     fsdept_b(ji,jj,jk) =      zcoef  * ( fsdepw_b(ji,jj,jk  ) + 0.5 * fse3w_b(ji,jj,jk))  &
929                         &                + (1-zcoef) * ( fsdept_b(ji,jj,jk-1) +       fse3w_b(ji,jj,jk)) 
930                  END DO
931               END DO
932            END DO
933            !
934         ENDIF       
935         !
936         ! 2) Updates at now time step:
937         ! ----------------------------
938         !
939         ! Update vertical scale factor at T-points:
940         fse3t_n(i1:i2,j1:j2,k1:k2) = ptab(i1:i2,j1:j2,k1:k2)
941         !
942         ! Update total depth:
943         ht(i1:i2,j1:j2) = 0._wp
944         DO jk = 1, jpkm1
945            ht(i1:i2,j1:j2) = ht(i1:i2,j1:j2) + fse3t_n(i1:i2,j1:j2,jk) * tmask(i1:i2,j1:j2,jk)
946         END DO
947         !
948         ! Update vertical scale factor at W-points and depths:
949         fse3w_n (i1:i2,j1:j2,1) = e3w_0(i1:i2,j1:j2,1) + fse3t_n(i1:i2,j1:j2,1) - e3t_0(i1:i2,j1:j2,1)
950         fsdept_n(i1:i2,j1:j2,1) = 0.5_wp * fse3w_n(i1:i2,j1:j2,1)
951         fsdepw_n(i1:i2,j1:j2,1) = 0.0_wp
952         fsde3w_n(i1:i2,j1:j2,1) = fsdept_n(i1:i2,j1:j2,1) - (ht(i1:i2,j1:j2)-ht_0(i1:i2,j1:j2)) ! Last term in the rhs is ssh
953         !
954         DO jk = 2, jpk
955            DO jj = j1,j2
956               DO ji = i1,i2           
957               zcoef = (tmask(ji,jj,jk) - wmask(ji,jj,jk))
958               fse3w_n(ji,jj,jk)  = e3w_0(ji,jj,jk) + ( 1.0_wp - 0.5_wp * tmask(ji,jj,jk) ) * ( fse3t_n(ji,jj,jk-1) - e3t_0(ji,jj,jk-1) )   &
959               &                                    +            0.5_wp * tmask(ji,jj,jk)   * ( fse3t_n(ji,jj,jk  ) - e3t_0(ji,jj,jk  ) )
960               fsdepw_n(ji,jj,jk) = fsdepw_n(ji,jj,jk-1) + fse3t_n(ji,jj,jk-1)
961               fsdept_n(ji,jj,jk) =      zcoef  * ( fsdepw_n(ji,jj,jk  ) + 0.5 * fse3w_n(ji,jj,jk))  &
962                   &                + (1-zcoef) * ( fsdept_n(ji,jj,jk-1) +       fse3w_n(ji,jj,jk)) 
963               fsde3w_n(ji,jj,jk) = fsdept_n(ji,jj,jk) - (ht(ji,jj)-ht_0(ji,jj)) ! Last term in the rhs is ssh
964               END DO
965            END DO
966         END DO
967         !
968         IF  ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN
969            fse3t_b (i1:i2,j1:j2,1:jpk)  = fse3t_n (i1:i2,j1:j2,1:jpk)
970            fse3w_b (i1:i2,j1:j2,1:jpk)  = fse3w_n (i1:i2,j1:j2,1:jpk)
971            fsdepw_b(i1:i2,j1:j2,1:jpk)  = fsdepw_n(i1:i2,j1:j2,1:jpk)
972            fsdept_b(i1:i2,j1:j2,1:jpk)  = fsdept_n(i1:i2,j1:j2,1:jpk)
973         ENDIF
974         !
975      ENDIF
976      !
977   END SUBROUTINE updatee3t
978
979#else
980CONTAINS
981   SUBROUTINE agrif_opa_update_empty
982      !!---------------------------------------------
983      !!   *** ROUTINE agrif_opa_update_empty ***
984      !!---------------------------------------------
985      WRITE(*,*)  'agrif_opa_update : You should not have seen this print! error?'
986   END SUBROUTINE agrif_opa_update_empty
987#endif
988END MODULE agrif_opa_update
989
Note: See TracBrowser for help on using the repository browser.