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_top_interp.F90 in branches/UKMO/dev_merge_2017_CICE_interface/NEMOGCM/NEMO/NST_SRC – NEMO

source: branches/UKMO/dev_merge_2017_CICE_interface/NEMOGCM/NEMO/NST_SRC/agrif_top_interp.F90 @ 9499

Last change on this file since 9499 was 9499, checked in by davestorkey, 6 years ago

branches/UKMO/dev_merge_2017_CICE_interface : clear SVN keywords.

File size: 11.5 KB
Line 
1MODULE agrif_top_interp
2   !!======================================================================
3   !!                   ***  MODULE  agrif_top_interp  ***
4   !! AGRIF: interpolation package for TOP
5   !!======================================================================
6   !! History :  2.0  !  ???
7   !!----------------------------------------------------------------------
8#if defined key_agrif && defined key_top
9   !!----------------------------------------------------------------------
10   !!   'key_agrif'                                              AGRIF zoom
11   !!   'key_top'                                           on-line tracers
12   !!----------------------------------------------------------------------
13   USE par_oce
14   USE oce
15   USE dom_oce     
16   USE agrif_oce
17   USE agrif_top_sponge
18   USE par_trc
19   USE trc
20   !
21   USE lib_mpp     ! MPP library
22
23   IMPLICIT NONE
24   PRIVATE
25
26   PUBLIC Agrif_trc, interptrn
27
28  !!----------------------------------------------------------------------
29   !! NEMO/NST 4.0 , NEMO Consortium (2017)
30   !! $Id$
31   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
32   !!----------------------------------------------------------------------
33CONTAINS
34
35   SUBROUTINE Agrif_trc
36      !!----------------------------------------------------------------------
37      !!                   ***  ROUTINE Agrif_trc  ***
38      !!----------------------------------------------------------------------
39      !
40      IF( Agrif_Root() )   RETURN
41      !
42      Agrif_SpecialValue    = 0._wp
43      Agrif_UseSpecialValue = .TRUE.
44      !
45      CALL Agrif_Bc_variable( trn_id, procname=interptrn )
46      Agrif_UseSpecialValue = .FALSE.
47      !
48   END SUBROUTINE Agrif_trc
49
50
51   SUBROUTINE interptrn( ptab, i1, i2, j1, j2, k1, k2, n1, n2, before, nb, ndir )
52      !!----------------------------------------------------------------------
53      !!                   *** ROUTINE interptrn ***
54      !!----------------------------------------------------------------------
55      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) ::   ptab
56      INTEGER                                     , INTENT(in   ) ::   i1, i2, j1, j2, k1, k2, n1, n2
57      LOGICAL                                     , INTENT(in   ) ::   before
58      INTEGER                                     , INTENT(in   ) ::   nb , ndir
59      !!
60      INTEGER ::   ji, jj, jk, jn   ! dummy loop indices
61      INTEGER ::   imin, imax, jmin, jmax
62      LOGICAL ::   ll_west, ll_east, ll_north, ll_south
63      REAL(wp) ::   zrhox, z1, z2, z3, z4, z5, z6, z7
64      !!----------------------------------------------------------------------
65      !
66      INTEGER  ::   ji, jj, jk, jn   ! dummy loop indices
67      INTEGER  ::   imin, imax, jmin, jmax, N_in, N_out
68      REAL(wp) ::   zrhox , zalpha1, zalpha2, zalpha3
69      REAL(wp) ::   zalpha4, zalpha5, zalpha6, zalpha7
70      LOGICAL :: western_side, eastern_side,northern_side,southern_side
71      ! vertical interpolation:
72      REAL(wp), DIMENSION(i1:i2,j1:j2,1:jpk,n1:n2) :: ptab_child
73      REAL(wp), DIMENSION(k1:k2,n1:n2-1) :: tabin
74      REAL(wp), DIMENSION(k1:k2) :: h_in
75      REAL(wp), DIMENSION(1:jpk) :: h_out(1:jpk)
76      REAL(wp) :: h_diff, zrhoxy
77
78      zrhoxy = Agrif_rhox()*Agrif_rhoy()
79      IF (before) THEN         
80         DO jn = 1,jpts
81            DO jk=k1,k2
82               DO jj=j1,j2
83                 DO ji=i1,i2
84                       ptab(ji,jj,jk,jn) = trn(ji,jj,jk,jn)
85                 END DO
86               END DO
87            END DO
88         END DO
89# if defined key_vertical
90        DO jk=k1,k2
91           DO jj=j1,j2
92              DO ji=i1,i2
93                 ptab(ji,jj,jk,jpts+1) = tmask(ji,jj,jk) * e3t_n(ji,jj,jk) 
94              END DO
95           END DO
96        END DO
97# endif
98
99      ELSE
100         !
101         western_side  = (nb == 1).AND.(ndir == 1)
102         eastern_side  = (nb == 1).AND.(ndir == 2)
103         southern_side = (nb == 2).AND.(ndir == 1)
104         northern_side = (nb == 2).AND.(ndir == 2)
105
106# if defined key_vertical             
107         DO jj=j1,j2
108            DO ji=i1,i2
109               iref = ji
110               jref = jj
111               if(western_side) iref=MAX(2,ji)
112               if(eastern_side) iref=MIN(nlci-1,ji)
113               if(southern_side) jref=MAX(2,jj)
114               if(northern_side) jref=MIN(nlcj-1,jj)
115               N_in = 0
116               DO jk=k1,k2 !k2 = jpk of parent grid
117                  IF (ptab(ji,jj,jk,n2) == 0) EXIT
118                  N_in = N_in + 1
119                  tabin(jk,:) = ptab(ji,jj,jk,n1:n2-1)
120                  h_in(N_in) = ptab(ji,jj,jk,n2)
121               END DO
122               N_out = 0
123               DO jk=1,jpk ! jpk of child grid
124                  IF (tmask(iref,jref,jk) == 0) EXIT
125                  N_out = N_out + 1
126                  h_out(jk) = e3t_n(iref,jref,jk)
127               ENDDO
128               IF (N_in > 0) THEN
129                  h_diff = sum(h_out(1:N_out))-sum(h_in(1:N_in))
130                  DO jn=1,jptra
131                     call reconstructandremap(tabin(1:N_in,jn),h_in,ptab_child(ji,jj,1:N_out,jn),h_out,N_in,N_out)
132                  ENDDO
133               ENDIF
134            ENDDO
135         ENDDO
136# else
137         ptab_child(i1:i2,j1:j2,1:jpk,1:jptra) = ptab(i1:i2,j1:j2,1:jpk,1:jptra)
138# endif
139
140         !
141         zrhox = Agrif_Rhox()
142         !
143         zalpha1 = ( zrhox - 1. ) * 0.5
144         zalpha2 = 1. - zalpha1
145         !
146         zalpha3 = ( zrhox - 1. ) / ( zrhox + 1. )
147         zalpha4 = 1. - zalpha3
148         !
149         zalpha6 = 2. * ( zrhox - 1. ) / ( zrhox + 1. )
150         zalpha7 =    - ( zrhox - 1. ) / ( zrhox + 3. )
151         zalpha5 = 1. - zalpha6 - zalpha7
152         !
153         imin = i1
154         imax = i2
155         jmin = j1
156         jmax = j2
157         !
158         ! Remove CORNERS
159         IF((nbondj == -1).OR.(nbondj == 2)) jmin = 3
160         IF((nbondj == +1).OR.(nbondj == 2)) jmax = nlcj-2
161         IF((nbondi == -1).OR.(nbondi == 2)) imin = 3
162         IF((nbondi == +1).OR.(nbondi == 2)) imax = nlci-2       
163         !
164         IF( eastern_side) THEN
165            DO jn = 1, jptra
166               tra(nlci,j1:j2,1:jpk,jn) = zalpha1 * ptab_child(nlci,j1:j2,1:jpk,jn) + zalpha2 * ptab_child(nlci-1,j1:j2,1:jpk,jn)
167               DO jk = 1, jpkm1
168                  DO jj = jmin,jmax
169                     IF( umask(nlci-2,jj,jk) == 0.e0 ) THEN
170                        tra(nlci-1,jj,jk,jn) = tra(nlci,jj,jk,jn) * tmask(nlci-1,jj,jk)
171                     ELSE
172                        tra(nlci-1,jj,jk,jn)=(zalpha4*tra(nlci,jj,jk,jn)+zalpha3*tra(nlci-2,jj,jk,jn))*tmask(nlci-1,jj,jk)
173                        IF( un(nlci-2,jj,jk) > 0.e0 ) THEN
174                           tra(nlci-1,jj,jk,jn)=( zalpha6*tra(nlci-2,jj,jk,jn)+zalpha5*tra(nlci,jj,jk,jn) & 
175                                 + zalpha7*tra(nlci-3,jj,jk,jn) ) * tmask(nlci-1,jj,jk)
176                        ENDIF
177                     END DO
178                  END DO
179               END DO
180            ENDDO
181         ENDIF
182         !
183         IF( northern_side ) THEN           
184            DO jn = 1, jptra
185               tra(i1:i2,nlcj,1:jpk,jn) = zalpha1 * ptab_child(i1:i2,nlcj,1:jpk,jn) + zalpha2 * ptab_child(i1:i2,nlcj-1,1:jpk,jn)
186               DO jk = 1, jpkm1
187                  DO ji = imin,imax
188                     IF( vmask(ji,nlcj-2,jk) == 0.e0 ) THEN
189                        tra(ji,nlcj-1,jk,jn) = tra(ji,nlcj,jk,jn) * tmask(ji,nlcj-1,jk)
190                     ELSE
191                        tra(ji,nlcj-1,jk,jn)=(zalpha4*tra(ji,nlcj,jk,jn)+zalpha3*tra(ji,nlcj-2,jk,jn))*tmask(ji,nlcj-1,jk)       
192                        IF (vn(ji,nlcj-2,jk) > 0.e0 ) THEN
193                           tra(ji,nlcj-1,jk,jn)=( zalpha6*tra(ji,nlcj-2,jk,jn)+zalpha5*tra(ji,nlcj,jk,jn)  &
194                                 + zalpha7*tra(ji,nlcj-3,jk,jn) ) * tmask(ji,nlcj-1,jk)
195                        ENDIF
196                     END DO
197                  END DO
198               END DO
199            ENDDO
200         ENDIF
201         !
202         IF( western_side) THEN           
203            DO jn = 1, jptra
204               tra(1,j1:j2,1:jpk,jn) = zalpha1 * ptab_child(1,j1:j2,1:jpk,jn) + zalpha2 * ptab_child(2,j1:j2,1:jpk,jn)
205               DO jk = 1, jpkm1
206                  DO jj = jmin,jmax
207                     IF( umask(2,jj,jk) == 0.e0 ) THEN
208                        tra(2,jj,jk,jn) = tra(1,jj,jk,jn) * tmask(2,jj,jk)
209                     ELSE
210                        tra(2,jj,jk,jn)=(zalpha4*tra(1,jj,jk,jn)+zalpha3*tra(3,jj,jk,jn))*tmask(2,jj,jk)       
211                        IF( un(2,jj,jk) < 0.e0 ) THEN
212                           tra(2,jj,jk,jn)=(zalpha6*tra(3,jj,jk,jn)+zalpha5*tra(1,jj,jk,jn)+zalpha7*tra(4,jj,jk,jn))*tmask(2,jj,jk)
213                        ENDIF
214                     END DO
215                  END DO
216               END DO
217            END DO
218         ENDIF
219         !
220         IF( southern_side ) THEN           
221            DO jn = 1, jptra
222               tra(i1:i2,1,1:jpk,jn) = zalpha1 * ptab_child(i1:i2,1,1:jpk,jn) + zalpha2 * ptab_child(i1:i2,2,1:jpk,jn)
223               DO jk=1,jpkm1
224                  DO ji=imin,imax
225                     IF( vmask(ji,2,jk) == 0.e0 ) THEN
226                        tra(ji,2,jk,jn)=tra(ji,1,jk,jn) * tmask(ji,2,jk)
227                     ELSE
228                        tra(ji,2,jk,jn)=(zalpha4*tra(ji,1,jk,jn)+zalpha3*tra(ji,3,jk,jn))*tmask(ji,2,jk)
229                        IF( vn(ji,2,jk) < 0.e0 ) THEN
230                           tra(ji,2,jk,jn)=(zalpha6*tra(ji,3,jk,jn)+zalpha5*tra(ji,1,jk,jn)+zalpha7*tra(ji,4,jk,jn))*tmask(ji,2,jk)
231                        ENDIF
232                     END DO
233                  END DO
234               END DO
235            ENDIF
236            !
237            ! Treatment of corners
238            IF( ll_east .AND.((nbondj == -1).OR.(nbondj == 2)) )   tra(nlci-1,   2  ,:,:) = ptab(nlci-1,   2  ,:,:)   ! East south
239            IF( ll_east .AND.((nbondj ==  1).OR.(nbondj == 2)) )   tra(nlci-1,nlcj-1,:,:) = ptab(nlci-1,nlcj-1,:,:)   ! East north
240            IF( ll_west .AND.((nbondj == -1).OR.(nbondj == 2)) )   tra(   2  ,   2  ,:,:) = ptab(   2  ,   2  ,:,:)   ! West south
241            IF( ll_west .AND.((nbondj ==  1).OR.(nbondj == 2)) )   tra(   2  ,nlcj-1,:,:) = ptab(   2  ,nlcj-1,:,:)   ! West north
242            !
243         ENDIF
244         !
245         ! Treatment of corners
246         !
247         ! East south
248         IF ((eastern_side).AND.((nbondj == -1).OR.(nbondj == 2))) THEN
249            tra(nlci-1,2,:,:) = ptab_child(nlci-1,2,:,:)
250         ENDIF
251         ! East north
252         IF ((eastern_side).AND.((nbondj == 1).OR.(nbondj == 2))) THEN
253            tra(nlci-1,nlcj-1,:,:) = ptab_child(nlci-1,nlcj-1,:,:)
254         ENDIF
255         ! West south
256         IF ((western_side).AND.((nbondj == -1).OR.(nbondj == 2))) THEN
257            tra(2,2,:,:) = ptab_child(2,2,:,:)
258         ENDIF
259         ! West north
260         IF ((western_side).AND.((nbondj == 1).OR.(nbondj == 2))) THEN
261            tra(2,nlcj-1,:,:) = ptab_child(2,nlcj-1,:,:)
262         ENDIF
263         !
264      ENDIF
265      !
266   END SUBROUTINE interptrn
267
268#else
269   !!----------------------------------------------------------------------
270   !!   Empty module                                           no TOP AGRIF
271   !!----------------------------------------------------------------------
272CONTAINS
273   SUBROUTINE Agrif_TOP_Interp_empty
274      !!---------------------------------------------
275      !!   *** ROUTINE agrif_Top_Interp_empty ***
276      !!---------------------------------------------
277      WRITE(*,*)  'agrif_top_interp : You should not have seen this print! error?'
278   END SUBROUTINE Agrif_TOP_Interp_empty
279#endif
280
281   !!======================================================================
282END MODULE agrif_top_interp
Note: See TracBrowser for help on using the repository browser.