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/2019/dev_r11233_AGRIF-05_jchanut_vert_coord_interp/src/NST – NEMO

source: NEMO/branches/2019/dev_r11233_AGRIF-05_jchanut_vert_coord_interp/src/NST/agrif_oce_interp.F90 @ 11625

Last change on this file since 11625 was 11625, checked in by jchanut, 5 years ago

#2222, add initialization to 0 of tracer open boundary data with vertical interpolation + various neutral optimizations
  • Property svn:keywords set to Id
File size: 44.8 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 
37   IMPLICIT NONE
38   PRIVATE
39
40   PUBLIC   Agrif_dyn, Agrif_ssh, Agrif_dyn_ts, Agrif_dyn_ts_flux, Agrif_ssh_ts, Agrif_dta_ts
41   PUBLIC   Agrif_tra, Agrif_avm
42   PUBLIC   interpun , interpvn
43   PUBLIC   interptsn, interpsshn, interpavm
44   PUBLIC   interpunb, interpvnb , interpub2b, interpvb2b
45   PUBLIC   interpe3t, interpumsk, interpvmsk
46
47   INTEGER ::   bdy_tinterp = 0
48
49#  include "vectopt_loop_substitute.h90"
50   !!----------------------------------------------------------------------
51   !! NEMO/NST 4.0 , NEMO Consortium (2018)
52   !! $Id$
53   !! Software governed by the CeCILL license (see ./LICENSE)
54   !!----------------------------------------------------------------------
55CONTAINS
56
57   SUBROUTINE Agrif_tra
58      !!----------------------------------------------------------------------
59      !!                  ***  ROUTINE Agrif_tra  ***
60      !!----------------------------------------------------------------------
61      !
62      IF( Agrif_Root() )   RETURN
63      !
64      Agrif_SpecialValue    = 0._wp
65      Agrif_UseSpecialValue = .TRUE.
66      !
67      CALL Agrif_Bc_variable( tsn_id, procname=interptsn )
68      !
69      Agrif_UseSpecialValue = .FALSE.
70      !
71   END SUBROUTINE Agrif_tra
72
73
74   SUBROUTINE Agrif_dyn( kt )
75      !!----------------------------------------------------------------------
76      !!                  ***  ROUTINE Agrif_DYN  ***
77      !!---------------------------------------------------------------------- 
78      INTEGER, INTENT(in) ::   kt
79      !
80      INTEGER ::   ji, jj, jk       ! dummy loop indices
81      INTEGER ::   ibdy1, jbdy1, ibdy2, jbdy2
82      REAL(wp), DIMENSION(jpi,jpj) ::   zub, zvb
83      !!---------------------------------------------------------------------- 
84      !
85      IF( Agrif_Root() )   RETURN
86      !
87      Agrif_SpecialValue    = 0._wp
88      Agrif_UseSpecialValue = ln_spc_dyn
89      !
90      CALL Agrif_Bc_variable( un_interp_id, procname=interpun )
91      CALL Agrif_Bc_variable( vn_interp_id, procname=interpvn )
92      !
93      Agrif_UseSpecialValue = .FALSE.
94      !
95      ! --- West --- !
96      ibdy1 = 2
97      ibdy2 = 1+nbghostcells 
98      !
99      IF( .NOT.ln_dynspg_ts ) THEN  ! Store transport
100         DO ji = mi0(ibdy1), mi1(ibdy2)
101            ua_b(ji,:) = 0._wp
102
103            DO jk = 1, jpkm1
104               DO jj = 1, jpj
105                  ua_b(ji,jj) = ua_b(ji,jj) + e3u_a(ji,jj,jk) * ua(ji,jj,jk) * umask(ji,jj,jk)
106               END DO
107            END DO
108
109            DO jj = 1, jpj
110               ua_b(ji,jj) = ua_b(ji,jj) * r1_hu_a(ji,jj)
111            END DO
112         END DO
113      ENDIF
114      !
115      DO ji = mi0(ibdy1), mi1(ibdy2)
116         zub(ji,:) = 0._wp    ! Correct transport
117         DO jk = 1, jpkm1
118            DO jj = 1, jpj
119               zub(ji,jj) = zub(ji,jj) & 
120                  & + e3u_a(ji,jj,jk)  * ua(ji,jj,jk)*umask(ji,jj,jk)
121            END DO
122         END DO
123         DO jj=1,jpj
124            zub(ji,jj) = zub(ji,jj) * r1_hu_a(ji,jj)
125         END DO
126           
127         DO jk = 1, jpkm1
128            DO jj = 1, jpj
129               ua(ji,jj,jk) = ( ua(ji,jj,jk) + ua_b(ji,jj)-zub(ji,jj)) * umask(ji,jj,jk)
130            END DO
131         END DO
132      END DO
133           
134      IF( ln_dynspg_ts ) THEN       ! Set tangential velocities to time splitting estimate
135         DO ji = mi0(ibdy1), mi1(ibdy2)
136            zvb(ji,:) = 0._wp
137            DO jk = 1, jpkm1
138               DO jj = 1, jpj
139                  zvb(ji,jj) = zvb(ji,jj) + e3v_a(ji,jj,jk) * va(ji,jj,jk) * vmask(ji,jj,jk)
140               END DO
141            END DO
142            DO jj = 1, jpj
143               zvb(ji,jj) = zvb(ji,jj) * r1_hv_a(ji,jj)
144            END DO
145            DO jk = 1, jpkm1
146               DO jj = 1, jpj
147                  va(ji,jj,jk) = ( va(ji,jj,jk) + va_b(ji,jj)-zvb(ji,jj))*vmask(ji,jj,jk)
148               END DO
149            END DO
150         END DO
151      ENDIF
152
153      ! --- East --- !
154      ibdy1 = jpiglo-1-nbghostcells
155      ibdy2 = jpiglo-2 
156      !
157      IF( .NOT.ln_dynspg_ts ) THEN  ! Store transport
158         DO ji = mi0(ibdy1), mi1(ibdy2)
159            ua_b(ji,:) = 0._wp
160            DO jk = 1, jpkm1
161               DO jj = 1, jpj
162                  ua_b(ji,jj) = ua_b(ji,jj) & 
163                      & + e3u_a(ji,jj,jk) * ua(ji,jj,jk) * umask(ji,jj,jk)
164               END DO
165            END DO
166            DO jj = 1, jpj
167               ua_b(ji,jj) = ua_b(ji,jj) * r1_hu_a(ji,jj)
168            END DO
169         END DO
170      ENDIF
171      !
172      DO ji = mi0(ibdy1), mi1(ibdy2)
173         zub(ji,:) = 0._wp    ! Correct transport
174         DO jk = 1, jpkm1
175            DO jj = 1, jpj
176               zub(ji,jj) = zub(ji,jj) & 
177                  & + e3u_a(ji,jj,jk)  * ua(ji,jj,jk) * umask(ji,jj,jk)
178            END DO
179         END DO
180         DO jj=1,jpj
181            zub(ji,jj) = zub(ji,jj) * r1_hu_a(ji,jj)
182         END DO
183           
184         DO jk = 1, jpkm1
185            DO jj = 1, jpj
186               ua(ji,jj,jk) = ( ua(ji,jj,jk) & 
187                 & + ua_b(ji,jj)-zub(ji,jj))*umask(ji,jj,jk)
188            END DO
189         END DO
190      END DO
191           
192      IF( ln_dynspg_ts ) THEN       ! Set tangential velocities to time splitting estimate
193         ibdy1 = jpiglo-nbghostcells
194         ibdy2 = jpiglo-1 
195         DO ji = mi0(ibdy1), mi1(ibdy2)
196            zvb(ji,:) = 0._wp
197            DO jk = 1, jpkm1
198               DO jj = 1, jpj
199                  zvb(ji,jj) = zvb(ji,jj) &
200                     & + e3v_a(ji,jj,jk) * va(ji,jj,jk) * vmask(ji,jj,jk)
201               END DO
202            END DO
203            DO jj = 1, jpj
204               zvb(ji,jj) = zvb(ji,jj) * r1_hv_a(ji,jj)
205            END DO
206            DO jk = 1, jpkm1
207               DO jj = 1, jpj
208                  va(ji,jj,jk) = ( va(ji,jj,jk) & 
209                      & + va_b(ji,jj)-zvb(ji,jj)) * vmask(ji,jj,jk)
210               END DO
211            END DO
212         END DO
213      ENDIF
214
215      ! --- South --- !
216      jbdy1 = 2
217      jbdy2 = 1+nbghostcells 
218      !
219      IF( .NOT.ln_dynspg_ts ) THEN  ! Store transport
220         DO jj = mj0(jbdy1), mj1(jbdy2)
221            va_b(:,jj) = 0._wp
222            DO jk = 1, jpkm1
223               DO ji = 1, jpi
224                  va_b(ji,jj) = va_b(ji,jj) & 
225                      & + e3v_a(ji,jj,jk) * va(ji,jj,jk) * vmask(ji,jj,jk)
226               END DO
227            END DO
228            DO ji=1,jpi
229               va_b(ji,jj) = va_b(ji,jj) * r1_hv_a(ji,jj)     
230            END DO
231         END DO
232      ENDIF
233      !
234      DO jj = mj0(jbdy1), mj1(jbdy2)
235         zvb(:,jj) = 0._wp    ! Correct transport
236         DO jk=1,jpkm1
237            DO ji=1,jpi
238               zvb(ji,jj) = zvb(ji,jj) & 
239                  & + e3v_a(ji,jj,jk) * va(ji,jj,jk) * vmask(ji,jj,jk)
240            END DO
241         END DO
242         DO ji = 1, jpi
243            zvb(ji,jj) = zvb(ji,jj) * r1_hv_a(ji,jj)
244         END DO
245
246         DO jk = 1, jpkm1
247            DO ji = 1, jpi
248               va(ji,jj,jk) = ( va(ji,jj,jk) & 
249                 & + va_b(ji,jj) - zvb(ji,jj) ) * vmask(ji,jj,jk)
250            END DO
251         END DO
252      END DO
253           
254      IF( ln_dynspg_ts ) THEN       ! Set tangential velocities to time splitting estimate
255         DO jj = mj0(jbdy1), mj1(jbdy2)
256            zub(:,jj) = 0._wp
257            DO jk = 1, jpkm1
258               DO ji = 1, jpi
259                  zub(ji,jj) = zub(ji,jj) & 
260                     & + e3u_a(ji,jj,jk) * ua(ji,jj,jk) * umask(ji,jj,jk)
261               END DO
262            END DO
263            DO ji = 1, jpi
264               zub(ji,jj) = zub(ji,jj) * r1_hu_a(ji,jj)
265            END DO
266               
267            DO jk = 1, jpkm1
268               DO ji = 1, jpi
269                  ua(ji,jj,jk) = ( ua(ji,jj,jk) & 
270                    & + ua_b(ji,jj) - zub(ji,jj) ) * umask(ji,jj,jk)
271               END DO
272            END DO
273         END DO
274      ENDIF
275
276      ! --- North --- !
277      jbdy1 = jpjglo-1-nbghostcells
278      jbdy2 = jpjglo-2 
279      !
280      IF( .NOT.ln_dynspg_ts ) THEN  ! Store transport
281         DO jj = mj0(jbdy1), mj1(jbdy2)
282            va_b(:,jj) = 0._wp
283            DO jk = 1, jpkm1
284               DO ji = 1, jpi
285                  va_b(ji,jj) = va_b(ji,jj) & 
286                      & + e3v_a(ji,jj,jk) * va(ji,jj,jk) * vmask(ji,jj,jk)
287               END DO
288            END DO
289            DO ji=1,jpi
290               va_b(ji,jj) = va_b(ji,jj) * r1_hv_a(ji,jj)
291            END DO
292         END DO
293      ENDIF
294      !
295      DO jj = mj0(jbdy1), mj1(jbdy2)
296         zvb(:,jj) = 0._wp    ! Correct transport
297         DO jk=1,jpkm1
298            DO ji=1,jpi
299               zvb(ji,jj) = zvb(ji,jj) & 
300                  & + e3v_a(ji,jj,jk) * va(ji,jj,jk) * vmask(ji,jj,jk)
301            END DO
302         END DO
303         DO ji = 1, jpi
304            zvb(ji,jj) = zvb(ji,jj) * r1_hv_a(ji,jj)
305         END DO
306
307         DO jk = 1, jpkm1
308            DO ji = 1, jpi
309               va(ji,jj,jk) = ( va(ji,jj,jk) & 
310                 & + va_b(ji,jj) - zvb(ji,jj) ) * vmask(ji,jj,jk)
311            END DO
312         END DO
313      END DO
314           
315      IF( ln_dynspg_ts ) THEN       ! Set tangential velocities to time splitting estimate
316         jbdy1 = jpjglo-nbghostcells
317         jbdy2 = jpjglo-1
318         DO jj = mj0(jbdy1), mj1(jbdy2)
319            zub(:,jj) = 0._wp
320            DO jk = 1, jpkm1
321               DO ji = 1, jpi
322                  zub(ji,jj) = zub(ji,jj) & 
323                     & + e3u_a(ji,jj,jk) * ua(ji,jj,jk) * umask(ji,jj,jk)
324               END DO
325            END DO
326            DO ji = 1, jpi
327               zub(ji,jj) = zub(ji,jj) * r1_hu_a(ji,jj)
328            END DO
329               
330            DO jk = 1, jpkm1
331               DO ji = 1, jpi
332                  ua(ji,jj,jk) = ( ua(ji,jj,jk) & 
333                    & + ua_b(ji,jj) - zub(ji,jj) ) * umask(ji,jj,jk)
334               END DO
335            END DO
336         END DO
337      ENDIF
338      !
339   END SUBROUTINE Agrif_dyn
340
341
342   SUBROUTINE Agrif_dyn_ts( jn )
343      !!----------------------------------------------------------------------
344      !!                  ***  ROUTINE Agrif_dyn_ts  ***
345      !!---------------------------------------------------------------------- 
346      INTEGER, INTENT(in) ::   jn
347      !!
348      INTEGER :: ji, jj
349      INTEGER :: istart, iend, jstart, jend
350      !!---------------------------------------------------------------------- 
351      !
352      IF( Agrif_Root() )   RETURN
353      !
354      !--- West ---!
355      istart = 2
356      iend   = nbghostcells+1
357      DO ji = mi0(istart), mi1(iend)
358         DO jj=1,jpj
359            va_e(ji,jj) = vbdy(ji,jj) * hvr_e(ji,jj)
360            ua_e(ji,jj) = ubdy(ji,jj) * hur_e(ji,jj)
361         END DO
362      END DO
363      !
364      !--- East ---!
365      istart = jpiglo-nbghostcells
366      iend   = jpiglo-1
367      DO ji = mi0(istart), mi1(iend)
368         DO jj=1,jpj
369            va_e(ji,jj) = vbdy(ji,jj) * hvr_e(ji,jj)
370         END DO
371      END DO
372      istart = jpiglo-nbghostcells-1
373      iend   = jpiglo-2
374      DO ji = mi0(istart), mi1(iend)
375         DO jj=1,jpj
376            ua_e(ji,jj) = ubdy(ji,jj) * hur_e(ji,jj)
377         END DO
378      END DO
379      !
380      !--- South ---!
381      jstart = 2
382      jend   = nbghostcells+1
383      DO jj = mj0(jstart), mj1(jend)
384         DO ji=1,jpi
385            ua_e(ji,jj) = ubdy(ji,jj) * hur_e(ji,jj)
386            va_e(ji,jj) = vbdy(ji,jj) * hvr_e(ji,jj)
387         END DO
388      END DO
389      !
390      !--- North ---!
391      jstart = jpjglo-nbghostcells
392      jend   = jpjglo-1
393      DO jj = mj0(jstart), mj1(jend)
394         DO ji=1,jpi
395            ua_e(ji,jj) = ubdy(ji,jj) * hur_e(ji,jj)
396         END DO
397      END DO
398      jstart = jpjglo-nbghostcells-1
399      jend   = jpjglo-2
400      DO jj = mj0(jstart), mj1(jend)
401         DO ji=1,jpi
402            va_e(ji,jj) = vbdy(ji,jj) * hvr_e(ji,jj)
403         END DO
404      END DO
405      !
406   END SUBROUTINE Agrif_dyn_ts
407
408   SUBROUTINE Agrif_dyn_ts_flux( jn, zu, zv )
409      !!----------------------------------------------------------------------
410      !!                  ***  ROUTINE Agrif_dyn_ts_flux  ***
411      !!---------------------------------------------------------------------- 
412      INTEGER, INTENT(in) ::   jn
413      REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) ::   zu, zv
414      !!
415      INTEGER :: ji, jj
416      INTEGER :: istart, iend, jstart, jend
417      !!---------------------------------------------------------------------- 
418      !
419      IF( Agrif_Root() )   RETURN
420      !
421      !--- West ---!
422      istart = 2
423      iend   = nbghostcells+1
424      DO ji = mi0(istart), mi1(iend)
425         DO jj=1,jpj
426            zv(ji,jj) = vbdy(ji,jj) * e1v(ji,jj)
427            zu(ji,jj) = ubdy(ji,jj) * e2u(ji,jj)
428         END DO
429      END DO
430      !
431      !--- East ---!
432      istart = jpiglo-nbghostcells
433      iend   = jpiglo-1
434      DO ji = mi0(istart), mi1(iend)
435         DO jj=1,jpj
436            zv(ji,jj) = vbdy(ji,jj) * e1v(ji,jj)
437         END DO
438      END DO
439      istart = jpiglo-nbghostcells-1
440      iend   = jpiglo-2
441      DO ji = mi0(istart), mi1(iend)
442         DO jj=1,jpj
443            zu(ji,jj) = ubdy(ji,jj) * e2u(ji,jj)
444         END DO
445      END DO
446      !
447      !--- South ---!
448      jstart = 2
449      jend   = nbghostcells+1
450      DO jj = mj0(jstart), mj1(jend)
451         DO ji=1,jpi
452            zu(ji,jj) = ubdy(ji,jj) * e2u(ji,jj)
453            zv(ji,jj) = vbdy(ji,jj) * e1v(ji,jj)
454         END DO
455      END DO
456      !
457      !--- North ---!
458      jstart = jpjglo-nbghostcells
459      jend   = jpjglo-1
460      DO jj = mj0(jstart), mj1(jend)
461         DO ji=1,jpi
462            zu(ji,jj) = ubdy(ji,jj) * e2u(ji,jj)
463         END DO
464      END DO
465      jstart = jpjglo-nbghostcells-1
466      jend   = jpjglo-2
467      DO jj = mj0(jstart), mj1(jend)
468         DO ji=1,jpi
469            zv(ji,jj) = vbdy(ji,jj) * e1v(ji,jj)
470         END DO
471      END DO
472      !
473   END SUBROUTINE Agrif_dyn_ts_flux
474
475   SUBROUTINE Agrif_dta_ts( kt )
476      !!----------------------------------------------------------------------
477      !!                  ***  ROUTINE Agrif_dta_ts  ***
478      !!---------------------------------------------------------------------- 
479      INTEGER, INTENT(in) ::   kt
480      !!
481      INTEGER :: ji, jj
482      LOGICAL :: ll_int_cons
483      !!---------------------------------------------------------------------- 
484      !
485      IF( Agrif_Root() )   RETURN
486      !
487      ll_int_cons = ln_bt_fw ! Assume conservative temporal integration in the forward case only
488      !
489      ! Enforce volume conservation if no time refinement: 
490      IF ( Agrif_rhot()==1 ) ll_int_cons=.TRUE. 
491      !
492      ! Interpolate barotropic fluxes
493      Agrif_SpecialValue = 0._wp
494      Agrif_UseSpecialValue = ln_spc_dyn
495      !
496      ! Set bdy time interpolation stage to 0 (latter incremented locally do deal with corners)
497      utint_stage(:,:) = 0
498      vtint_stage(:,:) = 0
499      !
500      IF( ll_int_cons ) THEN  ! Conservative interpolation
501         ! order matters here !!!!!!
502         CALL Agrif_Bc_variable( ub2b_interp_id, calledweight=1._wp, procname=interpub2b ) ! Time integrated
503         CALL Agrif_Bc_variable( vb2b_interp_id, calledweight=1._wp, procname=interpvb2b ) 
504         !
505         bdy_tinterp = 1
506         CALL Agrif_Bc_variable( unb_id        , calledweight=1._wp, procname=interpunb  ) ! After
507         CALL Agrif_Bc_variable( vnb_id        , calledweight=1._wp, procname=interpvnb  ) 
508         !
509         bdy_tinterp = 2
510         CALL Agrif_Bc_variable( unb_id        , calledweight=0._wp, procname=interpunb  ) ! Before
511         CALL Agrif_Bc_variable( vnb_id        , calledweight=0._wp, procname=interpvnb  )   
512      ELSE ! Linear interpolation
513         !
514         ubdy(:,:) = 0._wp   ;   vbdy(:,:) = 0._wp 
515         CALL Agrif_Bc_variable( unb_id, procname=interpunb )
516         CALL Agrif_Bc_variable( vnb_id, procname=interpvnb )
517      ENDIF
518      Agrif_UseSpecialValue = .FALSE.
519      !
520   END SUBROUTINE Agrif_dta_ts
521
522
523   SUBROUTINE Agrif_ssh( kt )
524      !!----------------------------------------------------------------------
525      !!                  ***  ROUTINE Agrif_ssh  ***
526      !!---------------------------------------------------------------------- 
527      INTEGER, INTENT(in) ::   kt
528      !
529      INTEGER  :: ji, jj
530      INTEGER  :: istart, iend, jstart, jend
531      !!---------------------------------------------------------------------- 
532      !
533      IF( Agrif_Root() )   RETURN
534      !     
535      ! Linear time interpolation of sea level
536      !
537      Agrif_SpecialValue    = 0._wp
538      Agrif_UseSpecialValue = .TRUE.
539      CALL Agrif_Bc_variable(sshn_id, procname=interpsshn )
540      Agrif_UseSpecialValue = .FALSE.
541      !
542      ! --- West --- !
543      istart = 2
544      iend   = 1 + nbghostcells
545      DO ji = mi0(istart), mi1(iend)
546         DO jj = 1, jpj
547            ssha(ji,jj) = hbdy(ji,jj)
548         ENDDO
549      ENDDO
550      !
551      ! --- East --- !
552      istart = jpiglo - nbghostcells
553      iend   = jpiglo - 1
554      DO ji = mi0(istart), mi1(iend)
555         DO jj = 1, jpj
556            ssha(ji,jj) = hbdy(ji,jj)
557         ENDDO
558      ENDDO
559      !
560      ! --- South --- !
561      jstart = 2
562      jend   = 1 + nbghostcells
563      DO jj = mj0(jstart), mj1(jend)
564         DO ji = 1, jpi
565            ssha(ji,jj) = hbdy(ji,jj)
566         ENDDO
567      ENDDO
568      !
569      ! --- North --- !
570      jstart = jpjglo - nbghostcells
571      jend   = jpjglo - 1
572      DO jj = mj0(jstart), mj1(jend)
573         DO ji = 1, jpi
574            ssha(ji,jj) = hbdy(ji,jj)
575         ENDDO
576      ENDDO
577      !
578   END SUBROUTINE Agrif_ssh
579
580
581   SUBROUTINE Agrif_ssh_ts( jn )
582      !!----------------------------------------------------------------------
583      !!                  ***  ROUTINE Agrif_ssh_ts  ***
584      !!---------------------------------------------------------------------- 
585      INTEGER, INTENT(in) ::   jn
586      !!
587      INTEGER :: ji, jj
588      INTEGER  :: istart, iend, jstart, jend
589      !!---------------------------------------------------------------------- 
590      !
591      IF( Agrif_Root() )   RETURN
592      !
593      ! --- West --- !
594      istart = 2
595      iend   = 1+nbghostcells
596      DO ji = mi0(istart), mi1(iend)
597         DO jj = 1, jpj
598            ssha_e(ji,jj) = hbdy(ji,jj)
599         ENDDO
600      ENDDO
601      !
602      ! --- East --- !
603      istart = jpiglo - nbghostcells
604      iend   = jpiglo - 1
605      DO ji = mi0(istart), mi1(iend)
606         DO jj = 1, jpj
607            ssha_e(ji,jj) = hbdy(ji,jj)
608         ENDDO
609      ENDDO
610      !
611      ! --- South --- !
612      jstart = 2
613      jend   = 1+nbghostcells
614      DO jj = mj0(jstart), mj1(jend)
615         DO ji = 1, jpi
616            ssha_e(ji,jj) = hbdy(ji,jj)
617         ENDDO
618      ENDDO
619      !
620      ! --- North --- !
621      jstart = jpjglo - nbghostcells
622      jend   = jpjglo - 1
623      DO jj = mj0(jstart), mj1(jend)
624         DO ji = 1, jpi
625            ssha_e(ji,jj) = hbdy(ji,jj)
626         ENDDO
627      ENDDO
628      !
629   END SUBROUTINE Agrif_ssh_ts
630
631   SUBROUTINE Agrif_avm
632      !!----------------------------------------------------------------------
633      !!                  ***  ROUTINE Agrif_avm  ***
634      !!---------------------------------------------------------------------- 
635      REAL(wp) ::   zalpha
636      !!---------------------------------------------------------------------- 
637      !
638      IF( Agrif_Root() )   RETURN
639      !
640      zalpha = 1._wp ! JC: proper time interpolation impossible 
641                     ! => use last available value from parent
642      !
643      Agrif_SpecialValue    = 0.e0
644      Agrif_UseSpecialValue = .TRUE.
645      !
646      CALL Agrif_Bc_variable( avm_id, calledweight=zalpha, procname=interpavm )       
647      !
648      Agrif_UseSpecialValue = .FALSE.
649      !
650   END SUBROUTINE Agrif_avm
651   
652
653   SUBROUTINE interptsn( ptab, i1, i2, j1, j2, k1, k2, n1, n2, before )
654      !!----------------------------------------------------------------------
655      !!                  *** ROUTINE interptsn ***
656      !!----------------------------------------------------------------------
657      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) ::   ptab
658      INTEGER                                     , INTENT(in   ) ::   i1, i2, j1, j2, k1, k2, n1, n2
659      LOGICAL                                     , INTENT(in   ) ::   before
660      !
661      INTEGER  ::   ji, jj, jk, jn  ! dummy loop indices
662      INTEGER  ::   N_in, N_out
663      ! vertical interpolation:
664      REAL(wp), DIMENSION(i1:i2,j1:j2,1:jpk,1:jpts) :: ptab_child
665      REAL(wp), DIMENSION(k1:k2,1:jpts) :: tabin
666      REAL(wp), DIMENSION(k1:k2) :: h_in
667      REAL(wp), DIMENSION(1:jpk) :: h_out
668      !!----------------------------------------------------------------------
669
670      IF( before ) THEN         
671         DO jn = 1,jpts
672            DO jk=k1,k2
673               DO jj=j1,j2
674                 DO ji=i1,i2
675                       ptab(ji,jj,jk,jn) = tsn(ji,jj,jk,jn)
676                 END DO
677              END DO
678           END DO
679        END DO
680
681# if defined key_vertical
682        DO jk=k1,k2
683           DO jj=j1,j2
684              DO ji=i1,i2
685                 ptab(ji,jj,jk,jpts+1) = tmask(ji,jj,jk) * e3t_n(ji,jj,jk) 
686              END DO
687           END DO
688        END DO
689# endif
690      ELSE 
691
692# if defined key_vertical             
693         DO jj=j1,j2
694            DO ji=i1,i2
695               ptab_child(ji,jj,:,:) = 0._wp
696               N_in = 0
697               DO jk=k1,k2 !k2 = jpk of parent grid
698                  IF (ptab(ji,jj,jk,n2) == 0) EXIT
699                  N_in = N_in + 1
700                  tabin(jk,:) = ptab(ji,jj,jk,n1:n2-1)
701                  h_in(N_in) = ptab(ji,jj,jk,n2)
702               END DO
703               N_out = 0
704               DO jk=1,jpk ! jpk of child grid
705                  IF (tmask(ji,jj,jk) == 0) EXIT
706                  N_out = N_out + 1
707                  h_out(jk) = e3t_a(ji,jj,jk)
708               ENDDO
709               IF (N_in*N_out > 0) THEN
710                  CALL reconstructandremap(tabin(1:N_in,1:jpts),h_in(1:N_in),ptab_child(ji,jj,1:N_out,1:jpts),h_out(1:N_out),N_in,N_out,jpts)
711               ENDIF
712            ENDDO
713         ENDDO
714# else
715         ptab_child(i1:i2,j1:j2,1:jpk,1:jpts) = ptab(i1:i2,j1:j2,1:jpk,1:jpts)
716# endif
717         !
718         DO jn=1, jpts
719            tsa(i1:i2,j1:j2,1:jpk,jn)=ptab_child(i1:i2,j1:j2,1:jpk,jn)*tmask(i1:i2,j1:j2,1:jpk) 
720         END DO
721
722      ENDIF
723      !
724   END SUBROUTINE interptsn
725
726   SUBROUTINE interpsshn( ptab, i1, i2, j1, j2, before )
727      !!----------------------------------------------------------------------
728      !!                  ***  ROUTINE interpsshn  ***
729      !!---------------------------------------------------------------------- 
730      INTEGER                         , INTENT(in   ) ::   i1, i2, j1, j2
731      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) ::   ptab
732      LOGICAL                         , INTENT(in   ) ::   before
733      !
734      !!---------------------------------------------------------------------- 
735      !
736      IF( before) THEN
737         ptab(i1:i2,j1:j2) = sshn(i1:i2,j1:j2)
738      ELSE
739         hbdy(i1:i2,j1:j2) = ptab(i1:i2,j1:j2) * tmask(i1:i2,j1:j2,1)
740      ENDIF
741      !
742   END SUBROUTINE interpsshn
743
744   SUBROUTINE interpun( ptab, i1, i2, j1, j2, k1, k2, m1, m2, before )
745      !!----------------------------------------------------------------------
746      !!                  *** ROUTINE interpun ***
747      !!---------------------------------------------   
748      !!
749      INTEGER, INTENT(in) :: i1,i2,j1,j2,k1,k2,m1,m2
750      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2,m1:m2), INTENT(inout) :: ptab
751      LOGICAL, INTENT(in) :: before
752      !!
753      INTEGER :: ji,jj,jk
754      REAL(wp) :: zrhoy
755      ! vertical interpolation:
756      REAL(wp), DIMENSION(k1:k2) :: tabin, h_in
757      REAL(wp), DIMENSION(1:jpk) :: h_out
758      INTEGER  :: N_in, N_out
759      REAL(wp) :: h_diff
760      !!---------------------------------------------   
761      !
762      IF (before) THEN
763         DO jk=1,jpk
764            DO jj=j1,j2
765               DO ji=i1,i2
766                  ptab(ji,jj,jk,1) = (e2u(ji,jj) * e3u_n(ji,jj,jk) * un(ji,jj,jk)*umask(ji,jj,jk)) 
767# if defined key_vertical
768                  ptab(ji,jj,jk,2) = (umask(ji,jj,jk) * e2u(ji,jj) * e3u_n(ji,jj,jk))
769# endif
770               END DO
771            END DO
772         END DO
773      ELSE
774         zrhoy = Agrif_rhoy()
775# if defined key_vertical
776! VERTICAL REFINEMENT BEGIN
777
778         DO ji=i1,i2
779            DO jj=j1,j2
780               N_in = 0
781               DO jk=k1,k2
782                  IF (ptab(ji,jj,jk,2) == 0) EXIT
783                  N_in = N_in + 1
784                  tabin(jk) = ptab(ji,jj,jk,1)/ptab(ji,jj,jk,2)
785                  h_in(N_in) = ptab(ji,jj,jk,2)/(e2u(ji,jj)*zrhoy) 
786              ENDDO
787         
788              IF (N_in == 0) THEN
789                 ua(ji,jj,:) = 0._wp
790                 CYCLE
791              ENDIF
792         
793              N_out = 0
794              DO jk=1,jpk
795                 if (umask(ji,jj,jk) == 0) EXIT
796                 N_out = N_out + 1
797                 h_out(N_out) = e3u_a(ji,jj,jk)
798              ENDDO
799         
800              IF (N_out == 0) THEN
801                 ua(ji,jj,:) = 0._wp
802                 CYCLE
803              ENDIF
804         
805              CALL reconstructandremap(tabin(1:N_in),h_in(1:N_in),ua(ji,jj,1:N_out),h_out(1:N_out),N_in,N_out,1)
806            ENDDO
807         ENDDO
808
809# else
810         DO jk = 1, jpkm1
811            DO jj=j1,j2
812               ua(i1:i2,jj,jk) = ptab(i1:i2,jj,jk,1) / ( zrhoy * e2u(i1:i2,jj) * e3u_a(i1:i2,jj,jk) )
813            END DO
814         END DO
815# endif
816
817      ENDIF
818      !
819   END SUBROUTINE interpun
820
821   SUBROUTINE interpvn( ptab, i1, i2, j1, j2, k1, k2, m1, m2, before )
822      !!----------------------------------------------------------------------
823      !!                  *** ROUTINE interpvn ***
824      !!----------------------------------------------------------------------
825      !
826      INTEGER, INTENT(in) :: i1,i2,j1,j2,k1,k2,m1,m2
827      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2,m1:m2), INTENT(inout) :: ptab
828      LOGICAL, INTENT(in) :: before
829      !
830      INTEGER :: ji,jj,jk
831      REAL(wp) :: zrhox
832      ! vertical interpolation:
833      REAL(wp), DIMENSION(k1:k2) :: tabin, h_in
834      REAL(wp), DIMENSION(1:jpk) :: h_out
835      INTEGER  :: N_in, N_out
836      REAL(wp) :: h_diff
837      !!---------------------------------------------   
838      !     
839      IF (before) THEN         
840         DO jk=k1,k2
841            DO jj=j1,j2
842               DO ji=i1,i2
843                  ptab(ji,jj,jk,1) = (e1v(ji,jj) * e3v_n(ji,jj,jk) * vn(ji,jj,jk)*vmask(ji,jj,jk))
844# if defined key_vertical
845                  ptab(ji,jj,jk,2) = vmask(ji,jj,jk) * e1v(ji,jj) * e3v_n(ji,jj,jk)
846# endif
847               END DO
848            END DO
849         END DO
850      ELSE       
851         zrhox = Agrif_rhox()
852# if defined key_vertical
853
854         DO jj=j1,j2
855            DO ji=i1,i2
856               N_in = 0
857               DO jk=k1,k2
858                  if (ptab(ji,jj,jk,2) == 0) EXIT
859                  N_in = N_in + 1
860                  tabin(jk) = ptab(ji,jj,jk,1)/ptab(ji,jj,jk,2)
861                  h_in(N_in) = ptab(ji,jj,jk,2)/(e1v(ji,jj)*zrhox)
862               END DO
863               IF (N_in == 0) THEN
864                  va(ji,jj,:) = 0._wp
865                  CYCLE
866               ENDIF
867         
868               N_out = 0
869               DO jk=1,jpk
870                  if (vmask(ji,jj,jk) == 0) EXIT
871                  N_out = N_out + 1
872                  h_out(N_out) = e3v_a(ji,jj,jk)
873               END DO
874               IF (N_out == 0) THEN
875                 va(ji,jj,:) = 0._wp
876                 CYCLE
877               ENDIF
878               call reconstructandremap(tabin(1:N_in),h_in(1:N_in),va(ji,jj,1:N_out),h_out(1:N_out),N_in,N_out,1)
879            END DO
880         END DO
881# else
882         DO jk = 1, jpkm1
883            va(i1:i2,j1:j2,jk) = ptab(i1:i2,j1:j2,jk,1) / ( zrhox * e1v(i1:i2,j1:j2) * e3v_a(i1:i2,j1:j2,jk) )
884         END DO
885# endif
886      ENDIF
887      !       
888   END SUBROUTINE interpvn
889
890   SUBROUTINE interpunb( ptab, i1, i2, j1, j2, before)
891      !!----------------------------------------------------------------------
892      !!                  ***  ROUTINE interpunb  ***
893      !!---------------------------------------------------------------------- 
894      INTEGER                         , INTENT(in   ) ::   i1, i2, j1, j2
895      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) ::   ptab
896      LOGICAL                         , INTENT(in   ) ::   before
897      !
898      INTEGER  ::   ji, jj
899      REAL(wp) ::   zrhoy, zrhot, zt0, zt1, ztcoeff
900      !!---------------------------------------------------------------------- 
901      !
902      IF( before ) THEN
903         ptab(i1:i2,j1:j2) = e2u(i1:i2,j1:j2) * hu_n(i1:i2,j1:j2) * un_b(i1:i2,j1:j2)
904      ELSE
905         zrhoy = Agrif_Rhoy()
906         zrhot = Agrif_rhot()
907         ! Time indexes bounds for integration
908         zt0 = REAL(Agrif_NbStepint()  , wp) / zrhot
909         zt1 = REAL(Agrif_NbStepint()+1, wp) / zrhot     
910         !
911         DO ji = i1, i2
912            DO jj = j1, j2
913               IF ( utint_stage(ji,jj) < (bdy_tinterp + 1) ) THEN
914                  IF    ( utint_stage(ji,jj) == 1  ) THEN
915                     ztcoeff = zrhot * (  zt1**2._wp * (       zt1 - 1._wp)        &
916                        &               - zt0**2._wp * (       zt0 - 1._wp)        )
917                  ELSEIF( utint_stage(ji,jj) == 2  ) THEN
918                     ztcoeff = zrhot * (  zt1        * (       zt1 - 1._wp)**2._wp &
919                        &               - zt0        * (       zt0 - 1._wp)**2._wp )
920                  ELSEIF( utint_stage(ji,jj) == 0  ) THEN               
921                     ztcoeff = 1._wp
922                  ELSE
923                     ztcoeff = 0._wp
924                  ENDIF
925                  !   
926                  ubdy(ji,jj) = ubdy(ji,jj) + ztcoeff * ptab(ji,jj)
927                  !           
928                  IF (( utint_stage(ji,jj) == 2 ).OR.( utint_stage(ji,jj) == 0 )) THEN
929                     ubdy(ji,jj) = ubdy(ji,jj) / (zrhoy*e2u(ji,jj)) * umask(ji,jj,1)
930                  ENDIF
931                  !
932                  utint_stage(ji,jj) = utint_stage(ji,jj) + 1
933               ENDIF
934            END DO
935         END DO
936      END IF
937      !
938   END SUBROUTINE interpunb
939
940
941   SUBROUTINE interpvnb( ptab, i1, i2, j1, j2, before )
942      !!----------------------------------------------------------------------
943      !!                  ***  ROUTINE interpvnb  ***
944      !!---------------------------------------------------------------------- 
945      INTEGER                         , INTENT(in   ) ::   i1, i2, j1, j2
946      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) ::   ptab
947      LOGICAL                         , INTENT(in   ) ::   before
948      !
949      INTEGER  ::   ji, jj
950      REAL(wp) ::   zrhox, zrhot, zt0, zt1, ztcoeff   
951      !!---------------------------------------------------------------------- 
952      !
953      IF( before ) THEN
954         ptab(i1:i2,j1:j2) = e1v(i1:i2,j1:j2) * hv_n(i1:i2,j1:j2) * vn_b(i1:i2,j1:j2)
955      ELSE
956         zrhox = Agrif_Rhox()
957         zrhot = Agrif_rhot()
958         ! Time indexes bounds for integration
959         zt0 = REAL(Agrif_NbStepint()  , wp) / zrhot
960         zt1 = REAL(Agrif_NbStepint()+1, wp) / zrhot 
961         !     
962         DO ji = i1, i2
963            DO jj = j1, j2
964               IF ( vtint_stage(ji,jj) < (bdy_tinterp + 1) ) THEN
965                  IF    ( vtint_stage(ji,jj) == 1  ) THEN
966                     ztcoeff = zrhot * (  zt1**2._wp * (       zt1 - 1._wp)        &
967                        &               - zt0**2._wp * (       zt0 - 1._wp)        )
968                  ELSEIF( vtint_stage(ji,jj) == 2  ) THEN
969                     ztcoeff = zrhot * (  zt1        * (       zt1 - 1._wp)**2._wp &
970                        &               - zt0        * (       zt0 - 1._wp)**2._wp )
971                  ELSEIF( vtint_stage(ji,jj) == 0  ) THEN               
972                     ztcoeff = 1._wp
973                  ELSE
974                     ztcoeff = 0._wp
975                  ENDIF
976                  !   
977                  vbdy(ji,jj) = vbdy(ji,jj) + ztcoeff * ptab(ji,jj)
978                  !           
979                  IF (( vtint_stage(ji,jj) == 2 ).OR.( vtint_stage(ji,jj) == 0 )) THEN
980                     vbdy(ji,jj) = vbdy(ji,jj) / (zrhox*e1v(ji,jj)) * vmask(ji,jj,1)
981                  ENDIF
982                  !
983                  vtint_stage(ji,jj) = vtint_stage(ji,jj) + 1
984               ENDIF
985            END DO
986         END DO         
987      ENDIF
988      !
989   END SUBROUTINE interpvnb
990
991
992   SUBROUTINE interpub2b( ptab, i1, i2, j1, j2, before )
993      !!----------------------------------------------------------------------
994      !!                  ***  ROUTINE interpub2b  ***
995      !!---------------------------------------------------------------------- 
996      INTEGER                         , INTENT(in   ) ::   i1, i2, j1, j2
997      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) ::   ptab
998      LOGICAL                         , INTENT(in   ) ::   before
999      !
1000      INTEGER  ::   ji,jj
1001      REAL(wp) ::   zrhot, zt0, zt1, zat
1002      !!---------------------------------------------------------------------- 
1003      IF( before ) THEN
1004         IF ( ln_bt_fw ) THEN
1005            ptab(i1:i2,j1:j2) = e2u(i1:i2,j1:j2) * ub2_b(i1:i2,j1:j2)
1006         ELSE
1007            ptab(i1:i2,j1:j2) = e2u(i1:i2,j1:j2) * un_adv(i1:i2,j1:j2)
1008         ENDIF
1009      ELSE
1010         zrhot = Agrif_rhot()
1011         ! Time indexes bounds for integration
1012         zt0 = REAL(Agrif_NbStepint()  , wp) / zrhot
1013         zt1 = REAL(Agrif_NbStepint()+1, wp) / zrhot
1014         ! Polynomial interpolation coefficients:
1015         zat = zrhot * (  zt1**2._wp * (-2._wp*zt1 + 3._wp)    &
1016            &           - zt0**2._wp * (-2._wp*zt0 + 3._wp)    ) 
1017         !
1018         ubdy(i1:i2,j1:j2) = zat * ptab(i1:i2,j1:j2) 
1019         !
1020         ! Update interpolation stage:
1021         utint_stage(i1:i2,j1:j2) = 1
1022      ENDIF
1023      !
1024   END SUBROUTINE interpub2b
1025   
1026
1027   SUBROUTINE interpvb2b( ptab, i1, i2, j1, j2, before )
1028      !!----------------------------------------------------------------------
1029      !!                  ***  ROUTINE interpvb2b  ***
1030      !!---------------------------------------------------------------------- 
1031      INTEGER                         , INTENT(in   ) ::   i1, i2, j1, j2
1032      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) ::   ptab
1033      LOGICAL                         , INTENT(in   ) ::   before
1034      !
1035      INTEGER ::   ji,jj
1036      REAL(wp) ::   zrhot, zt0, zt1, zat
1037      !!---------------------------------------------------------------------- 
1038      !
1039      IF( before ) THEN
1040         IF ( ln_bt_fw ) THEN
1041            ptab(i1:i2,j1:j2) = e1v(i1:i2,j1:j2) * vb2_b(i1:i2,j1:j2)
1042         ELSE
1043            ptab(i1:i2,j1:j2) = e1v(i1:i2,j1:j2) * vn_adv(i1:i2,j1:j2)
1044         ENDIF
1045      ELSE     
1046         zrhot = Agrif_rhot()
1047         ! Time indexes bounds for integration
1048         zt0 = REAL(Agrif_NbStepint()  , wp) / zrhot
1049         zt1 = REAL(Agrif_NbStepint()+1, wp) / zrhot
1050         ! Polynomial interpolation coefficients:
1051         zat = zrhot * (  zt1**2._wp * (-2._wp*zt1 + 3._wp)    &
1052            &           - zt0**2._wp * (-2._wp*zt0 + 3._wp)    ) 
1053         !
1054         vbdy(i1:i2,j1:j2) = zat * ptab(i1:i2,j1:j2)
1055         !
1056         ! update interpolation stage:
1057         vtint_stage(i1:i2,j1:j2) = 1
1058      ENDIF
1059      !     
1060   END SUBROUTINE interpvb2b
1061
1062
1063   SUBROUTINE interpe3t( ptab, i1, i2, j1, j2, k1, k2, before, nb, ndir )
1064      !!----------------------------------------------------------------------
1065      !!                  ***  ROUTINE interpe3t  ***
1066      !!---------------------------------------------------------------------- 
1067      INTEGER                              , INTENT(in   ) :: i1, i2, j1, j2, k1, k2
1068      REAL(wp),DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: ptab
1069      LOGICAL                              , INTENT(in   ) :: before
1070      INTEGER                              , INTENT(in   ) :: nb , ndir
1071      !
1072      INTEGER :: ji, jj, jk
1073      LOGICAL :: western_side, eastern_side, northern_side, southern_side
1074      !!---------------------------------------------------------------------- 
1075      !   
1076      IF( before ) THEN
1077         ptab(i1:i2,j1:j2,k1:k2) = tmask(i1:i2,j1:j2,k1:k2) * e3t_0(i1:i2,j1:j2,k1:k2)
1078      ELSE
1079         western_side  = (nb == 1).AND.(ndir == 1)
1080         eastern_side  = (nb == 1).AND.(ndir == 2)
1081         southern_side = (nb == 2).AND.(ndir == 1)
1082         northern_side = (nb == 2).AND.(ndir == 2)
1083         !
1084         DO jk = k1, k2
1085            DO jj = j1, j2
1086               DO ji = i1, i2
1087                  !
1088                  IF( ABS( ptab(ji,jj,jk) - tmask(ji,jj,jk) * e3t_0(ji,jj,jk) ) > 1.D-2) THEN
1089                     IF (western_side.AND.(ptab(i1+nbghostcells-1,jj,jk)>0._wp)) THEN
1090                        WRITE(numout,*) 'ERROR bathymetry merge at the western border ji,jj,jk ', ji+nimpp-1,jj+njmpp-1,jk
1091                        WRITE(numout,*)  ptab(ji,jj,jk), e3t_0(ji,jj,jk) 
1092                        kindic_agr = kindic_agr + 1
1093                     ELSEIF (eastern_side.AND.(ptab(i2-nbghostcells+1,jj,jk)>0._wp)) THEN
1094                        WRITE(numout,*) 'ERROR bathymetry merge at the eastern border ji,jj,jk ', ji+nimpp-1,jj+njmpp-1,jk
1095                        WRITE(numout,*)  ptab(ji,jj,jk), e3t_0(ji,jj,jk)
1096                        kindic_agr = kindic_agr + 1
1097                     ELSEIF (southern_side.AND.(ptab(ji,j1+nbghostcells-1,jk)>0._wp)) THEN
1098                        WRITE(numout,*) 'ERROR bathymetry merge at the southern border ji,jj,jk', ji+nimpp-1,jj+njmpp-1,jk
1099                        WRITE(numout,*)  ptab(ji,jj,jk), e3t_0(ji,jj,jk)
1100                        kindic_agr = kindic_agr + 1
1101                     ELSEIF (northern_side.AND.(ptab(ji,j2-nbghostcells+1,jk)>0._wp)) THEN
1102                        WRITE(numout,*) 'ERROR bathymetry merge at the northen border ji,jj,jk', ji+nimpp-1,jj+njmpp-1,jk
1103                        WRITE(numout,*)  ptab(ji,jj,jk), e3t_0(ji,jj,jk)
1104                        kindic_agr = kindic_agr + 1
1105                     ENDIF
1106                  ENDIF
1107               END DO
1108            END DO
1109         END DO
1110         !
1111      ENDIF
1112      !
1113   END SUBROUTINE interpe3t
1114
1115
1116   SUBROUTINE interpumsk( ptab, i1, i2, j1, j2, k1, k2, before, nb, ndir )
1117      !!----------------------------------------------------------------------
1118      !!                  ***  ROUTINE interpumsk  ***
1119      !!---------------------------------------------------------------------- 
1120      INTEGER                              , INTENT(in   ) ::   i1, i2, j1, j2, k1, k2
1121      REAL(wp),DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) ::   ptab
1122      LOGICAL                              , INTENT(in   ) ::   before
1123      INTEGER                              , INTENT(in   ) ::   nb , ndir
1124      !
1125      INTEGER ::   ji, jj, jk
1126      LOGICAL ::   western_side, eastern_side   
1127      !!---------------------------------------------------------------------- 
1128      !   
1129      IF( before ) THEN
1130         ptab(i1:i2,j1:j2,k1:k2) = umask(i1:i2,j1:j2,k1:k2)
1131      ELSE
1132         western_side = (nb == 1).AND.(ndir == 1)
1133         eastern_side = (nb == 1).AND.(ndir == 2)
1134         DO jk = k1, k2
1135            DO jj = j1, j2
1136               DO ji = i1, i2
1137                   ! Velocity mask at boundary edge points:
1138                  IF (ABS(ptab(ji,jj,jk) - umask(ji,jj,jk)) > 1.D-2) THEN
1139                     IF (western_side) THEN
1140                        WRITE(numout,*) 'ERROR with umask at the western border ji,jj,jk ', ji+nimpp-1,jj+njmpp-1,jk
1141                        WRITE(numout,*) '      masks: parent, child ', ptab(ji,jj,jk), umask(ji,jj,jk)
1142                        kindic_agr = kindic_agr + 1
1143                     ELSEIF (eastern_side) THEN
1144                        WRITE(numout,*) 'ERROR with umask at the eastern border ji,jj,jk ', ji+nimpp-1,jj+njmpp-1,jk
1145                        WRITE(numout,*) '      masks: parent, child ', ptab(ji,jj,jk), umask(ji,jj,jk)
1146                        kindic_agr = kindic_agr + 1
1147                     ENDIF
1148                  ENDIF
1149               END DO
1150            END DO
1151         END DO
1152         !
1153      ENDIF
1154      !
1155   END SUBROUTINE interpumsk
1156
1157
1158   SUBROUTINE interpvmsk( ptab, i1, i2, j1, j2, k1, k2, before, nb, ndir )
1159      !!----------------------------------------------------------------------
1160      !!                  ***  ROUTINE interpvmsk  ***
1161      !!---------------------------------------------------------------------- 
1162      INTEGER                              , INTENT(in   ) ::   i1,i2,j1,j2,k1,k2
1163      REAL(wp),DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) ::   ptab
1164      LOGICAL                              , INTENT(in   ) ::   before
1165      INTEGER                              , INTENT(in   ) :: nb , ndir
1166      !
1167      INTEGER ::   ji, jj, jk
1168      LOGICAL ::   northern_side, southern_side     
1169      !!---------------------------------------------------------------------- 
1170      !   
1171      IF( before ) THEN
1172         ptab(i1:i2,j1:j2,k1:k2) = vmask(i1:i2,j1:j2,k1:k2)
1173      ELSE
1174         southern_side = (nb == 2).AND.(ndir == 1)
1175         northern_side = (nb == 2).AND.(ndir == 2)
1176         DO jk = k1, k2
1177            DO jj = j1, j2
1178               DO ji = i1, i2
1179                   ! Velocity mask at boundary edge points:
1180                  IF (ABS(ptab(ji,jj,jk) - vmask(ji,jj,jk)) > 1.D-2) THEN
1181                     IF (southern_side) THEN
1182                        WRITE(numout,*) 'ERROR with vmask at the southern border ji,jj,jk ', ji+nimpp-1,jj+njmpp-1,jk
1183                        WRITE(numout,*) '      masks: parent, child ', ptab(ji,jj,jk), vmask(ji,jj,jk)
1184                        kindic_agr = kindic_agr + 1
1185                     ELSEIF (northern_side) THEN
1186                        WRITE(numout,*) 'ERROR with vmask at the northern border ji,jj,jk ', ji+nimpp-1,jj+njmpp-1,jk
1187                        WRITE(numout,*) '      masks: parent, child ', ptab(ji,jj,jk), vmask(ji,jj,jk)
1188                        kindic_agr = kindic_agr + 1
1189                     ENDIF
1190                  ENDIF
1191               END DO
1192            END DO
1193         END DO
1194         !
1195      ENDIF
1196      !
1197   END SUBROUTINE interpvmsk
1198
1199
1200   SUBROUTINE interpavm( ptab, i1, i2, j1, j2, k1, k2, m1, m2, before )
1201      !!----------------------------------------------------------------------
1202      !!                  ***  ROUTINE interavm  ***
1203      !!---------------------------------------------------------------------- 
1204      INTEGER                                    , INTENT(in   ) ::   i1, i2, j1, j2, k1, k2, m1, m2
1205      REAL(wp),DIMENSION(i1:i2,j1:j2,k1:k2,m1:m2), INTENT(inout) ::   ptab
1206      LOGICAL                                    , INTENT(in   ) ::   before
1207      REAL(wp), DIMENSION(k1:k2) :: tabin, h_in
1208      REAL(wp), DIMENSION(1:jpk) :: h_out
1209      INTEGER  :: N_in, N_out, ji, jj, jk
1210      !!---------------------------------------------------------------------- 
1211      !     
1212      IF (before) THEN         
1213         DO jk=k1,k2
1214            DO jj=j1,j2
1215              DO ji=i1,i2
1216                    ptab(ji,jj,jk,1) = avm_k(ji,jj,jk)
1217              END DO
1218           END DO
1219        END DO
1220#ifdef key_vertical         
1221        DO jk=k1,k2
1222           DO jj=j1,j2
1223              DO ji=i1,i2
1224                 ptab(ji,jj,jk,2) = wmask(ji,jj,jk) * e3w_n(ji,jj,jk) 
1225              END DO
1226           END DO
1227        END DO
1228#endif
1229      ELSE 
1230#ifdef key_vertical         
1231         avm_k(i1:i2,j1:j2,1:jpk) = 0.
1232         DO jj=j1,j2
1233            DO ji=i1,i2
1234               N_in = 0
1235               DO jk=k1,k2 !k2 = jpk of parent grid
1236                  IF (ptab(ji,jj,jk,2) == 0) EXIT
1237                  N_in = N_in + 1
1238                  tabin(jk) = ptab(ji,jj,jk,1)
1239                  h_in(N_in) = ptab(ji,jj,jk,2)
1240               END DO
1241               N_out = 0
1242               DO jk=1,jpk ! jpk of child grid
1243                  IF (wmask(ji,jj,jk) == 0) EXIT
1244                  N_out = N_out + 1
1245                  h_out(jk) = e3t_n(ji,jj,jk)
1246               ENDDO
1247               IF (N_in > 0) THEN
1248                  CALL reconstructandremap(tabin(1:N_in),h_in,avm_k(ji,jj,1:N_out),h_out,N_in,N_out,1)
1249               ENDIF
1250            ENDDO
1251         ENDDO
1252#else
1253         avm_k(i1:i2,j1:j2,k1:k2) = ptab (i1:i2,j1:j2,k1:k2,1)
1254#endif
1255      ENDIF
1256      !
1257   END SUBROUTINE interpavm
1258
1259#else
1260   !!----------------------------------------------------------------------
1261   !!   Empty module                                          no AGRIF zoom
1262   !!----------------------------------------------------------------------
1263CONTAINS
1264   SUBROUTINE Agrif_OCE_Interp_empty
1265      WRITE(*,*)  'agrif_oce_interp : You should not have seen this print! error?'
1266   END SUBROUTINE Agrif_OCE_Interp_empty
1267#endif
1268
1269   !!======================================================================
1270END MODULE agrif_oce_interp
Note: See TracBrowser for help on using the repository browser.