New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
agrif_oce_interp.F90 in NEMO/branches/2020/dev_r13312_AGRIF-03-04_jchanut_vinterp_tstep/src/NST – NEMO

source: NEMO/branches/2020/dev_r13312_AGRIF-03-04_jchanut_vinterp_tstep/src/NST/agrif_oce_interp.F90 @ 14013

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

Manually merge with Christian's branch agrif top sponge/interp #2129

  • Property svn:keywords set to Id
File size: 73.6 KB
Line 
1MODULE agrif_oce_interp
2   !!======================================================================
3   !!                   ***  MODULE  agrif_oce_interp  ***
4   !! AGRIF: interpolation package for the ocean dynamics (OPA)
5   !!======================================================================
6   !! History :  2.0  !  2002-06  (L. Debreu)  Original cade
7   !!            3.2  !  2009-04  (R. Benshila)
8   !!            3.6  !  2014-09  (R. Benshila)
9   !!----------------------------------------------------------------------
10#if defined key_agrif
11   !!----------------------------------------------------------------------
12   !!   'key_agrif'                                              AGRIF zoom
13   !!----------------------------------------------------------------------
14   !!   Agrif_tra     :
15   !!   Agrif_dyn     :
16   !!   Agrif_ssh     :
17   !!   Agrif_dyn_ts  :
18   !!   Agrif_dta_ts  :
19   !!   Agrif_ssh_ts  :
20   !!   Agrif_avm     :
21   !!   interpu       :
22   !!   interpv       :
23   !!----------------------------------------------------------------------
24   USE par_oce
25   USE oce
26   USE dom_oce     
27   USE zdf_oce
28   USE agrif_oce
29   USE phycst
30   USE dynspg_ts, ONLY: un_adv, vn_adv
31   !
32   USE in_out_manager
33   USE agrif_oce_sponge
34   USE lib_mpp
35   USE vremap
36   USE lbclnk
37 
38   IMPLICIT NONE
39   PRIVATE
40
41   PUBLIC   Agrif_dyn, Agrif_ssh, Agrif_dyn_ts, Agrif_dyn_ts_flux, Agrif_ssh_ts, Agrif_dta_ts
42   PUBLIC   Agrif_tra, Agrif_avm
43   PUBLIC   interpun , interpvn
44   PUBLIC   interptsn, interpsshn, interpavm
45   PUBLIC   interpunb, interpvnb , interpub2b, interpvb2b
46   PUBLIC   interpe3t, interpglamt, interpgphit
47   PUBLIC   interpht0, interpmbkt, interpe3t0_vremap
48   PUBLIC   agrif_istate_oce, agrif_istate_ssh   ! called by icestate.F90 and domvvl.F90
49   PUBLIC   agrif_check_bat
50
51   INTEGER ::   bdy_tinterp = 0
52
53   !!----------------------------------------------------------------------
54   !! NEMO/NST 4.0 , NEMO Consortium (2018)
55   !! $Id$
56   !! Software governed by the CeCILL license (see ./LICENSE)
57   !!----------------------------------------------------------------------
58CONTAINS
59
60   SUBROUTINE Agrif_istate_oce( Kbb, Kmm, Kaa )
61      !!----------------------------------------------------------------------
62      !!                 *** ROUTINE agrif_istate_oce ***
63      !!
64      !!                 set initial t, s, u, v, ssh from parent
65      !!----------------------------------------------------------------------
66      !
67      IMPLICIT NONE
68      !
69      INTEGER, INTENT(in)  :: Kbb, Kmm, Kaa
70      INTEGER :: jn
71      !!----------------------------------------------------------------------
72      IF(lwp) WRITE(numout,*) ' '
73      IF(lwp) WRITE(numout,*) 'Agrif_istate_oce : interp child initial state from parent'
74      IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~~~~'
75      IF(lwp) WRITE(numout,*) ' '
76
77      IF ( ln_rstart ) & 
78         & CALL ctl_stop('AGRIF hot start should be desactivated in restarting mode')
79
80      IF ( .NOT.Agrif_Parent(l_1st_euler) ) & 
81         & CALL ctl_stop('AGRIF hot start requires to force Euler first step on parent')
82
83      l_ini_child           = .TRUE.
84      Agrif_SpecialValue    = 0.0_wp
85      Agrif_UseSpecialValue = .TRUE.
86
87      ts(:,:,:,:,:) = 0.0_wp
88      uu(:,:,:,:)   = 0.0_wp
89      vv(:,:,:,:)   = 0.0_wp 
90      ssh(:,:,:)    = 0._wp
91       
92      Krhs_a = Kbb   ;   Kmm_a = Kbb
93
94      CALL Agrif_Init_Variable(tsini_id, procname=interptsn)
95      CALL Agrif_Init_Variable(sshini_id, procname=interpsshn)
96
97      Agrif_UseSpecialValue = ln_spc_dyn
98      use_sign_north = .TRUE.
99      sign_north = -1._wp
100      CALL Agrif_Init_Variable(uini_id , procname=interpun )
101      CALL Agrif_Init_Variable(vini_id , procname=interpvn )
102      use_sign_north = .FALSE.
103
104      Agrif_UseSpecialValue = .FALSE.
105      l_ini_child           = .FALSE.
106
107      Krhs_a = Kaa   ;   Kmm_a = Kmm
108
109      ssh(:,:,Kbb) = ssh(:,:,Kbb) * tmask(:,:,1)
110
111      DO jn = 1, jpts
112         ts(:,:,:,jn,Kbb) = ts(:,:,:,jn,Kbb) * tmask(:,:,:)
113      END DO
114      uu(:,:,:,Kbb) = uu(:,:,:,Kbb) * umask(:,:,:)     
115      vv(:,:,:,Kbb) = vv(:,:,:,Kbb) * vmask(:,:,:) 
116
117      CALL lbc_lnk_multi( 'agrif_istate_oce', uu(:,:,:  ,Kbb), 'U', -1.0_wp , vv(:,:,:,Kbb), 'V', -1.0_wp )
118      CALL lbc_lnk( 'agrif_istate_oce', ts(:,:,:,:,Kbb), 'T',  1.0_wp )
119      CALL lbc_lnk( 'agrif_istate_oce', ssh(:,:,Kbb), 'T', 1.0_wp )
120
121   END SUBROUTINE Agrif_istate_oce
122
123
124   SUBROUTINE Agrif_istate_ssh( Kbb, Kmm )
125      !!----------------------------------------------------------------------
126      !!                 *** ROUTINE agrif_istate_ssh ***
127      !!
128      !!                    set initial ssh from parent
129      !!----------------------------------------------------------------------
130      !
131      IMPLICIT NONE
132      !
133      INTEGER, INTENT(in)  :: Kbb, Kmm 
134      !!----------------------------------------------------------------------
135      IF(lwp) WRITE(numout,*) ' '
136      IF(lwp) WRITE(numout,*) 'Agrif_istate_ssh : interp child ssh from parent'
137      IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~~~~'
138      IF(lwp) WRITE(numout,*) ' '
139
140      IF ( ln_rstart ) & 
141         & CALL ctl_stop('AGRIF hot start should be desactivated in restarting mode')
142
143      IF ( .NOT.Agrif_Parent(l_1st_euler) ) & 
144         & CALL ctl_stop('AGRIF hot start requires to force Euler first step on parent')
145
146      Kmm_a = Kmm
147      ssh(:,:,Kmm) = 0._wp
148
149      Agrif_SpecialValue    = 0._wp
150      Agrif_UseSpecialValue = .TRUE.
151      l_ini_child           = .TRUE.
152      !
153      CALL Agrif_Init_Variable(sshini_id, procname=interpsshn)
154      !
155      Agrif_UseSpecialValue = .FALSE.
156      l_ini_child           = .FALSE.
157      CALL lbc_lnk( 'dom_vvl_rst', ssh(:,:,Kmm), 'T', 1._wp )
158
159   END SUBROUTINE Agrif_istate_ssh
160
161
162   SUBROUTINE Agrif_tra
163      !!----------------------------------------------------------------------
164      !!                  ***  ROUTINE Agrif_tra  ***
165      !!----------------------------------------------------------------------
166      !
167      IF( Agrif_Root() )   RETURN
168      !
169      Agrif_SpecialValue    = 0._wp
170      Agrif_UseSpecialValue = .TRUE.
171      l_vremap           = ln_vert_remap
172      !
173      CALL Agrif_Bc_variable( ts_interp_id, procname=interptsn )
174      !
175      Agrif_UseSpecialValue = .FALSE.
176      l_vremap              = .FALSE.
177      !
178   END SUBROUTINE Agrif_tra
179
180
181   SUBROUTINE Agrif_dyn( kt )
182      !!----------------------------------------------------------------------
183      !!                  ***  ROUTINE Agrif_DYN  ***
184      !!---------------------------------------------------------------------- 
185      INTEGER, INTENT(in) ::   kt
186      !
187      INTEGER ::   ji, jj, jk       ! dummy loop indices
188      INTEGER ::   ibdy1, jbdy1, ibdy2, jbdy2
189      REAL(wp), DIMENSION(jpi,jpj) ::   zub, zvb
190      !!---------------------------------------------------------------------- 
191      !
192      IF( Agrif_Root() )   RETURN
193      !
194      Agrif_SpecialValue    = 0.0_wp
195      Agrif_UseSpecialValue = ln_spc_dyn
196      l_vremap              = ln_vert_remap
197      !
198      use_sign_north = .TRUE.
199      sign_north = -1.0_wp
200      CALL Agrif_Bc_variable( un_interp_id, procname=interpun )
201      CALL Agrif_Bc_variable( vn_interp_id, procname=interpvn )
202
203      IF( .NOT.ln_dynspg_ts ) THEN ! Get transports
204         ubdy(:,:) = 0._wp   ;   vbdy(:,:) = 0._wp 
205         CALL Agrif_Bc_variable( unb_interp_id, procname=interpunb )
206         CALL Agrif_Bc_variable( vnb_interp_id, procname=interpvnb )
207      ENDIF
208
209      use_sign_north = .FALSE.
210      !
211      Agrif_UseSpecialValue = .FALSE.
212      l_vremap              = .FALSE.
213      !
214      ! Ensure below that vertically integrated transports match
215      ! either transports out of time splitting procedure (ln_dynspg_ts=.TRUE.)
216      ! or parent grid transports (ln_dynspg_ts=.FALSE.)
217      !
218      ! --- West --- !
219      IF( lk_west ) THEN
220         ibdy1 = nn_hls + 2                  ! halo + land + 1
221         ibdy2 = nn_hls + 1 + nbghostcells + nn_shift_bar*Agrif_Rhox()   ! halo + land + nbghostcells
222         !
223         IF( .NOT.ln_dynspg_ts ) THEN  ! Store transport
224            DO ji = mi0(ibdy1), mi1(ibdy2)
225               DO jj = 1, jpj
226                  uu_b(ji,jj,Krhs_a) = ubdy(ji,jj) * r1_hu(ji,jj,Krhs_a)
227                  vv_b(ji,jj,Krhs_a) = vbdy(ji,jj) * r1_hv(ji,jj,Krhs_a)
228               END DO
229            END DO
230         ENDIF
231         !
232         DO ji = mi0(ibdy1), mi1(ibdy2)
233            zub(ji,:) = 0._wp 
234            DO jk = 1, jpkm1
235               DO jj = 1, jpj
236                  zub(ji,jj) = zub(ji,jj) + e3u(ji,jj,jk,Krhs_a)  * uu(ji,jj,jk,Krhs_a) * umask(ji,jj,jk)
237               END DO
238            END DO
239            DO jj=1,jpj
240               zub(ji,jj) = zub(ji,jj) * r1_hu(ji,jj,Krhs_a)
241            END DO
242            DO jk = 1, jpkm1
243               DO jj = 1, jpj
244                  uu(ji,jj,jk,Krhs_a) = ( uu(ji,jj,jk,Krhs_a) + uu_b(ji,jj,Krhs_a) - zub(ji,jj) ) * umask(ji,jj,jk)
245               END DO
246            END DO
247         END DO
248         !   
249         DO ji = mi0(ibdy1), mi1(ibdy2)
250            zvb(ji,:) = 0._wp
251            DO jk = 1, jpkm1
252               DO jj = 1, jpj
253                  zvb(ji,jj) = zvb(ji,jj) + e3v(ji,jj,jk,Krhs_a) * vv(ji,jj,jk,Krhs_a) * vmask(ji,jj,jk)
254               END DO
255            END DO
256            DO jj = 1, jpj
257               zvb(ji,jj) = zvb(ji,jj) * r1_hv(ji,jj,Krhs_a)
258            END DO
259            DO jk = 1, jpkm1
260               DO jj = 1, jpj
261                  vv(ji,jj,jk,Krhs_a) = ( vv(ji,jj,jk,Krhs_a) + vv_b(ji,jj,Krhs_a) - zvb(ji,jj) )*vmask(ji,jj,jk)
262               END DO
263            END DO
264         END DO
265         !
266      ENDIF
267
268      ! --- East --- !
269      IF( lk_east) THEN
270         ibdy1 = jpiglo - ( nn_hls + nbghostcells + 1) - nn_shift_bar*Agrif_Rhox()   
271         ibdy2 = jpiglo - ( nn_hls + 2 )                 
272         !
273         IF( .NOT.ln_dynspg_ts ) THEN
274            DO ji = mi0(ibdy1), mi1(ibdy2)
275               uu_b(ji,:,Krhs_a) = 0._wp
276               DO jk = 1, jpkm1
277                  DO jj = 1, jpj
278                     uu_b(ji,jj,Krhs_a) = uu_b(ji,jj,Krhs_a) + e3u(ji,jj,jk,Krhs_a) * uu(ji,jj,jk,Krhs_a) * umask(ji,jj,jk)
279                  END DO
280               END DO
281               DO jj = 1, jpj
282                  uu_b(ji,jj,Krhs_a) = ubdy(ji,jj) * r1_hu(ji,jj,Krhs_a)
283               END DO
284            END DO
285         ENDIF
286         !
287         DO ji = mi0(ibdy1), mi1(ibdy2)
288            zub(ji,:) = 0._wp   
289            DO jk = 1, jpkm1
290               DO jj = 1, jpj
291                  zub(ji,jj) = zub(ji,jj) + e3u(ji,jj,jk,Krhs_a)  * uu(ji,jj,jk,Krhs_a) * umask(ji,jj,jk)
292               END DO
293            END DO
294            DO jj=1,jpj
295               zub(ji,jj) = zub(ji,jj) * r1_hu(ji,jj,Krhs_a)
296            END DO
297            DO jk = 1, jpkm1
298               DO jj = 1, jpj
299                  uu(ji,jj,jk,Krhs_a) = ( uu(ji,jj,jk,Krhs_a) + uu_b(ji,jj,Krhs_a) - zub(ji,jj) ) * umask(ji,jj,jk)
300               END DO
301            END DO
302         END DO
303         !
304         ibdy1 = jpiglo - ( nn_hls + nbghostcells ) - nn_shift_bar*Agrif_Rhox() 
305         ibdy2 = jpiglo - ( nn_hls + 1 )              !
306         IF( .NOT.ln_dynspg_ts ) THEN
307            DO ji = mi0(ibdy1), mi1(ibdy2)
308               vv_b(ji,:,Krhs_a) = 0._wp
309               DO jk = 1, jpkm1
310                  DO jj = 1, jpj
311                     vv_b(ji,jj,Krhs_a) = vv_b(ji,jj,Krhs_a) + e3v(ji,jj,jk,Krhs_a) * vv(ji,jj,jk,Krhs_a) * vmask(ji,jj,jk)
312                  END DO
313               END DO
314               DO jj = 1, jpj
315                  vv_b(ji,jj,Krhs_a) = vbdy(ji,jj) * r1_hv(ji,jj,Krhs_a)
316               END DO
317            END DO
318         ENDIF
319
320         DO ji = mi0(ibdy1), mi1(ibdy2)
321            zvb(ji,:) = 0._wp
322            DO jk = 1, jpkm1
323               DO jj = 1, jpj
324                     zvb(ji,jj) = zvb(ji,jj) + e3v(ji,jj,jk,Krhs_a) * vv(ji,jj,jk,Krhs_a) * vmask(ji,jj,jk)
325               END DO
326            END DO
327            DO jj = 1, jpj
328               zvb(ji,jj) = zvb(ji,jj) * r1_hv(ji,jj,Krhs_a)
329            END DO
330            DO jk = 1, jpkm1
331               DO jj = 1, jpj
332                     vv(ji,jj,jk,Krhs_a) = ( vv(ji,jj,jk,Krhs_a) + vv_b(ji,jj,Krhs_a) - zvb(ji,jj) ) * vmask(ji,jj,jk)
333               END DO
334            END DO
335         END DO
336         !
337      ENDIF
338
339      ! --- South --- !
340      IF( lk_south ) THEN
341         jbdy1 = nn_hls + 2                 
342         jbdy2 = nn_hls + 1 + nbghostcells + nn_shift_bar*Agrif_Rhoy()   
343         !
344         IF( .NOT.ln_dynspg_ts ) THEN
345            DO jj = mj0(jbdy1), mj1(jbdy2)
346               vv_b(:,jj,Krhs_a) = 0._wp
347               DO jk = 1, jpkm1
348                  DO ji = 1, jpi
349                     vv_b(ji,jj,Krhs_a) = vv_b(ji,jj,Krhs_a) + e3v(ji,jj,jk,Krhs_a) * vv(ji,jj,jk,Krhs_a) * vmask(ji,jj,jk)
350                  END DO
351               END DO
352               DO ji=1,jpi
353                  vv_b(ji,jj,Krhs_a) = vv_b(ji,jj,Krhs_a) * r1_hv(ji,jj,Krhs_a) 
354                  uu_b(ji,jj,Krhs_a) = uu_b(ji,jj,Krhs_a) * r1_hu(ji,jj,Krhs_a)
355               END DO
356            END DO
357         ENDIF
358         !
359         DO jj = mj0(jbdy1), mj1(jbdy2)
360            zvb(:,jj) = 0._wp
361            DO jk=1,jpkm1
362               DO ji=1,jpi
363                  zvb(ji,jj) = zvb(ji,jj) + e3v(ji,jj,jk,Krhs_a) * vv(ji,jj,jk,Krhs_a) * vmask(ji,jj,jk)
364               END DO
365            END DO
366            DO ji = 1, jpi
367               zvb(ji,jj) = zvb(ji,jj) * r1_hv(ji,jj,Krhs_a)
368            END DO
369            DO jk = 1, jpkm1
370               DO ji = 1, jpi
371                  vv(ji,jj,jk,Krhs_a) = ( vv(ji,jj,jk,Krhs_a) + vv_b(ji,jj,Krhs_a) - zvb(ji,jj) ) * vmask(ji,jj,jk)
372               END DO
373            END DO
374         END DO
375         !
376         DO jj = mj0(jbdy1), mj1(jbdy2)
377            zub(:,jj) = 0._wp
378            DO jk = 1, jpkm1
379               DO ji = 1, jpi
380                  zub(ji,jj) = zub(ji,jj) + e3u(ji,jj,jk,Krhs_a) * uu(ji,jj,jk,Krhs_a) * umask(ji,jj,jk)
381               END DO
382            END DO
383            DO ji = 1, jpi
384               zub(ji,jj) = zub(ji,jj) * r1_hu(ji,jj,Krhs_a)
385            END DO
386            DO jk = 1, jpkm1
387               DO ji = 1, jpi
388                  uu(ji,jj,jk,Krhs_a) = ( uu(ji,jj,jk,Krhs_a) + uu_b(ji,jj,Krhs_a) - zub(ji,jj) ) * umask(ji,jj,jk)
389               END DO
390            END DO
391         END DO
392         !
393      ENDIF
394
395      ! --- North --- !
396      IF( lk_north ) THEN
397         jbdy1 = jpjglo - ( nn_hls + nbghostcells + 1) - nn_shift_bar*Agrif_Rhoy() 
398         jbdy2 = jpjglo - ( nn_hls + 2 )
399         !
400         IF( .NOT.ln_dynspg_ts ) THEN
401            DO jj = mj0(jbdy1), mj1(jbdy2)
402               vv_b(:,jj,Krhs_a) = 0._wp
403               DO jk = 1, jpkm1
404                  DO ji = 1, jpi
405                     vv_b(ji,jj,Krhs_a) = vv_b(ji,jj,Krhs_a) + e3v(ji,jj,jk,Krhs_a) * vv(ji,jj,jk,Krhs_a) * vmask(ji,jj,jk)
406                  END DO
407               END DO
408               DO ji=1,jpi
409                  vv_b(ji,jj,Krhs_a) = vbdy(ji,jj) * r1_hv(ji,jj,Krhs_a)
410               END DO
411            END DO
412         ENDIF
413         !
414         DO jj = mj0(jbdy1), mj1(jbdy2)
415            zvb(:,jj) = 0._wp 
416            DO jk=1,jpkm1
417               DO ji=1,jpi
418                  zvb(ji,jj) = zvb(ji,jj) + e3v(ji,jj,jk,Krhs_a) * vv(ji,jj,jk,Krhs_a) * vmask(ji,jj,jk)
419               END DO
420            END DO
421            DO ji = 1, jpi
422               zvb(ji,jj) = zvb(ji,jj) * r1_hv(ji,jj,Krhs_a)
423            END DO
424            DO jk = 1, jpkm1
425               DO ji = 1, jpi
426                  vv(ji,jj,jk,Krhs_a) = ( vv(ji,jj,jk,Krhs_a) + vv_b(ji,jj,Krhs_a) - zvb(ji,jj) ) * vmask(ji,jj,jk)
427               END DO
428            END DO
429         END DO
430         !
431         jbdy1 = jpjglo - ( nn_hls + nbghostcells ) - nn_shift_bar*Agrif_Rhoy() 
432         jbdy2 = jpjglo - ( nn_hls + 1 )
433         IF( .NOT.ln_dynspg_ts ) THEN
434            DO jj = mj0(jbdy1), mj1(jbdy2)
435               uu_b(:,jj,Krhs_a) = 0._wp
436               DO jk = 1, jpkm1
437                  DO ji = 1, jpi
438                     uu_b(ji,jj,Krhs_a) = uu_b(ji,jj,Krhs_a) + e3u(ji,jj,jk,Krhs_a) * uu(ji,jj,jk,Krhs_a) * umask(ji,jj,jk)
439                  END DO
440               END DO
441               DO ji=1,jpi
442                  uu_b(ji,jj,Krhs_a) = ubdy(ji,jj) * r1_hu(ji,jj,Krhs_a)
443               END DO
444            END DO
445         ENDIF
446
447         DO jj = mj0(jbdy1), mj1(jbdy2)
448            zub(:,jj) = 0._wp
449            DO jk = 1, jpkm1
450               DO ji = 1, jpi
451                  zub(ji,jj) = zub(ji,jj) + e3u(ji,jj,jk,Krhs_a) * uu(ji,jj,jk,Krhs_a) * umask(ji,jj,jk)
452               END DO
453            END DO
454            DO ji = 1, jpi
455               zub(ji,jj) = zub(ji,jj) * r1_hu(ji,jj,Krhs_a)
456            END DO
457            DO jk = 1, jpkm1
458               DO ji = 1, jpi
459                     uu(ji,jj,jk,Krhs_a) = ( uu(ji,jj,jk,Krhs_a) + uu_b(ji,jj,Krhs_a) - zub(ji,jj) ) * umask(ji,jj,jk)
460               END DO
461            END DO
462         END DO
463         !
464      ENDIF
465      !
466   END SUBROUTINE Agrif_dyn
467
468
469   SUBROUTINE Agrif_dyn_ts( jn )
470      !!----------------------------------------------------------------------
471      !!                  ***  ROUTINE Agrif_dyn_ts  ***
472      !!---------------------------------------------------------------------- 
473      INTEGER, INTENT(in) ::   jn
474      !!
475      INTEGER :: ji, jj
476      INTEGER :: istart, iend, jstart, jend
477      !!---------------------------------------------------------------------- 
478      !
479      IF( Agrif_Root() )   RETURN
480      !
481      !--- West ---!
482      IF( lk_west ) THEN
483         istart = nn_hls + 2                              ! halo + land + 1
484         iend   = nn_hls + 1 + nbghostcells  + nn_shift_bar*Agrif_Rhox()              ! halo + land + nbghostcells
485         DO ji = mi0(istart), mi1(iend)
486            DO jj=1,jpj
487               va_e(ji,jj) = vbdy(ji,jj) * hvr_e(ji,jj)
488               ua_e(ji,jj) = ubdy(ji,jj) * hur_e(ji,jj)
489            END DO
490         END DO
491      ENDIF
492      !
493      !--- East ---!
494      IF( lk_east ) THEN
495         istart = jpiglo - ( nn_hls + nbghostcells ) - nn_shift_bar*Agrif_Rhox() 
496         iend   = jpiglo - ( nn_hls + 1 )               
497         DO ji = mi0(istart), mi1(iend)
498
499            DO jj=1,jpj
500               va_e(ji,jj) = vbdy(ji,jj) * hvr_e(ji,jj)
501            END DO
502         END DO
503         istart = jpiglo - ( nn_hls + nbghostcells + 1) - nn_shift_bar*Agrif_Rhox() 
504         iend   = jpiglo - ( nn_hls + 2 )               
505         DO ji = mi0(istart), mi1(iend)
506            DO jj=1,jpj
507               ua_e(ji,jj) = ubdy(ji,jj) * hur_e(ji,jj)
508            END DO
509         END DO
510      ENDIF 
511      !
512      !--- South ---!
513      IF( lk_south ) THEN
514         jstart = nn_hls + 2                             
515         jend   = nn_hls + 1 + nbghostcells + nn_shift_bar*Agrif_Rhoy()           
516         DO jj = mj0(jstart), mj1(jend)
517
518            DO ji=1,jpi
519               ua_e(ji,jj) = ubdy(ji,jj) * hur_e(ji,jj)
520               va_e(ji,jj) = vbdy(ji,jj) * hvr_e(ji,jj)
521            END DO
522         END DO
523      ENDIF       
524      !
525      !--- North ---!
526      IF( lk_north ) THEN
527         jstart = jpjglo - ( nn_hls + nbghostcells ) - nn_shift_bar*Agrif_Rhoy()     
528         jend   = jpjglo - ( nn_hls + 1 )               
529         DO jj = mj0(jstart), mj1(jend)
530            DO ji=1,jpi
531               ua_e(ji,jj) = ubdy(ji,jj) * hur_e(ji,jj)
532            END DO
533         END DO
534         jstart = jpjglo - ( nn_hls + nbghostcells + 1) - nn_shift_bar*Agrif_Rhoy() 
535         jend   = jpjglo - ( nn_hls + 2 )               
536         DO jj = mj0(jstart), mj1(jend)
537            DO ji=1,jpi
538               va_e(ji,jj) = vbdy(ji,jj) * hvr_e(ji,jj)
539            END DO
540         END DO
541      ENDIF 
542      !
543   END SUBROUTINE Agrif_dyn_ts
544
545   
546   SUBROUTINE Agrif_dyn_ts_flux( jn, zu, zv )
547      !!----------------------------------------------------------------------
548      !!                  ***  ROUTINE Agrif_dyn_ts_flux  ***
549      !!---------------------------------------------------------------------- 
550      INTEGER, INTENT(in) ::   jn
551      REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) ::   zu, zv
552      !!
553      INTEGER :: ji, jj
554      INTEGER :: istart, iend, jstart, jend
555      !!---------------------------------------------------------------------- 
556      !
557      IF( Agrif_Root() )   RETURN
558      !
559      !--- West ---!
560      IF( lk_west ) THEN
561         istart = nn_hls + 2                             
562         iend   = nn_hls + 1 + nbghostcells + nn_shift_bar*Agrif_Rhox() 
563         DO ji = mi0(istart), mi1(iend)
564            DO jj=1,jpj
565               zv(ji,jj) = vbdy(ji,jj) * e1v(ji,jj)
566               zu(ji,jj) = ubdy(ji,jj) * e2u(ji,jj)
567            END DO
568         END DO
569      ENDIF
570      !
571      !--- East ---!
572      IF( lk_east ) THEN
573         istart = jpiglo - ( nn_hls + nbghostcells ) - nn_shift_bar*Agrif_Rhox()
574         iend   = jpiglo - ( nn_hls + 1 )                 
575         DO ji = mi0(istart), mi1(iend)
576            DO jj=1,jpj
577               zv(ji,jj) = vbdy(ji,jj) * e1v(ji,jj)
578            END DO
579         END DO
580         istart = jpiglo - ( nn_hls + nbghostcells + 1) - nn_shift_bar*Agrif_Rhox() 
581         iend   = jpiglo - ( nn_hls + 2 )                 
582         DO ji = mi0(istart), mi1(iend)
583            DO jj=1,jpj
584               zu(ji,jj) = ubdy(ji,jj) * e2u(ji,jj)
585            END DO
586         END DO
587      ENDIF
588      !
589      !--- South ---!
590      IF( lk_south ) THEN
591         jstart = nn_hls + 2                             
592         jend   = nn_hls + 1 + nbghostcells + nn_shift_bar*Agrif_Rhoy() 
593         DO jj = mj0(jstart), mj1(jend)
594            DO ji=1,jpi
595               zu(ji,jj) = ubdy(ji,jj) * e2u(ji,jj)
596               zv(ji,jj) = vbdy(ji,jj) * e1v(ji,jj)
597            END DO
598         END DO
599      ENDIF
600      !
601      !--- North ---!
602      IF( lk_north ) THEN
603         jstart = jpjglo - ( nn_hls + nbghostcells ) - nn_shift_bar*Agrif_Rhoy() 
604         jend   = jpjglo - ( nn_hls + 1 )               
605         DO jj = mj0(jstart), mj1(jend)
606            DO ji=1,jpi
607               zu(ji,jj) = ubdy(ji,jj) * e2u(ji,jj)
608            END DO
609         END DO
610         jstart = jpjglo - ( nn_hls + nbghostcells + 1) - nn_shift_bar*Agrif_Rhoy() 
611         jend   = jpjglo - ( nn_hls + 2 )               
612         DO jj = mj0(jstart), mj1(jend)
613            DO ji=1,jpi
614               zv(ji,jj) = vbdy(ji,jj) * e1v(ji,jj)
615            END DO
616         END DO
617      ENDIF
618      !
619   END SUBROUTINE Agrif_dyn_ts_flux
620
621   
622   SUBROUTINE Agrif_dta_ts( kt )
623      !!----------------------------------------------------------------------
624      !!                  ***  ROUTINE Agrif_dta_ts  ***
625      !!---------------------------------------------------------------------- 
626      INTEGER, INTENT(in) ::   kt
627      !!
628      LOGICAL :: ll_int_cons
629      !!---------------------------------------------------------------------- 
630      !
631      IF( Agrif_Root() )   RETURN
632      !
633      ll_int_cons = ln_bt_fw ! Assume conservative temporal integration in the forward case only
634      !
635      ! Enforce volume conservation if no time refinement: 
636      IF ( Agrif_rhot()==1 ) ll_int_cons=.TRUE. 
637      !
638      ! Interpolate barotropic fluxes
639      Agrif_SpecialValue = 0._wp
640      Agrif_UseSpecialValue = ln_spc_dyn
641
642      use_sign_north = .TRUE.
643      sign_north = -1.
644
645      !
646      ! Set bdy time interpolation stage to 0 (latter incremented locally do deal with corners)
647      utint_stage(:,:) = 0
648      vtint_stage(:,:) = 0
649      !
650      IF( ll_int_cons ) THEN  ! Conservative interpolation
651         IF ( lk_tint2d_notinterp ) THEN
652            CALL Agrif_Bc_variable( ub2b_interp_id, calledweight=1._wp, procname=interpub2b_const )
653            CALL Agrif_Bc_variable( vb2b_interp_id, calledweight=1._wp, procname=interpvb2b_const ) 
654            ! Divergence conserving correction terms:
655            CALL Agrif_Bc_variable(    ub2b_cor_id, calledweight=1._wp, procname=        ub2b_cor )
656            CALL Agrif_Bc_variable(    vb2b_cor_id, calledweight=1._wp, procname=        vb2b_cor )
657         ELSE
658            ! order matters here !!!!!!
659            CALL Agrif_Bc_variable( ub2b_interp_id, calledweight=1._wp, procname=interpub2b ) ! Time integrated
660            CALL Agrif_Bc_variable( vb2b_interp_id, calledweight=1._wp, procname=interpvb2b ) 
661            !
662            bdy_tinterp = 1
663            CALL Agrif_Bc_variable( unb_interp_id , calledweight=1._wp, procname=interpunb  ) ! After
664            CALL Agrif_Bc_variable( vnb_interp_id , calledweight=1._wp, procname=interpvnb  ) 
665            !
666            bdy_tinterp = 2
667            CALL Agrif_Bc_variable( unb_interp_id , calledweight=0._wp, procname=interpunb  ) ! Before
668            CALL Agrif_Bc_variable( vnb_interp_id , calledweight=0._wp, procname=interpvnb  )   
669         ENDIF
670      ELSE ! Linear interpolation
671         !
672         ubdy(:,:) = 0._wp   ;   vbdy(:,:) = 0._wp 
673         CALL Agrif_Bc_variable( unb_interp_id, procname=interpunb )
674         CALL Agrif_Bc_variable( vnb_interp_id, procname=interpvnb )
675      ENDIF
676      Agrif_UseSpecialValue = .FALSE.
677      use_sign_north = .FALSE.
678      !
679   END SUBROUTINE Agrif_dta_ts
680
681
682   SUBROUTINE Agrif_ssh( kt )
683      !!----------------------------------------------------------------------
684      !!                  ***  ROUTINE Agrif_ssh  ***
685      !!---------------------------------------------------------------------- 
686      INTEGER, INTENT(in) ::   kt
687      !
688      INTEGER  :: ji, jj
689      INTEGER  :: istart, iend, jstart, jend
690      !!---------------------------------------------------------------------- 
691      !
692      IF( Agrif_Root() )   RETURN
693      !     
694      ! Linear time interpolation of sea level
695      !
696      Agrif_SpecialValue    = 0._wp
697      Agrif_UseSpecialValue = .TRUE.
698      CALL Agrif_Bc_variable(sshn_id, procname=interpsshn )
699      Agrif_UseSpecialValue = .FALSE.
700      !
701      ! --- West --- !
702      IF(lk_west) THEN
703         istart = nn_hls + 2                                                          ! halo + land + 1
704         iend   = nn_hls + 1 + nbghostcells + nn_shift_bar*Agrif_Rhox()               ! halo + land + nbghostcells
705         DO ji = mi0(istart), mi1(iend)
706            DO jj = 1, jpj
707               ssh(ji,jj,Krhs_a) = hbdy(ji,jj)
708            END DO
709         END DO
710      ENDIF
711      !
712      ! --- East --- !
713      IF(lk_east) THEN
714         istart = jpiglo - ( nn_hls + nbghostcells ) - nn_shift_bar*Agrif_Rhox()       ! halo + land + nbghostcells - 1
715         iend   = jpiglo - ( nn_hls + 1 )                                              ! halo + land + 1            - 1
716         DO ji = mi0(istart), mi1(iend)
717            DO jj = 1, jpj
718               ssh(ji,jj,Krhs_a) = hbdy(ji,jj)
719            END DO
720         END DO
721      ENDIF
722      !
723      ! --- South --- !
724      IF(lk_south) THEN
725         jstart = nn_hls + 2                                                          ! halo + land + 1
726         jend   = nn_hls + 1 + nbghostcells + nn_shift_bar*Agrif_Rhoy()               ! halo + land + nbghostcells
727         DO jj = mj0(jstart), mj1(jend)
728            DO ji = 1, jpi
729               ssh(ji,jj,Krhs_a) = hbdy(ji,jj)
730            END DO
731         END DO
732      ENDIF
733      !
734      ! --- North --- !
735      IF(lk_north) THEN
736         jstart = jpjglo - ( nn_hls + nbghostcells ) - nn_shift_bar*Agrif_Rhoy()     ! halo + land + nbghostcells - 1
737         jend   = jpjglo - ( nn_hls + 1 )                                            ! halo + land + 1            - 1
738         DO jj = mj0(jstart), mj1(jend)
739            DO ji = 1, jpi
740               ssh(ji,jj,Krhs_a) = hbdy(ji,jj)
741            END DO
742         END DO
743      ENDIF
744      !
745   END SUBROUTINE Agrif_ssh
746
747
748   SUBROUTINE Agrif_ssh_ts( jn )
749      !!----------------------------------------------------------------------
750      !!                  ***  ROUTINE Agrif_ssh_ts  ***
751      !!---------------------------------------------------------------------- 
752      INTEGER, INTENT(in) ::   jn
753      !!
754      INTEGER :: ji, jj
755      INTEGER  :: istart, iend, jstart, jend
756      !!---------------------------------------------------------------------- 
757      !
758      IF( Agrif_Root() )   RETURN
759      !
760      ! --- West --- !
761      IF(lk_west) THEN
762         istart = nn_hls + 2                                                        ! halo + land + 1
763         iend   = nn_hls + 1 + nbghostcells + nn_shift_bar*Agrif_Rhox()             ! halo + land + nbghostcells
764         DO ji = mi0(istart), mi1(iend)
765            DO jj = 1, jpj
766               ssha_e(ji,jj) = hbdy(ji,jj)
767            END DO
768         END DO
769      ENDIF
770      !
771      ! --- East --- !
772      IF(lk_east) THEN
773         istart = jpiglo - ( nn_hls + nbghostcells ) - nn_shift_bar*Agrif_Rhox()    ! halo + land + nbghostcells - 1
774         iend   = jpiglo - ( nn_hls + 1 )                                           ! halo + land + 1            - 1
775         DO ji = mi0(istart), mi1(iend)
776            DO jj = 1, jpj
777               ssha_e(ji,jj) = hbdy(ji,jj)
778            END DO
779         END DO
780      ENDIF
781      !
782      ! --- South --- !
783      IF(lk_south) THEN
784         jstart = nn_hls + 2                                                        ! halo + land + 1
785         jend   = nn_hls + 1 + nbghostcells + nn_shift_bar*Agrif_Rhoy()             ! halo + land + nbghostcells
786         DO jj = mj0(jstart), mj1(jend)
787            DO ji = 1, jpi
788               ssha_e(ji,jj) = hbdy(ji,jj)
789            END DO
790         END DO
791      ENDIF
792      !
793      ! --- North --- !
794      IF(lk_north) THEN
795         jstart = jpjglo - ( nn_hls + nbghostcells ) - nn_shift_bar*Agrif_Rhoy()    ! halo + land + nbghostcells - 1
796         jend   = jpjglo - ( nn_hls + 1 )                                           ! halo + land + 1            - 1
797         DO jj = mj0(jstart), mj1(jend)
798            DO ji = 1, jpi
799               ssha_e(ji,jj) = hbdy(ji,jj)
800            END DO
801         END DO
802      ENDIF
803      !
804   END SUBROUTINE Agrif_ssh_ts
805
806   
807   SUBROUTINE Agrif_avm
808      !!----------------------------------------------------------------------
809      !!                  ***  ROUTINE Agrif_avm  ***
810      !!---------------------------------------------------------------------- 
811      REAL(wp) ::   zalpha
812      !!---------------------------------------------------------------------- 
813      !
814      IF( Agrif_Root() )   RETURN
815      !
816      zalpha = 1._wp ! JC: proper time interpolation impossible 
817                     ! => use last available value from parent
818      !
819      Agrif_SpecialValue    = 0.e0
820      Agrif_UseSpecialValue = .TRUE.
821      l_vremap              = ln_vert_remap
822      !
823      CALL Agrif_Bc_variable( avm_id, calledweight=zalpha, procname=interpavm )       
824      !
825      Agrif_UseSpecialValue = .FALSE.
826      l_vremap              = .FALSE.
827      !
828   END SUBROUTINE Agrif_avm
829
830
831   SUBROUTINE interptsn( ptab, i1, i2, j1, j2, k1, k2, n1, n2, before )
832      !!----------------------------------------------------------------------
833      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) ::   ptab
834      INTEGER                                     , INTENT(in   ) ::   i1, i2, j1, j2, k1, k2, n1, n2
835      LOGICAL                                     , INTENT(in   ) ::   before
836      !
837      INTEGER  ::   ji, jj, jk, jn  ! dummy loop indices
838      INTEGER  ::   N_in, N_out
839      INTEGER  :: item
840      ! vertical interpolation:
841      REAL(wp) :: zhtot, zwgt
842      REAL(wp), DIMENSION(k1:k2,1:jpts) :: tabin, tabin_i
843      REAL(wp), DIMENSION(k1:k2) :: z_in, h_in_i, z_in_i
844      REAL(wp), DIMENSION(1:jpk) :: h_out, z_out
845      !!----------------------------------------------------------------------
846
847      IF( before ) THEN
848
849         item = Kmm_a
850         IF( l_ini_child )   Kmm_a = Kbb_a 
851
852         DO jn = 1,jpts
853            DO jk=k1,k2
854               DO jj=j1,j2
855                 DO ji=i1,i2
856                       ptab(ji,jj,jk,jn) = ts(ji,jj,jk,jn,Kmm_a)
857                 END DO
858              END DO
859           END DO
860         END DO
861
862         IF( l_vremap .OR. l_ini_child .OR. ln_zps ) THEN
863
864            ! Fill cell depths (i.e. gdept) to be interpolated
865            ! Warning: these are masked, hence extrapolated prior interpolation.
866            DO jj=j1,j2
867               DO ji=i1,i2
868                  ptab(ji,jj,k1,jpts+1) = 0.5_wp * tmask(ji,jj,k1) * e3t(ji,jj,k1,Kmm_a)
869                  DO jk=k1+1,k2
870                     ptab(ji,jj,jk,jpts+1) = tmask(ji,jj,jk) * &
871                        & ( ptab(ji,jj,jk-1,jpts+1) + 0.5_wp * (e3t(ji,jj,jk-1,Kmm_a)+e3t(ji,jj,jk,Kmm_a)) )
872                  END DO
873               END DO
874            END DO
875         
876            ! Save ssh at last level:
877            IF (.NOT.ln_linssh) THEN
878               ptab(i1:i2,j1:j2,k2,jpts+1) = ssh(i1:i2,j1:j2,Kmm_a)*tmask(i1:i2,j1:j2,1) 
879            END IF     
880         ENDIF
881         Kmm_a = item
882
883      ELSE
884         item = Krhs_a
885         IF( l_ini_child )   Krhs_a = Kbb_a 
886
887         IF( l_vremap .OR. l_ini_child ) THEN
888            IF (ln_linssh) ptab(i1:i2,j1:j2,k2,n2) = 0._wp 
889            DO jj=j1,j2
890               DO ji=i1,i2
891                  ts(ji,jj,:,:,Krhs_a) = 0. 
892                  !
893                  ! Build vertical grids:
894                  N_in = mbkt_parent(ji,jj)
895                  ! Input grid (account for partial cells if any):
896                  DO jk=1,N_in
897                     z_in(jk) = ptab(ji,jj,jk,n2) - ptab(ji,jj,k2,n2)
898                     tabin(jk,1:jpts) = ptab(ji,jj,jk,1:jpts)
899                  END DO
900                 
901                  ! Intermediate grid:
902                  IF ( l_vremap ) THEN
903                     DO jk = 1, N_in
904                        h_in_i(jk) = e3t0_parent(ji,jj,jk) * & 
905                          &       (1._wp + ptab(ji,jj,k2,n2)/(ht0_parent(ji,jj)*ssmask(ji,jj) + 1._wp - ssmask(ji,jj)))
906                     END DO
907                     z_in_i(1) = 0.5_wp * h_in_i(1)
908                     DO jk=2,N_in
909                        z_in_i(jk) = z_in_i(jk-1) + 0.5_wp * ( h_in_i(jk) + h_in_i(jk-1) )
910                     END DO
911                     z_in_i(1:N_in) = z_in_i(1:N_in)  - ptab(ji,jj,k2,n2)
912                  ENDIF                             
913
914                  ! Output (Child) grid:
915                  N_out = mbkt(ji,jj)
916                  DO jk=1,N_out
917                     h_out(jk) = e3t(ji,jj,jk,Krhs_a)
918                  END DO
919                  z_out(1) = 0.5_wp * h_out(1)
920                  DO jk=2,N_out
921                     z_out(jk) = z_out(jk-1) + 0.5_wp * ( h_out(jk)+h_out(jk-1) )
922                  END DO
923                  IF (.NOT.ln_linssh) z_out(1:N_out) = z_out(1:N_out)  - ssh(ji,jj,Krhs_a)
924
925                  IF (N_in*N_out > 0) THEN
926                     IF( l_ini_child ) THEN
927                        CALL remap_linear(tabin(1:N_in,1:jpts),z_in(1:N_in),ts(ji,jj,1:N_out,1:jpts,Krhs_a),              &
928                                      &   z_out(1:N_out),N_in,N_out,jpts) 
929                     ELSE
930                        CALL remap_linear(tabin(1:N_in,1:jpts),z_in(1:N_in),tabin_i(1:N_in,1:jpts),                       &
931                                     &   z_in_i(1:N_in),N_in,N_in,jpts)
932                        CALL reconstructandremap(tabin_i(1:N_in,1:jpts),h_in_i(1:N_in),ts(ji,jj,1:N_out,1:jpts,Krhs_a),   &
933                                      &   h_out(1:N_out),N_in,N_out,jpts) 
934                     ENDIF
935                  ENDIF
936               END DO
937            END DO
938            Krhs_a = item
939 
940         ELSE
941         
942            IF ( Agrif_Parent(ln_zps) ) THEN ! Account for partial cells
943                                             ! linear vertical interpolation
944               DO jj=j1,j2
945                  DO ji=i1,i2
946                     !
947                     N_in  = mbkt(ji,jj)
948                     N_out = mbkt(ji,jj)
949                     z_in(1) = ptab(ji,jj,1,n2)
950                     tabin(1,1:jpts) = ptab(ji,jj,1,1:jpts)
951                     DO jk=2, N_in
952                        z_in(jk) = ptab(ji,jj,jk,n2)
953                        tabin(jk,1:jpts) = ptab(ji,jj,jk,1:jpts)
954                     END DO
955                     IF (.NOT.ln_linssh) z_in(1:N_in) = z_in(1:N_in) - ptab(ji,jj,k2,n2)
956                     z_out(1) = 0.5_wp * e3t(ji,jj,1,Krhs_a)
957                     DO jk=2, N_out
958                        z_out(jk) = z_out(jk-1) + 0.5_wp * (e3t(ji,jj,jk-1,Krhs_a) + e3t(ji,jj,jk,Krhs_a))
959                     END DO
960                     IF (.NOT.ln_linssh) z_out(1:N_out) = z_out(1:N_out) - ssh(ji,jj,Krhs_a)
961                     CALL remap_linear(tabin(1:N_in,1:jpts),z_in(1:N_in),ptab(ji,jj,1:N_out,1:jpts), &
962                                   &   z_out(1:N_out),N_in,N_out,jpts) 
963                  END DO
964               END DO
965            ENDIF
966
967            DO jn =1, jpts
968               ts(i1:i2,j1:j2,1:jpk,jn,Krhs_a) = ptab(i1:i2,j1:j2,1:jpk,jn)*tmask(i1:i2,j1:j2,1:jpk)
969            END DO
970         ENDIF
971
972      ENDIF
973      !
974   END SUBROUTINE interptsn
975
976   
977   SUBROUTINE interpsshn( ptab, i1, i2, j1, j2, before )
978      !!----------------------------------------------------------------------
979      !!                  ***  ROUTINE interpsshn  ***
980      !!---------------------------------------------------------------------- 
981      INTEGER                         , INTENT(in   ) ::   i1, i2, j1, j2
982      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) ::   ptab
983      LOGICAL                         , INTENT(in   ) ::   before
984      !
985      !!---------------------------------------------------------------------- 
986      !
987      IF( before) THEN
988         ptab(i1:i2,j1:j2) = ssh(i1:i2,j1:j2,Kmm_a)
989      ELSE
990         IF( l_ini_child ) THEN
991            ssh(i1:i2,j1:j2,Kmm_a) = ptab(i1:i2,j1:j2) * tmask(i1:i2,j1:j2,1)
992         ELSE
993            hbdy(i1:i2,j1:j2) = ptab(i1:i2,j1:j2) * tmask(i1:i2,j1:j2,1)
994         ENDIF
995      ENDIF
996      !
997   END SUBROUTINE interpsshn
998
999   
1000   SUBROUTINE interpun( ptab, i1, i2, j1, j2, k1, k2, m1, m2, before )
1001      !!----------------------------------------------------------------------
1002      !!                  *** ROUTINE interpun ***
1003      !!---------------------------------------------   
1004      !!
1005      INTEGER, INTENT(in) :: i1,i2,j1,j2,k1,k2,m1,m2
1006      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2,m1:m2), INTENT(inout) :: ptab
1007      LOGICAL, INTENT(in) :: before
1008      !!
1009      INTEGER :: ji,jj,jk
1010      REAL(wp) :: zrhoy, zhtot
1011      ! vertical interpolation:
1012      REAL(wp), DIMENSION(k1:k2) :: tabin, h_in, z_in
1013      REAL(wp), DIMENSION(1:jpk) :: h_out, z_out
1014      INTEGER  :: N_in, N_out,item
1015      REAL(wp) :: h_diff
1016      !!---------------------------------------------   
1017      !
1018      IF (before) THEN
1019
1020         item = Kmm_a
1021         IF( l_ini_child )   Kmm_a = Kbb_a     
1022
1023         DO jk=1,jpk
1024            DO jj=j1,j2
1025               DO ji=i1,i2
1026                  ptab(ji,jj,jk,1) = (e2u(ji,jj) * e3u(ji,jj,jk,Kmm_a) * uu(ji,jj,jk,Kmm_a)*umask(ji,jj,jk)) 
1027                  IF( l_vremap .OR. l_ini_child) THEN
1028                     ! Interpolate thicknesses (masked for subsequent extrapolation)
1029                     ptab(ji,jj,jk,2) = umask(ji,jj,jk) * e2u(ji,jj) * e3u(ji,jj,jk,Kmm_a)
1030                  ENDIF
1031               END DO
1032            END DO
1033         END DO
1034
1035        IF( l_vremap .OR. l_ini_child ) THEN
1036         ! Extrapolate thicknesses in partial bottom cells:
1037         ! Set them to Agrif_SpecialValue (0.). Correct bottom thicknesses are retrieved later on
1038            IF (ln_zps) THEN
1039               DO jj=j1,j2
1040                  DO ji=i1,i2
1041                     jk = mbku(ji,jj)
1042                     ptab(ji,jj,jk,2) = 0._wp
1043                  END DO
1044               END DO           
1045            END IF
1046
1047           ! Save ssh at last level:
1048           ptab(i1:i2,j1:j2,k2,2) = 0._wp
1049           IF (.NOT.ln_linssh) THEN
1050              ! This vertical sum below should be replaced by the sea-level at U-points (optimization):
1051              DO jk=1,jpk
1052                 ptab(i1:i2,j1:j2,k2,2) = ptab(i1:i2,j1:j2,k2,2) + e3u(i1:i2,j1:j2,jk,Kmm_a) * umask(i1:i2,j1:j2,jk)
1053              END DO
1054              ptab(i1:i2,j1:j2,k2,2) = ptab(i1:i2,j1:j2,k2,2) - hu_0(i1:i2,j1:j2)
1055           END IF
1056        ENDIF
1057
1058         Kmm_a = item
1059         !
1060      ELSE
1061         zrhoy = Agrif_rhoy()
1062
1063        IF( l_vremap .OR. l_ini_child) THEN
1064! VERTICAL REFINEMENT BEGIN
1065
1066            IF (ln_linssh) ptab(i1:i2,j1:j2,k2,2) = 0._wp 
1067
1068            DO ji=i1,i2
1069               DO jj=j1,j2
1070                  uu(ji,jj,:,Krhs_a) = 0._wp
1071                  N_in = mbku_parent(ji,jj)
1072                  zhtot = 0._wp
1073                  DO jk=1,N_in
1074                     !IF (jk==N_in) THEN
1075                     !   h_in(jk) = hu0_parent(ji,jj) + ptab(ji,jj,k2,2) - zhtot
1076                     !ELSE
1077                     !   h_in(jk) = ptab(ji,jj,jk,2)/(e2u(ji,jj)*zrhoy)
1078                     !ENDIF
1079                     IF ( l_vremap ) THEN
1080                        h_in(jk) = e3u0_parent(ji,jj,jk) 
1081                        zhtot = zhtot + h_in(jk)
1082                     ENDIF
1083                     IF( h_in(jk) .GT. 0. ) THEN
1084                     tabin(jk) = ptab(ji,jj,jk,1)/(e2u(ji,jj)*zrhoy*h_in(jk))
1085                     ELSE
1086                     tabin(jk) = 0.
1087                     ENDIF
1088                 END DO
1089                 z_in(1) = 0.5_wp * h_in(1) - zhtot + hu0_parent(ji,jj) 
1090                 DO jk=2,N_in
1091                    z_in(jk) = z_in(jk-1) + 0.5_wp * (h_in(jk)+h_in(jk-1))
1092                 END DO
1093                     
1094                 N_out = 0
1095                 DO jk=1,jpk
1096                    IF (umask(ji,jj,jk) == 0) EXIT
1097                    N_out = N_out + 1
1098                    h_out(N_out) = e3u(ji,jj,jk,Krhs_a)
1099                 END DO
1100
1101                 z_out(1) = 0.5_wp * h_out(1) - SUM(h_out(1:N_out)) + hu_0(ji,jj)
1102                 DO jk=2,N_out
1103                    z_out(jk) = z_out(jk-1) + 0.5_wp * (h_out(jk-1) + h_out(jk)) 
1104                 END DO 
1105
1106                 IF (N_in*N_out > 0) THEN
1107                     IF( l_ini_child ) THEN
1108                        CALL remap_linear       (tabin(1:N_in),z_in(1:N_in),uu(ji,jj,1:N_out,Krhs_a),z_out(1:N_out),N_in,N_out,1)
1109                     ELSE
1110                        CALL reconstructandremap(tabin(1:N_in),h_in(1:N_in),uu(ji,jj,1:N_out,Krhs_a),h_out(1:N_out),N_in,N_out,1)
1111                     ENDIF   
1112                 ENDIF
1113               END DO
1114            END DO
1115         ELSE
1116            DO jk = 1, jpkm1
1117               DO jj=j1,j2
1118                  uu(i1:i2,jj,jk,Krhs_a) = ptab(i1:i2,jj,jk,1) / ( zrhoy * e2u(i1:i2,jj) * e3u(i1:i2,jj,jk,Krhs_a) )
1119               END DO
1120            END DO
1121         ENDIF
1122
1123      ENDIF
1124      !
1125   END SUBROUTINE interpun
1126
1127   
1128   SUBROUTINE interpvn( ptab, i1, i2, j1, j2, k1, k2, m1, m2, before )
1129      !!----------------------------------------------------------------------
1130      !!                  *** ROUTINE interpvn ***
1131      !!----------------------------------------------------------------------
1132      !
1133      INTEGER, INTENT(in) :: i1,i2,j1,j2,k1,k2,m1,m2
1134      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2,m1:m2), INTENT(inout) :: ptab
1135      LOGICAL, INTENT(in) :: before
1136      !
1137      INTEGER :: ji,jj,jk
1138      REAL(wp) :: zrhox
1139      ! vertical interpolation:
1140      REAL(wp), DIMENSION(k1:k2) :: tabin, h_in, z_in
1141      REAL(wp), DIMENSION(1:jpk) :: h_out, z_out
1142      INTEGER  :: N_in, N_out, item
1143      REAL(wp) :: h_diff, zhtot
1144      !!---------------------------------------------   
1145      !     
1146      IF (before) THEN   
1147
1148         item = Kmm_a
1149         IF( l_ini_child )   Kmm_a = Kbb_a     
1150       
1151         DO jk=k1,k2
1152            DO jj=j1,j2
1153               DO ji=i1,i2
1154                  ptab(ji,jj,jk,1) = (e1v(ji,jj) * e3v(ji,jj,jk,Kmm_a) * vv(ji,jj,jk,Kmm_a)*vmask(ji,jj,jk))
1155                  IF( l_vremap .OR. l_ini_child) THEN
1156                     ! Interpolate thicknesses (masked for subsequent extrapolation)
1157                     ptab(ji,jj,jk,2) = vmask(ji,jj,jk) * e1v(ji,jj) * e3v(ji,jj,jk,Kmm_a)
1158                  ENDIF
1159               END DO
1160            END DO
1161         END DO
1162
1163         IF( l_vremap .OR. l_ini_child) THEN
1164         ! Extrapolate thicknesses in partial bottom cells:
1165         ! Set them to Agrif_SpecialValue (0.). Correct bottom thicknesses are retrieved later on
1166            IF (ln_zps) THEN
1167               DO jj=j1,j2
1168                  DO ji=i1,i2
1169                     jk = mbkv(ji,jj)
1170                     ptab(ji,jj,jk,2) = 0._wp
1171                  END DO
1172               END DO           
1173            END IF
1174            ! Save ssh at last level:
1175            ptab(i1:i2,j1:j2,k2,2) = 0._wp
1176            IF (.NOT.ln_linssh) THEN
1177               ! This vertical sum below should be replaced by the sea-level at V-points (optimization):
1178               DO jk=1,jpk
1179                  ptab(i1:i2,j1:j2,k2,2) = ptab(i1:i2,j1:j2,k2,2) + e3v(i1:i2,j1:j2,jk,Kmm_a) * vmask(i1:i2,j1:j2,jk)
1180               END DO
1181               ptab(i1:i2,j1:j2,k2,2) = ptab(i1:i2,j1:j2,k2,2) - hv_0(i1:i2,j1:j2)
1182            END IF
1183         ENDIF
1184         item = Kmm_a
1185
1186      ELSE       
1187         zrhox = Agrif_rhox()
1188
1189         IF( l_vremap .OR. l_ini_child ) THEN
1190
1191            IF (ln_linssh) ptab(i1:i2,j1:j2,k2,2) = 0._wp 
1192
1193            DO jj=j1,j2
1194               DO ji=i1,i2
1195                  vv(ji,jj,:,Krhs_a) = 0._wp
1196                  N_in = mbkv_parent(ji,jj)
1197                  zhtot = 0._wp
1198                  DO jk=1,N_in
1199                     !IF (jk==N_in) THEN
1200                     !   h_in(jk) = hv0_parent(ji,jj) + ptab(ji,jj,k2,2) - zhtot
1201                     !ELSE
1202                     !   h_in(jk) = ptab(ji,jj,jk,2)/(e1v(ji,jj)*zrhox)
1203                     !ENDIF
1204                     IF (l_vremap) THEN
1205                        h_in(jk) = e3v0_parent(ji,jj,jk)
1206                        zhtot = zhtot + h_in(jk)
1207                     ENDIF
1208                     IF( h_in(jk) .GT. 0. ) THEN
1209                       tabin(jk) = ptab(ji,jj,jk,1)/(e1v(ji,jj)*zrhox*h_in(jk))
1210                     ELSE
1211                       tabin(jk)  = 0.
1212                     ENDIF
1213                  END DO
1214
1215                  z_in(1) = 0.5_wp * h_in(1) - zhtot + hv0_parent(ji,jj)
1216                  DO jk=2,N_in
1217                     z_in(jk) = z_in(jk-1) + 0.5_wp * (h_in(jk-1)+h_in(jk))
1218                  END DO
1219
1220                  N_out = 0
1221                  DO jk=1,jpk
1222                     IF (vmask(ji,jj,jk) == 0) EXIT
1223                     N_out = N_out + 1
1224                     h_out(N_out) = e3v(ji,jj,jk,Krhs_a)
1225                  END DO
1226
1227                  z_out(1) = 0.5_wp * h_out(1) - SUM(h_out(1:N_out)) + hv_0(ji,jj)
1228                  DO jk=2,N_out
1229                     z_out(jk) = z_out(jk-1) + 0.5_wp * (h_out(jk-1)+h_out(jk))
1230                  END DO
1231 
1232                  IF (N_in*N_out > 0) THEN
1233                     IF( l_ini_child ) THEN
1234                        CALL remap_linear       (tabin(1:N_in),z_in(1:N_in),vv(ji,jj,1:N_out,Krhs_a),z_out(1:N_out),N_in,N_out,1)
1235                     ELSE
1236                        CALL reconstructandremap(tabin(1:N_in),h_in(1:N_in),vv(ji,jj,1:N_out,Krhs_a),h_out(1:N_out),N_in,N_out,1)
1237                     ENDIF   
1238                  ENDIF
1239               END DO
1240            END DO
1241         ELSE
1242            DO jk = 1, jpkm1
1243               vv(i1:i2,j1:j2,jk,Krhs_a) = ptab(i1:i2,j1:j2,jk,1) / ( zrhox * e1v(i1:i2,j1:j2) * e3v(i1:i2,j1:j2,jk,Krhs_a) )
1244            END DO
1245         ENDIF
1246      ENDIF
1247      !       
1248   END SUBROUTINE interpvn
1249
1250   SUBROUTINE interpunb( ptab, i1, i2, j1, j2, before)
1251      !!----------------------------------------------------------------------
1252      !!                  ***  ROUTINE interpunb  ***
1253      !!---------------------------------------------------------------------- 
1254      INTEGER                         , INTENT(in   ) ::   i1, i2, j1, j2
1255      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) ::   ptab
1256      LOGICAL                         , INTENT(in   ) ::   before
1257      !
1258      INTEGER  ::   ji, jj
1259      REAL(wp) ::   zrhoy, zrhot, zt0, zt1, ztcoeff
1260      !!---------------------------------------------------------------------- 
1261      !
1262      IF( before ) THEN
1263         ptab(i1:i2,j1:j2) = e2u(i1:i2,j1:j2) * hu(i1:i2,j1:j2,Kmm_a) * uu_b(i1:i2,j1:j2,Kmm_a)
1264      ELSE
1265         zrhoy = Agrif_Rhoy()
1266         zrhot = Agrif_rhot()
1267         ! Time indexes bounds for integration
1268         zt0 = REAL(Agrif_NbStepint()  , wp) / zrhot
1269         zt1 = REAL(Agrif_NbStepint()+1, wp) / zrhot     
1270         !
1271         DO ji = i1, i2
1272            DO jj = j1, j2
1273               IF ( utint_stage(ji,jj) < (bdy_tinterp + 1) ) THEN
1274                  IF    ( utint_stage(ji,jj) == 1  ) THEN
1275                     ztcoeff = zrhot * (  zt1**2._wp * (       zt1 - 1._wp)        &
1276                        &               - zt0**2._wp * (       zt0 - 1._wp)        )
1277                  ELSEIF( utint_stage(ji,jj) == 2  ) THEN
1278                     ztcoeff = zrhot * (  zt1        * (       zt1 - 1._wp)**2._wp &
1279                        &               - zt0        * (       zt0 - 1._wp)**2._wp )
1280                  ELSEIF( utint_stage(ji,jj) == 0  ) THEN               
1281                     ztcoeff = 1._wp
1282                  ELSE
1283                     ztcoeff = 0._wp
1284                  ENDIF
1285                  !   
1286                  ubdy(ji,jj) = ubdy(ji,jj) + ztcoeff * ptab(ji,jj)
1287                  !           
1288                  IF (( utint_stage(ji,jj) == 2 ).OR.( utint_stage(ji,jj) == 0 )) THEN
1289                     ubdy(ji,jj) = ubdy(ji,jj) / (zrhoy*e2u(ji,jj)) * umask(ji,jj,1)
1290                  ENDIF
1291                  !
1292                  utint_stage(ji,jj) = utint_stage(ji,jj) + 1
1293               ENDIF
1294            END DO
1295         END DO
1296      END IF
1297      !
1298   END SUBROUTINE interpunb
1299
1300
1301   SUBROUTINE interpvnb( ptab, i1, i2, j1, j2, before )
1302      !!----------------------------------------------------------------------
1303      !!                  ***  ROUTINE interpvnb  ***
1304      !!---------------------------------------------------------------------- 
1305      INTEGER                         , INTENT(in   ) ::   i1, i2, j1, j2
1306      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) ::   ptab
1307      LOGICAL                         , INTENT(in   ) ::   before
1308      !
1309      INTEGER  ::   ji, jj
1310      REAL(wp) ::   zrhox, zrhot, zt0, zt1, ztcoeff   
1311      !!---------------------------------------------------------------------- 
1312      !
1313      IF( before ) THEN
1314         ptab(i1:i2,j1:j2) = e1v(i1:i2,j1:j2) * hv(i1:i2,j1:j2,Kmm_a) * vv_b(i1:i2,j1:j2,Kmm_a)
1315      ELSE
1316         zrhox = Agrif_Rhox()
1317         zrhot = Agrif_rhot()
1318         ! Time indexes bounds for integration
1319         zt0 = REAL(Agrif_NbStepint()  , wp) / zrhot
1320         zt1 = REAL(Agrif_NbStepint()+1, wp) / zrhot 
1321         !     
1322         DO ji = i1, i2
1323            DO jj = j1, j2
1324               IF ( vtint_stage(ji,jj) < (bdy_tinterp + 1) ) THEN
1325                  IF    ( vtint_stage(ji,jj) == 1  ) THEN
1326                     ztcoeff = zrhot * (  zt1**2._wp * (       zt1 - 1._wp)        &
1327                        &               - zt0**2._wp * (       zt0 - 1._wp)        )
1328                  ELSEIF( vtint_stage(ji,jj) == 2  ) THEN
1329                     ztcoeff = zrhot * (  zt1        * (       zt1 - 1._wp)**2._wp &
1330                        &               - zt0        * (       zt0 - 1._wp)**2._wp )
1331                  ELSEIF( vtint_stage(ji,jj) == 0  ) THEN               
1332                     ztcoeff = 1._wp
1333                  ELSE
1334                     ztcoeff = 0._wp
1335                  ENDIF
1336                  !   
1337                  vbdy(ji,jj) = vbdy(ji,jj) + ztcoeff * ptab(ji,jj)
1338                  !           
1339                  IF (( vtint_stage(ji,jj) == 2 ).OR.( vtint_stage(ji,jj) == 0 )) THEN
1340                     vbdy(ji,jj) = vbdy(ji,jj) / (zrhox*e1v(ji,jj)) * vmask(ji,jj,1)
1341                  ENDIF
1342                  !
1343                  vtint_stage(ji,jj) = vtint_stage(ji,jj) + 1
1344               ENDIF
1345            END DO
1346         END DO         
1347      ENDIF
1348      !
1349   END SUBROUTINE interpvnb
1350
1351
1352   SUBROUTINE interpub2b( ptab, i1, i2, j1, j2, before )
1353      !!----------------------------------------------------------------------
1354      !!                  ***  ROUTINE interpub2b  ***
1355      !!---------------------------------------------------------------------- 
1356      INTEGER                         , INTENT(in   ) ::   i1, i2, j1, j2
1357      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) ::   ptab
1358      LOGICAL                         , INTENT(in   ) ::   before
1359      !
1360      INTEGER  ::   ji,jj
1361      REAL(wp) ::   zrhot, zt0, zt1, zat
1362      !!---------------------------------------------------------------------- 
1363      IF( before ) THEN
1364         IF ( ln_bt_fw ) THEN
1365            ptab(i1:i2,j1:j2) = e2u(i1:i2,j1:j2) * ub2_b(i1:i2,j1:j2)
1366         ELSE
1367            ptab(i1:i2,j1:j2) = e2u(i1:i2,j1:j2) * un_adv(i1:i2,j1:j2)
1368         ENDIF
1369      ELSE
1370         zrhot = Agrif_rhot()
1371         ! Time indexes bounds for integration
1372         zt0 = REAL(Agrif_NbStepint()  , wp) / zrhot
1373         zt1 = REAL(Agrif_NbStepint()+1, wp) / zrhot
1374         ! Polynomial interpolation coefficients:
1375         zat = zrhot * (  zt1**2._wp * (-2._wp*zt1 + 3._wp)    &
1376            &           - zt0**2._wp * (-2._wp*zt0 + 3._wp)    ) 
1377         !
1378         ubdy(i1:i2,j1:j2) = zat * ptab(i1:i2,j1:j2) 
1379         !
1380         ! Update interpolation stage:
1381         utint_stage(i1:i2,j1:j2) = 1
1382      ENDIF
1383      !
1384   END SUBROUTINE interpub2b
1385   
1386   SUBROUTINE interpub2b_const( ptab, i1, i2, j1, j2, before )
1387      !!----------------------------------------------------------------------
1388      !!                  ***  ROUTINE interpub2b_const  ***
1389      !!---------------------------------------------------------------------- 
1390      INTEGER                         , INTENT(in   ) ::   i1, i2, j1, j2
1391      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) ::   ptab
1392      LOGICAL                         , INTENT(in   ) ::   before
1393      !
1394      REAL(wp) :: zrhoy
1395      !!---------------------------------------------------------------------- 
1396      IF( before ) THEN
1397         IF ( ln_bt_fw ) THEN
1398            ptab(i1:i2,j1:j2) = e2u(i1:i2,j1:j2) * ub2_b(i1:i2,j1:j2)
1399         ELSE
1400            ptab(i1:i2,j1:j2) = e2u(i1:i2,j1:j2) * un_adv(i1:i2,j1:j2)
1401         ENDIF
1402      ELSE
1403         zrhoy = Agrif_Rhoy()
1404         !
1405         ubdy(i1:i2,j1:j2) = ptab(i1:i2,j1:j2) & 
1406                           & / (zrhoy*e2u(i1:i2,j1:j2)) * umask(i1:i2,j1:j2,1)
1407         !
1408      ENDIF
1409      !
1410   END SUBROUTINE interpub2b_const
1411
1412
1413   SUBROUTINE ub2b_cor( ptab, i1, i2, j1, j2, before )
1414      !!----------------------------------------------------------------------
1415      !!                  ***  ROUTINE ub2b_cor  ***
1416      !!---------------------------------------------------------------------- 
1417      INTEGER                         , INTENT(in   ) ::   i1, i2, j1, j2
1418      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) ::   ptab
1419      LOGICAL                         , INTENT(in   ) ::   before
1420      !
1421      INTEGER  :: ji, jj
1422      REAL(wp) :: zrhox, zrhoy, zx
1423      !!---------------------------------------------------------------------- 
1424      IF( before ) THEN
1425         ptab(:,:) = 0._wp
1426         DO ji=i1+1,i2-1
1427            DO jj=j1+1,j2
1428               ptab(ji,jj) = 0.25_wp*( ( vb2_b(ji+1,jj  )*e1v(ji+1,jj  )   & 
1429                           &            -vb2_b(ji-1,jj  )*e1v(ji-1,jj  ) ) &
1430                           &          -( vb2_b(ji+1,jj-1)*e1v(ji+1,jj-1)   &
1431                           &            -vb2_b(ji-1,jj-1)*e1v(ji-1,jj-1) ) )
1432            END DO
1433         END DO
1434      ELSE
1435         !
1436         zrhox = Agrif_Rhox() 
1437         zrhoy = Agrif_Rhoy()
1438         DO ji=i1,i2
1439            DO jj=j1,j2
1440               IF (utint_stage(ji,jj)==0) THEN
1441                  zx = 2._wp*MOD(ABS(mig0(ji)-nbghostcells-1), INT(Agrif_Rhox()))/zrhox - 1._wp 
1442                  ubdy(ji,jj) = ubdy(ji,jj) + 0.25_wp*(1._wp-zx*zx) * ptab(ji,jj) & 
1443                              &         / zrhoy *r1_e2u(ji,jj) * umask(ji,jj,1) 
1444                  utint_stage(ji,jj) = 1 
1445               ENDIF
1446            END DO
1447         END DO 
1448         !
1449      ENDIF
1450      !
1451   END SUBROUTINE ub2b_cor
1452
1453
1454   SUBROUTINE interpvb2b( ptab, i1, i2, j1, j2, before )
1455      !!----------------------------------------------------------------------
1456      !!                  ***  ROUTINE interpvb2b  ***
1457      !!---------------------------------------------------------------------- 
1458      INTEGER                         , INTENT(in   ) ::   i1, i2, j1, j2
1459      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) ::   ptab
1460      LOGICAL                         , INTENT(in   ) ::   before
1461      !
1462      INTEGER ::   ji,jj
1463      REAL(wp) ::   zrhot, zt0, zt1, zat
1464      !!---------------------------------------------------------------------- 
1465      !
1466      IF( before ) THEN
1467         IF ( ln_bt_fw ) THEN
1468            ptab(i1:i2,j1:j2) = e1v(i1:i2,j1:j2) * vb2_b(i1:i2,j1:j2)
1469         ELSE
1470            ptab(i1:i2,j1:j2) = e1v(i1:i2,j1:j2) * vn_adv(i1:i2,j1:j2)
1471         ENDIF
1472      ELSE     
1473         zrhot = Agrif_rhot()
1474         ! Time indexes bounds for integration
1475         zt0 = REAL(Agrif_NbStepint()  , wp) / zrhot
1476         zt1 = REAL(Agrif_NbStepint()+1, wp) / zrhot
1477         ! Polynomial interpolation coefficients:
1478         zat = zrhot * (  zt1**2._wp * (-2._wp*zt1 + 3._wp)    &
1479            &           - zt0**2._wp * (-2._wp*zt0 + 3._wp)    ) 
1480         !
1481         vbdy(i1:i2,j1:j2) = zat * ptab(i1:i2,j1:j2)
1482         !
1483         ! update interpolation stage:
1484         vtint_stage(i1:i2,j1:j2) = 1
1485      ENDIF
1486      !     
1487   END SUBROUTINE interpvb2b
1488
1489
1490   SUBROUTINE interpvb2b_const( ptab, i1, i2, j1, j2, before )
1491      !!----------------------------------------------------------------------
1492      !!                  ***  ROUTINE interpub2b_const  ***
1493      !!---------------------------------------------------------------------- 
1494      INTEGER                         , INTENT(in   ) ::   i1, i2, j1, j2
1495      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) ::   ptab
1496      LOGICAL                         , INTENT(in   ) ::   before
1497      !
1498      REAL(wp) :: zrhox
1499      !!---------------------------------------------------------------------- 
1500      IF( before ) THEN
1501         IF ( ln_bt_fw ) THEN
1502            ptab(i1:i2,j1:j2) = e1v(i1:i2,j1:j2) * vb2_b(i1:i2,j1:j2)
1503         ELSE
1504            ptab(i1:i2,j1:j2) = e1v(i1:i2,j1:j2) * vn_adv(i1:i2,j1:j2)
1505         ENDIF
1506      ELSE
1507         zrhox = Agrif_Rhox()
1508         !
1509         vbdy(i1:i2,j1:j2) = ptab(i1:i2,j1:j2) &
1510                           & / (zrhox*e1v(i1:i2,j1:j2)) * vmask(i1:i2,j1:j2,1)
1511         !
1512      ENDIF
1513      !
1514   END SUBROUTINE interpvb2b_const
1515
1516 
1517   SUBROUTINE vb2b_cor( ptab, i1, i2, j1, j2, before )
1518      !!----------------------------------------------------------------------
1519      !!                  ***  ROUTINE vb2b_cor  ***
1520      !!---------------------------------------------------------------------- 
1521      INTEGER                         , INTENT(in   ) ::   i1, i2, j1, j2
1522      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) ::   ptab
1523      LOGICAL                         , INTENT(in   ) ::   before
1524      !
1525      INTEGER  :: ji, jj
1526      REAL(wp) :: zrhox, zrhoy, zy
1527      !!---------------------------------------------------------------------- 
1528      IF( before ) THEN
1529         ptab(:,:) = 0._wp
1530         DO ji=i1+1,i2-1
1531            DO jj=j1+1,j2
1532               ptab(ji,jj) = 0.25_wp*( ( ub2_b(ji  ,jj+1)*e2u(ji  ,jj+1)   & 
1533                           &            -ub2_b(ji  ,jj-1)*e2u(ji  ,jj-1) ) &
1534                           &          -( ub2_b(ji-1,jj+1)*e2u(ji-1,jj+1)   &
1535                           &            -ub2_b(ji-1,jj-1)*e2u(ji-1,jj-1) ) )
1536            END DO
1537         END DO
1538      ELSE
1539         !
1540         zrhox = Agrif_Rhox() 
1541         zrhoy = Agrif_Rhoy()
1542         DO ji=i1,i2
1543            DO jj=j1,j2
1544               IF (vtint_stage(ji,jj)==0) THEN
1545                  zy = 2._wp*MOD(ABS(mjg0(jj)-nbghostcells-1), INT(Agrif_Rhoy()))/zrhoy - 1._wp 
1546                  vbdy(ji,jj) = vbdy(ji,jj) + 0.25_wp*(1._wp-zy*zy) * ptab(ji,jj) & 
1547                              &         / zrhox * r1_e1v(ji,jj) * vmask(ji,jj,1) 
1548                  vtint_stage(ji,jj) = 1 
1549               ENDIF
1550            END DO
1551         END DO 
1552         !
1553      ENDIF
1554      !
1555   END SUBROUTINE vb2b_cor
1556
1557
1558   SUBROUTINE interpe3t( ptab, i1, i2, j1, j2, k1, k2, before )
1559      !!----------------------------------------------------------------------
1560      !!                  ***  ROUTINE interpe3t  ***
1561      !!---------------------------------------------------------------------- 
1562      INTEGER                              , INTENT(in   ) :: i1, i2, j1, j2, k1, k2
1563      REAL(wp),DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: ptab
1564      LOGICAL                              , INTENT(in   ) :: before
1565      !
1566      INTEGER :: ji, jj, jk
1567      !!---------------------------------------------------------------------- 
1568      !   
1569      IF( before ) THEN
1570         ptab(i1:i2,j1:j2,k1:k2) = tmask(i1:i2,j1:j2,k1:k2) * e3t_0(i1:i2,j1:j2,k1:k2)
1571      ELSE
1572         !
1573         DO jk = k1, k2
1574            DO jj = j1, j2
1575               DO ji = i1, i2
1576                  IF( ABS( ptab(ji,jj,jk) - tmask(ji,jj,jk) * e3t_0(ji,jj,jk) ) > 1.D-2) THEN
1577                     WRITE(numout,*) ' Agrif error for e3t_0: parent , child, i, j, k ',  & 
1578                     &                 ptab(ji,jj,jk), tmask(ji,jj,jk) * e3t_0(ji,jj,jk), &
1579                     &                 mig0(ji), mjg0(jj), jk
1580                !     kindic_agr = kindic_agr + 1
1581                  ENDIF
1582               END DO
1583            END DO
1584         END DO
1585         !
1586      ENDIF
1587      !
1588   END SUBROUTINE interpe3t
1589
1590
1591   SUBROUTINE interpe3t0_vremap( ptab, i1, i2, j1, j2, k1, k2, before )
1592      !!----------------------------------------------------------------------
1593      !!                  ***  ROUTINE interpe3t0_vremap  ***
1594      !!---------------------------------------------------------------------- 
1595      INTEGER                              , INTENT(in   ) :: i1, i2, j1, j2, k1, k2
1596      REAL(wp),DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: ptab
1597      LOGICAL                              , INTENT(in   ) :: before
1598      !
1599      INTEGER :: ji, jj, jk
1600      REAL(wp) :: zh
1601      !!---------------------------------------------------------------------- 
1602      !   
1603      IF( before ) THEN
1604         IF ( ln_zps ) THEN
1605            DO jk = k1, k2
1606               DO jj = j1, j2
1607                  DO ji = i1, i2
1608                     ptab(ji, jj, jk) = e3t_1d(jk)
1609                  END DO
1610               END DO
1611            END DO
1612         ELSE
1613            ptab(i1:i2,j1:j2,k1:k2) = e3t_0(i1:i2,j1:j2,k1:k2)
1614         ENDIF
1615      ELSE
1616         !
1617         DO jk = k1, k2
1618            DO jj = j1, j2
1619               DO ji = i1, i2
1620                  e3t0_parent(ji,jj,jk) = ptab(ji,jj,jk)
1621               END DO
1622            END DO
1623         END DO
1624
1625         ! Retrieve correct scale factor at the bottom:
1626         DO jj = j1, j2
1627            DO ji = i1, i2
1628               zh = 0._wp
1629               DO jk = 1, mbkt_parent(ji, jj)-1
1630                  zh = zh + e3t0_parent(ji,jj,jk)
1631               END DO
1632               e3t0_parent(ji,jj,mbkt_parent(ji,jj)) = ht0_parent(ji, jj) - zh
1633            END DO
1634         END DO
1635         
1636      ENDIF
1637      !
1638   END SUBROUTINE interpe3t0_vremap
1639
1640
1641   SUBROUTINE interpglamt( ptab, i1, i2, j1, j2, before )
1642      !!----------------------------------------------------------------------
1643      !!                  ***  ROUTINE interpglamt  ***
1644      !!---------------------------------------------------------------------- 
1645      INTEGER                        , INTENT(in   ) :: i1, i2, j1, j2
1646      REAL(wp),DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab
1647      LOGICAL                        , INTENT(in   ) :: before
1648      !
1649      INTEGER :: ji, jj, jk
1650      REAL(wp):: ztst
1651      !!---------------------------------------------------------------------- 
1652      !   
1653      IF( before ) THEN
1654         ptab(i1:i2,j1:j2) = glamt(i1:i2,j1:j2)
1655      ELSE
1656         ztst = MAXVAL(ABS(glamt(i1:i2,j1:j2)))*1.e-4
1657         DO jj = j1, j2
1658            DO ji = i1, i2
1659               IF( ABS( ptab(ji,jj) - glamt(ji,jj) ) > ztst ) THEN
1660                  WRITE(numout,*) ' Agrif error for glamt: parent, child, i, j ', ptab(ji,jj), glamt(ji,jj), mig0(ji), mig0(jj)
1661!                  kindic_agr = kindic_agr + 1
1662               ENDIF
1663            END DO
1664         END DO
1665      ENDIF
1666      !
1667   END SUBROUTINE interpglamt
1668
1669
1670   SUBROUTINE interpgphit( ptab, i1, i2, j1, j2, before )
1671      !!----------------------------------------------------------------------
1672      !!                  ***  ROUTINE interpgphit  ***
1673      !!---------------------------------------------------------------------- 
1674      INTEGER                        , INTENT(in   ) :: i1, i2, j1, j2
1675      REAL(wp),DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab
1676      LOGICAL                        , INTENT(in   ) :: before
1677      !
1678      INTEGER :: ji, jj, jk
1679      REAL(wp):: ztst
1680      !!---------------------------------------------------------------------- 
1681      !   
1682      IF( before ) THEN
1683         ptab(i1:i2,j1:j2) = gphit(i1:i2,j1:j2)
1684      ELSE
1685         ztst = MAXVAL(ABS(gphit(i1:i2,j1:j2)))*1.e-4
1686         DO jj = j1, j2
1687            DO ji = i1, i2
1688               IF( ABS( ptab(ji,jj) - gphit(ji,jj) ) > ztst ) THEN
1689                  WRITE(numout,*) ' Agrif error for gphit: parent, child, i, j ', ptab(ji,jj), gphit(ji,jj), mig0(ji), mig0(jj)
1690!                  kindic_agr = kindic_agr + 1
1691               ENDIF
1692            END DO
1693         END DO
1694      ENDIF
1695      !
1696   END SUBROUTINE interpgphit
1697
1698
1699   SUBROUTINE interpavm( ptab, i1, i2, j1, j2, k1, k2, m1, m2, before )
1700      !!----------------------------------------------------------------------
1701      !!                  ***  ROUTINE interavm  ***
1702      !!---------------------------------------------------------------------- 
1703      INTEGER                                    , INTENT(in   ) ::   i1, i2, j1, j2, k1, k2, m1, m2
1704      REAL(wp),DIMENSION(i1:i2,j1:j2,k1:k2,m1:m2), INTENT(inout) ::   ptab
1705      LOGICAL                                    , INTENT(in   ) ::   before
1706      !
1707      INTEGER  :: ji, jj, jk
1708      INTEGER  :: N_in, N_out
1709      REAL(wp), DIMENSION(k1:k2) :: tabin, z_in
1710      REAL(wp), DIMENSION(1:jpk) :: z_out
1711      !!---------------------------------------------------------------------- 
1712      !     
1713      IF (before) THEN         
1714         DO jk=k1,k2
1715            DO jj=j1,j2
1716              DO ji=i1,i2
1717                    ptab(ji,jj,jk,1) = avm_k(ji,jj,jk)
1718              END DO
1719           END DO
1720         END DO
1721
1722         IF( l_vremap ) THEN
1723            ! Interpolate thicknesses
1724            ! Warning: these are masked, hence extrapolated prior interpolation.
1725            DO jk=k1,k2
1726               DO jj=j1,j2
1727                  DO ji=i1,i2
1728                      ptab(ji,jj,jk,2) = tmask(ji,jj,jk) * e3t(ji,jj,jk,Kmm_a)
1729                  END DO
1730               END DO
1731            END DO
1732
1733            ! Extrapolate thicknesses in partial bottom cells:
1734            ! Set them to Agrif_SpecialValue (0.). Correct bottom thicknesses are retrieved later on
1735            IF (ln_zps) THEN
1736               DO jj=j1,j2
1737                  DO ji=i1,i2
1738                      jk = mbkt(ji,jj)
1739                      ptab(ji,jj,jk,2) = 0._wp
1740                  END DO
1741               END DO           
1742            END IF
1743       
1744           ! Save ssh at last level:
1745            IF (.NOT.ln_linssh) THEN
1746               ptab(i1:i2,j1:j2,k2,2) = ssh(i1:i2,j1:j2,Kmm_a)*tmask(i1:i2,j1:j2,1) 
1747            ELSE
1748               ptab(i1:i2,j1:j2,k2,2) = 0._wp
1749            END IF     
1750          ENDIF
1751
1752      ELSE
1753
1754         IF( l_vremap ) THEN
1755            IF (ln_linssh) ptab(i1:i2,j1:j2,k2,2) = 0._wp 
1756            avm_k(i1:i2,j1:j2,k1:k2) = 0._wp
1757               
1758            DO jj = j1, j2
1759               DO ji =i1, i2
1760                  N_in = mbkt_parent(ji,jj)
1761                  IF ( tmask(ji,jj,1) == 0._wp) N_in = 0
1762                  z_in(N_in+1) = ht0_parent(ji,jj) + ptab(ji,jj,k2,2)
1763                  DO jk = N_in, 1, -1  ! Parent vertical grid               
1764                        z_in(jk) = z_in(jk+1) - ptab(ji,jj,jk,2)
1765                       tabin(jk) = ptab(ji,jj,jk,1)
1766                  END DO
1767                  N_out = mbkt(ji,jj) 
1768                  DO jk = 1, N_out        ! Child vertical grid
1769                     z_out(jk) = gdepw(ji,jj,jk,Kmm_a)
1770                  END DO
1771                  IF (N_in*N_out > 0) THEN
1772                     CALL remap_linear(tabin(1:N_in),z_in(1:N_in),avm_k(ji,jj,1:N_out),z_out(1:N_out),N_in,N_out,1)
1773                  ENDIF
1774               END DO
1775            END DO
1776         ELSE
1777            avm_k(i1:i2,j1:j2,k1:k2) = ptab (i1:i2,j1:j2,k1:k2,1)
1778         ENDIF
1779      ENDIF
1780      !
1781   END SUBROUTINE interpavm
1782
1783   
1784   SUBROUTINE interpmbkt( ptab, i1, i2, j1, j2, before )
1785      !!----------------------------------------------------------------------
1786      !!                  ***  ROUTINE interpmbkt  ***
1787      !!---------------------------------------------------------------------- 
1788      INTEGER                         , INTENT(in   ) ::   i1, i2, j1, j2
1789      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) ::   ptab
1790      LOGICAL                         , INTENT(in   ) ::   before
1791      !
1792      !!---------------------------------------------------------------------- 
1793      !
1794      IF( before) THEN
1795         ptab(i1:i2,j1:j2) = REAL(mbkt(i1:i2,j1:j2),wp)
1796      ELSE
1797         mbkt_parent(i1:i2,j1:j2) = NINT(ptab(i1:i2,j1:j2))
1798      ENDIF
1799      !
1800   END SUBROUTINE interpmbkt
1801
1802   
1803   SUBROUTINE interpht0( ptab, i1, i2, j1, j2, before )
1804      !!----------------------------------------------------------------------
1805      !!                  ***  ROUTINE interpht0  ***
1806      !!---------------------------------------------------------------------- 
1807      INTEGER                         , INTENT(in   ) ::   i1, i2, j1, j2
1808      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) ::   ptab
1809      LOGICAL                         , INTENT(in   ) ::   before
1810      !
1811      !!---------------------------------------------------------------------- 
1812      !
1813      IF( before) THEN
1814         ptab(i1:i2,j1:j2) = ht_0(i1:i2,j1:j2)
1815      ELSE
1816         ht0_parent(i1:i2,j1:j2) = ptab(i1:i2,j1:j2)
1817      ENDIF
1818      !
1819   END SUBROUTINE interpht0
1820
1821   SUBROUTINE Agrif_check_bat( iindic )
1822      !!----------------------------------------------------------------------
1823      !!                  ***  ROUTINE Agrif_check_bat  ***
1824      !!---------------------------------------------------------------------- 
1825      INTEGER, INTENT(inout) ::   iindic
1826      !!
1827      INTEGER :: ji, jj
1828      INTEGER  :: istart, iend, jstart, jend, ispon
1829      !!---------------------------------------------------------------------- 
1830      !
1831      !
1832      ! --- West --- !
1833      IF(lk_west) THEN
1834         ispon  = nn_sponge_len * Agrif_irhox()
1835         istart = nn_hls + 2                                  ! halo + land + 1
1836         iend   = nn_hls + 1 + nbghostcells + ispon           ! halo + land + nbghostcells + sponge
1837         jstart = nn_hls + 2
1838         jend   = jpjglo - nn_hls - 1
1839         DO ji = mi0(istart), mi1(iend)
1840            DO jj = mj0(jstart), mj1(jend)
1841               IF ( ABS(ht0_parent(ji,jj)-ht_0(ji,jj)) > 1.e-3 ) iindic = iindic + 1
1842            END DO
1843            DO jj = mj0(jstart), mj1(jend-1)
1844               IF ( ABS(hv0_parent(ji,jj)-hv_0(ji,jj)) > 1.e-3 ) iindic = iindic + 1
1845            END DO
1846         END DO
1847         DO ji = mi0(istart), mi1(iend-1)
1848            DO jj = mj0(jstart), mj1(jend)
1849               IF ( ABS(hu0_parent(ji,jj)-hu_0(ji,jj)) > 1.e-3 ) iindic = iindic + 1
1850            END DO
1851         END DO
1852      ENDIF
1853      !
1854      ! --- East --- !
1855      IF(lk_east) THEN
1856         ispon  = nn_sponge_len * Agrif_irhox() 
1857         istart = jpiglo - ( nn_hls + nbghostcells + ispon )  ! halo + land + nbghostcells + sponge - 1
1858         iend   = jpiglo - ( nn_hls + 1 )                     ! halo + land + 1                     - 1
1859         jstart = nn_hls + 2
1860         jend   = jpjglo - nn_hls - 1 
1861         DO ji = mi0(istart), mi1(iend)
1862            DO jj = mj0(jstart), mj1(jend)
1863               IF ( ABS(ht0_parent(ji,jj)-ht_0(ji,jj)) > 1.e-3 ) iindic = iindic + 1
1864            END DO
1865            DO jj = mj0(jstart), mj1(jend-1)
1866               IF ( ABS(hv0_parent(ji,jj)-hv_0(ji,jj)) > 1.e-3 ) iindic = iindic + 1
1867            END DO
1868         END DO
1869         DO ji = mi0(istart+1), mi1(iend-1)
1870            DO jj = mj0(jstart), mj1(jend)
1871               IF ( ABS(hu0_parent(ji,jj)-hu_0(ji,jj)) > 1.e-3 ) iindic = iindic + 1
1872            END DO
1873         END DO
1874      ENDIF
1875      !
1876      ! --- South --- !
1877      IF(lk_south) THEN
1878         ispon  = nn_sponge_len * Agrif_irhoy() 
1879         jstart = nn_hls + 2                                 ! halo + land + 1
1880         jend   = nn_hls + 1 + nbghostcells + ispon          ! halo + land + nbghostcells + sponge
1881         istart = nn_hls + 2
1882         iend   = jpiglo - nn_hls - 1
1883         DO jj = mj0(jstart), mj1(jend)
1884            DO ji = mi0(istart), mi1(iend)
1885               IF ( ABS(ht0_parent(ji,jj)-ht_0(ji,jj)) > 1.e-3 ) iindic = iindic + 1
1886            END DO
1887            DO ji = mi0(istart), mi1(iend-1)
1888               IF ( ABS(hu0_parent(ji,jj)-hu_0(ji,jj)) > 1.e-3 ) iindic = iindic + 1
1889            END DO
1890         END DO
1891         DO jj = mj0(jstart), mj1(jend-1)
1892            DO ji = mi0(istart), mi1(iend)
1893               IF ( ABS(hv0_parent(ji,jj)-hv_0(ji,jj)) > 1.e-3 ) iindic = iindic + 1
1894            END DO
1895         END DO
1896      ENDIF
1897      !
1898      ! --- North --- !
1899      IF(lk_north) THEN
1900         ispon  = nn_sponge_len * Agrif_irhoy() 
1901         jstart = jpjglo - ( nn_hls + nbghostcells + ispon)  ! halo + land + nbghostcells +sponge - 1
1902         jend   = jpjglo - ( nn_hls + 1 )                    ! halo + land + 1            - 1
1903         istart = nn_hls + 2
1904         iend   = jpiglo - nn_hls - 1
1905         DO jj = mj0(jstart), mj1(jend)
1906            DO ji = mi0(istart), mi1(iend)
1907               IF ( ABS(ht0_parent(ji,jj)-ht_0(ji,jj)) > 1.e-3 ) iindic = iindic + 1
1908            END DO
1909            DO ji = mi0(istart), mi1(iend-1)
1910               IF ( ABS(hu0_parent(ji,jj)-hu_0(ji,jj)) > 1.e-3 ) iindic = iindic + 1
1911            END DO
1912         END DO
1913         DO jj = mj0(jstart+1), mj1(jend-1)
1914            DO ji = mi0(istart), mi1(iend)
1915               IF ( ABS(hv0_parent(ji,jj)-hv_0(ji,jj)) > 1.e-3 ) iindic = iindic + 1
1916            END DO
1917         END DO
1918      ENDIF
1919      !
1920   END SUBROUTINE Agrif_check_bat
1921   
1922#else
1923   !!----------------------------------------------------------------------
1924   !!   Empty module                                          no AGRIF zoom
1925   !!----------------------------------------------------------------------
1926CONTAINS
1927   SUBROUTINE Agrif_OCE_Interp_empty
1928      WRITE(*,*)  'agrif_oce_interp : You should not have seen this print! error?'
1929   END SUBROUTINE Agrif_OCE_Interp_empty
1930#endif
1931
1932   !!======================================================================
1933END MODULE agrif_oce_interp
Note: See TracBrowser for help on using the repository browser.