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_r14116_HPC-04_mcastril_Mixed_Precision_implementation_final/src/NST – NEMO

source: NEMO/branches/2020/dev_r14116_HPC-04_mcastril_Mixed_Precision_implementation_final/src/NST/agrif_oce_interp.F90 @ 14200

Last change on this file since 14200 was 14200, checked in by mcastril, 4 years ago

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