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_r10973_AGRIF-01_jchanut_small_jpi_jpj/src/NST – NEMO

source: NEMO/branches/2019/dev_r10973_AGRIF-01_jchanut_small_jpi_jpj/src/NST/agrif_oce_interp.F90 @ 11226

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

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