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 NEMO/branches/2020/dev_r13312_AGRIF-03-04_jchanut_vinterp_tstep/src/NST – NEMO

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

Last change on this file since 13521 was 13352, checked in by jchanut, 4 years ago

#2222, #2129: correct depths used in linear vertical interpolation. NB: these are also used for initial state interpolation/extrapolation with identical vertical grids

  • Property svn:keywords set to Id
File size: 7.1 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   USE vremap
21   !
22   USE lib_mpp     ! MPP library
23
24   IMPLICIT NONE
25   PRIVATE
26
27   PUBLIC Agrif_trc, interptrn
28
29  !!----------------------------------------------------------------------
30   !! NEMO/NST 4.0 , NEMO Consortium (2018)
31   !! $Id$
32   !! Software governed by the CeCILL license (see ./LICENSE)
33   !!----------------------------------------------------------------------
34CONTAINS
35
36   SUBROUTINE Agrif_trc
37      !!----------------------------------------------------------------------
38      !!                   ***  ROUTINE Agrif_trc  ***
39      !!----------------------------------------------------------------------
40      !
41      IF( Agrif_Root() )   RETURN
42      !
43      Agrif_SpecialValue    = 0._wp
44      Agrif_UseSpecialValue = .TRUE.
45      l_vremap = ln_vert_remap
46      !
47      CALL Agrif_Bc_variable( trn_id, procname=interptrn )
48      !
49      Agrif_UseSpecialValue = .FALSE.
50      l_vremap              = .FALSE.
51      !
52   END SUBROUTINE Agrif_trc
53
54   SUBROUTINE interptrn( ptab, i1, i2, j1, j2, k1, k2, n1, n2, before )
55      !!----------------------------------------------------------------------
56      !!                  *** ROUTINE interptrn ***
57      !!----------------------------------------------------------------------
58      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) ::   ptab
59      INTEGER                                     , INTENT(in   ) ::   i1, i2, j1, j2, k1, k2, n1, n2
60      LOGICAL                                     , INTENT(in   ) ::   before
61      !
62      INTEGER  ::   ji, jj, jk, jn  ! dummy loop indices
63      INTEGER  ::   N_in, N_out
64      INTEGER  :: item
65      ! vertical interpolation:
66      REAL(wp) :: zhtot
67      REAL(wp), DIMENSION(k1:k2,1:jptra) :: tabin
68      REAL(wp), DIMENSION(k1:k2) :: h_in, z_in
69      REAL(wp), DIMENSION(1:jpk) :: h_out, z_out
70      !!----------------------------------------------------------------------
71
72      IF( before ) THEN
73
74         item = Kmm_a
75         IF( l_ini_child )   Kmm_a = Kbb_a 
76
77         DO jn = 1,jptra
78            DO jk=k1,k2
79               DO jj=j1,j2
80                 DO ji=i1,i2
81                       ptab(ji,jj,jk,jn) = tr(ji,jj,jk,jn,Kmm_a)
82                 END DO
83              END DO
84           END DO
85         END DO
86
87         IF( l_vremap .OR. l_ini_child) THEN
88            ! Interpolate thicknesses
89            ! Warning: these are masked, hence extrapolated prior interpolation.
90            DO jk=k1,k2
91               DO jj=j1,j2
92                  DO ji=i1,i2
93                      ptab(ji,jj,jk,jptra+1) = tmask(ji,jj,jk) * e3t(ji,jj,jk,Kmm_a)
94
95                  END DO
96               END DO
97            END DO
98
99            ! Extrapolate thicknesses in partial bottom cells:
100            ! Set them to Agrif_SpecialValue (0.). Correct bottom thicknesses are retrieved later on
101            IF (ln_zps) THEN
102               DO jj=j1,j2
103                  DO ji=i1,i2
104                      jk = mbkt(ji,jj)
105                      ptab(ji,jj,jk,jptra+1) = 0._wp
106                  END DO
107               END DO           
108            END IF
109       
110            ! Save ssh at last level:
111            IF (.NOT.ln_linssh) THEN
112               ptab(i1:i2,j1:j2,k2,jptra+1) = ssh(i1:i2,j1:j2,Kmm_a)*tmask(i1:i2,j1:j2,1) 
113            ELSE
114               ptab(i1:i2,j1:j2,k2,jptra+1) = 0._wp
115            END IF     
116         ENDIF
117         Kmm_a = item
118
119      ELSE
120         item = Krhs_a
121         IF( l_ini_child )   Krhs_a = Kbb_a 
122
123         IF( l_vremap .OR. l_ini_child ) THEN
124            IF (ln_linssh) ptab(i1:i2,j1:j2,k2,n2) = 0._wp 
125               
126            DO jj=j1,j2
127               DO ji=i1,i2
128                  tr(ji,jj,:,:,Krhs_a) = 0.                 
129                  N_in = mbkt_parent(ji,jj)
130                  zhtot = 0._wp
131                  DO jk=1,N_in !k2 = jpk of parent grid
132                     IF (jk==N_in) THEN
133                        h_in(jk) = ht0_parent(ji,jj) + ptab(ji,jj,k2,n2) - zhtot
134                     ELSE
135                        h_in(jk) = ptab(ji,jj,jk,n2)
136                     ENDIF
137                     zhtot = zhtot + h_in(jk)
138                     tabin(jk,:) = ptab(ji,jj,jk,n1:n2-1)
139                  END DO
140                  z_in(1) = 0.5_wp * h_in(1) - zhtot + ht0_parent(ji,jj)
141                  DO jk=2,N_in
142                     z_in(jk) = z_in(jk-1) + 0.5_wp * (h_in(jk-1)+h_in(jk))
143                  END DO
144
145                  N_out = 0
146                  DO jk=1,jpk ! jpk of child grid
147                     IF (tmask(ji,jj,jk) == 0._wp) EXIT
148                     N_out = N_out + 1
149                     h_out(jk) = e3t(ji,jj,jk,Krhs_a)
150                  END DO
151
152                  z_out(1) = 0.5_wp * h_out(1) - SUM(h_out(1:N_out)) + ht_0(ji,jj)
153                  DO jk=2,N_out
154                     z_out(jk) = z_out(jk-1) + 0.5_wp * (h_out(jk-1)+h_out(jk))
155                  END DO
156
157                  IF (N_in*N_out > 0) THEN
158                     IF( l_ini_child ) THEN
159                        CALL remap_linear(tabin(1:N_in,1:jptra),z_in(1:N_in),tr(ji,jj,1:N_out,1:jptra,Krhs_a),          &
160                                      &   z_out(1:N_out),N_in,N_out,jptra) 
161                     ELSE
162                        CALL reconstructandremap(tabin(1:N_in,1:jptra),h_in(1:N_in),tr(ji,jj,1:N_out,1:jptra,Krhs_a),   &
163                                      &   h_out(1:N_out),N_in,N_out,jptra) 
164                     ENDIF
165                  ENDIF
166               END DO
167            END DO
168            Krhs_a = item
169 
170         ELSE
171         
172            DO jn=1, jptra
173                tr(i1:i2,j1:j2,1:jpk,jn,Krhs_a)=ptab(i1:i2,j1:j2,1:jpk,jn)*tmask(i1:i2,j1:j2,1:jpk) 
174            END DO
175         ENDIF
176
177      ENDIF
178      !
179   END SUBROUTINE interptrn
180
181#else
182   !!----------------------------------------------------------------------
183   !!   Empty module                                           no TOP AGRIF
184   !!----------------------------------------------------------------------
185CONTAINS
186   SUBROUTINE Agrif_TOP_Interp_empty
187      !!---------------------------------------------
188      !!   *** ROUTINE agrif_Top_Interp_empty ***
189      !!---------------------------------------------
190      WRITE(*,*)  'agrif_top_interp : You should not have seen this print! error?'
191   END SUBROUTINE Agrif_TOP_Interp_empty
192#endif
193
194   !!======================================================================
195END MODULE agrif_top_interp
Note: See TracBrowser for help on using the repository browser.