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_update.F90 in NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/NST – NEMO

source: NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/NST/agrif_top_update.F90 @ 11036

Last change on this file since 11036 was 10989, checked in by acc, 5 years ago

2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps : Convert NST routines in preparation for getting AGRIF back up and running. AGRIF conv stage now works but requires some renaming of recently changes DIU modules (included in this commit). AGRIF compile and link stage not yet working (agrif routines need to be passed the time-level indices) but non-AGRIF SETTE tests are all OK
  • Property svn:keywords set to Id
File size: 9.2 KB
Line 
1#define TWO_WAY
2#undef DECAL_FEEDBACK
3
4MODULE agrif_top_update
5   !!======================================================================
6   !!                ***  MODULE agrif_top_update  ***
7   !! AGRIF :   update package for passive tracers (TOP)
8   !!======================================================================
9   !! History : 
10   !!----------------------------------------------------------------------
11#if defined key_agrif && defined key_top
12   !!----------------------------------------------------------------------
13   !!   'key_agrif'                                              AGRIF zoom
14   !!   'key_TOP'                                           on-line tracers
15   !!----------------------------------------------------------------------
16   USE par_oce
17   USE oce
18   USE dom_oce
19   USE agrif_oce
20   USE par_trc
21   USE trc
22
23   IMPLICIT NONE
24   PRIVATE
25
26   PUBLIC Agrif_Update_Trc
27
28   !!----------------------------------------------------------------------
29   !! NEMO/NST 4.0 , NEMO Consortium (2018)
30   !! $Id$
31   !! Software governed by the CeCILL license (see ./LICENSE)
32   !!----------------------------------------------------------------------
33CONTAINS
34
35   SUBROUTINE Agrif_Update_Trc( )
36      !!----------------------------------------------------------------------
37      !!                   *** ROUTINE Agrif_Update_Trc ***
38      !!----------------------------------------------------------------------
39      !
40      IF (Agrif_Root()) RETURN 
41      !
42#if defined TWO_WAY   
43      Agrif_UseSpecialValueInUpdate = .TRUE.
44      Agrif_SpecialValueFineGrid    = 0._wp
45      !
46# if ! defined DECAL_FEEDBACK
47      CALL Agrif_Update_Variable(trn_id, procname=updateTRC )
48!      CALL Agrif_Update_Variable( trn_id, locupdate=(/0,2/), procname=updateTRC )
49# else
50      CALL Agrif_Update_Variable(trn_id, locupdate=(/1,0/),procname=updateTRC )
51!      CALL Agrif_Update_Variable( trn_id, locupdate=(/1,2/), procname=updateTRC )
52# endif
53      !
54      Agrif_UseSpecialValueInUpdate = .FALSE.
55      !
56#endif
57      !
58   END SUBROUTINE Agrif_Update_Trc
59
60#ifdef key_vertical
61   SUBROUTINE updateTRC( tabres, i1, i2, j1, j2, k1, k2, n1, n2, before )
62      !!---------------------------------------------
63      !!           *** ROUTINE updateT ***
64      !!---------------------------------------------
65      INTEGER, INTENT(in) :: i1,i2,j1,j2,k1,k2,n1,n2
66      REAL(wp),DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) :: tabres
67      LOGICAL, INTENT(in) :: before
68      !!
69      INTEGER :: ji,jj,jk,jn
70      REAL(wp), DIMENSION(i1:i2,j1:j2,1:jpk,n1:n2) :: tabres_child
71      REAL(wp) :: h_in(k1:k2)
72      REAL(wp) :: h_out(1:jpk)
73      INTEGER  :: N_in, N_out
74      REAL(wp) :: h_diff
75      REAL(wp) :: zrho_xy
76      REAL(wp) :: tabin(k1:k2,n1:n2)
77      !!---------------------------------------------
78      !
79      IF (before) THEN
80         AGRIF_SpecialValue = -999._wp
81         zrho_xy = Agrif_rhox() * Agrif_rhoy() 
82         DO jn = n1,n2-1
83            DO jk=k1,k2
84               DO jj=j1,j2
85                  DO ji=i1,i2
86                     tabres(ji,jj,jk,jn) = (tr(ji,jj,jk,jn,Kmm) * e3t(ji,jj,jk,Kmm) ) &
87                                           * tmask(ji,jj,jk) + (tmask(ji,jj,jk)-1)*999._wp
88                  END DO
89               END DO
90            END DO
91         END DO
92         DO jk=k1,k2
93            DO jj=j1,j2
94               DO ji=i1,i2
95                  tabres(ji,jj,jk,n2) = tmask(ji,jj,jk) * e3t(ji,jj,jk,Kmm) &
96                                           + (tmask(ji,jj,jk)-1)*999._wp
97               END DO
98            END DO
99         END DO
100      ELSE
101         tabres_child(:,:,:,:) = 0.
102         AGRIF_SpecialValue = 0._wp
103         DO jj=j1,j2
104            DO ji=i1,i2
105               N_in = 0
106               DO jk=k1,k2 !k2 = jpk of child grid
107                  IF (tabres(ji,jj,jk,n2) == 0  ) EXIT
108                  N_in = N_in + 1
109                  tabin(jk,:) = tabres(ji,jj,jk,n1:n2-1)/tabres(ji,jj,jk,n2)
110                  h_in(N_in) = tabres(ji,jj,jk,n2)
111               ENDDO
112               N_out = 0
113               DO jk=1,jpk ! jpk of parent grid
114                  IF (tmask(ji,jj,jk) < -900) EXIT ! TODO: Will not work with ISF
115                  N_out = N_out + 1
116                  h_out(N_out) = e3t(ji,jj,jk,Kmm) !Parent grid scale factors. Could multiply by e1e2t here instead of division above
117               ENDDO
118               IF (N_in > 0) THEN !Remove this?
119                  h_diff = sum(h_out(1:N_out))-sum(h_in(1:N_in))
120                  IF (h_diff < -1.e-4) THEN
121                     print *,'CHECK YOUR bathy T points ...',ji,jj,h_diff,sum(h_in(1:N_in)),sum(h_out(1:N_out))
122                     print *,h_in(1:N_in)
123                     print *,h_out(1:N_out)
124                     STOP
125                  ENDIF
126                  DO jn=1,jptra
127                     CALL reconstructandremap(tabin(1:N_in,jn),h_in(1:N_in),tabres_child(ji,jj,1:N_out,jn),h_out(1:N_out),N_in,N_out)
128                  ENDDO
129               ENDIF
130            ENDDO
131         ENDDO
132
133         IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN
134            ! Add asselin part
135            DO jn = 1,jptra
136               DO jk=1,jpk
137                  DO jj=j1,j2
138                     DO ji=i1,i2
139                        IF( tabres_child(ji,jj,jk,jn) .NE. 0. ) THEN
140                           tr(ji,jj,jk,jn,Kbb) = ts(ji,jj,jk,jn,Kbb) & 
141                                 & + atfp * ( tabres_child(ji,jj,jk,jn) &
142                                 &          - tr(ji,jj,jk,jn,Kmm) ) * tmask(ji,jj,jk)
143                        ENDIF
144                     ENDDO
145                  ENDDO
146               ENDDO
147            ENDDO
148         ENDIF
149         DO jn = 1,jptra
150            DO jk=1,jpk
151               DO jj=j1,j2
152                  DO ji=i1,i2
153                     IF( tabres_child(ji,jj,jk,jn) .NE. 0. ) THEN
154                        tr(ji,jj,jk,jn,Kmm) = tabres_child(ji,jj,jk,jn) * tmask(ji,jj,jk)
155                     END IF
156                  END DO
157               END DO
158            END DO
159         END DO
160      ENDIF
161      !
162   END SUBROUTINE updateTRC
163
164
165#else
166   SUBROUTINE updateTRC( tabres, i1, i2, j1, j2, k1, k2, n1, n2, before )
167      !!----------------------------------------------------------------------
168      !!                      *** ROUTINE updateTRC ***
169      !!----------------------------------------------------------------------
170      INTEGER                                    , INTENT(in   ) ::   i1, i2, j1, j2, k1, k2, n1, n2
171      REAL(wp),DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) ::   tabres
172      LOGICAL                                    , INTENT(in   ) ::   before
173      !!
174      INTEGER :: ji,jj,jk,jn
175      REAL(wp) :: ztb, ztnu, ztno
176      !!----------------------------------------------------------------------
177      !
178      !
179      IF (before) THEN
180         DO jn = n1,n2
181            DO jk=k1,k2
182               DO jj=j1,j2
183                  DO ji=i1,i2
184!> jc tmp
185                     tabres(ji,jj,jk,jn) = tr(ji,jj,jk,jn,Kmm)  * e3t(ji,jj,jk,Kmm) / e3t_0(ji,jj,jk)
186!                     tabres(ji,jj,jk,jn) = tr(ji,jj,jk,jn,Kmm)  * e3t(ji,jj,jk,Kmm)
187!< jc tmp
188                  END DO
189               END DO
190            END DO
191         END DO
192      ELSE
193!> jc tmp
194         DO jn = n1,n2
195            tabres(i1:i2,j1:j2,k1:k2,jn) =  tabres(i1:i2,j1:j2,k1:k2,jn) * e3t_0(i1:i2,j1:j2,k1:k2) &
196                                         & * tmask(i1:i2,j1:j2,k1:k2)
197         ENDDO
198!< jc tmp
199         IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN
200            ! Add asselin part
201            DO jn = n1,n2
202               DO jk=k1,k2
203                  DO jj=j1,j2
204                     DO ji=i1,i2
205                        IF( tabres(ji,jj,jk,jn) .NE. 0. ) THEN
206                           ztb  = tr(ji,jj,jk,jn,Kbb) * e3t(ji,jj,jk,Kbb) ! fse3t_b prior update should be used
207                           ztnu = tabres(ji,jj,jk,jn)
208                           ztno = tr(ji,jj,jk,jn,Kmm) * e3t(ji,jj,jk,Krhs)
209                           tr(ji,jj,jk,jn,Kbb) = ( ztb + atfp * ( ztnu - ztno) )  & 
210                                     &        * tmask(ji,jj,jk) / e3t(ji,jj,jk,Kbb)
211                        ENDIF
212                     ENDDO
213                  ENDDO
214               ENDDO
215            ENDDO
216         ENDIF
217         DO jn = n1,n2
218            DO jk=k1,k2
219               DO jj=j1,j2
220                  DO ji=i1,i2
221                     IF( tabres(ji,jj,jk,jn) .NE. 0. ) THEN
222                        tr(ji,jj,jk,jn,Kmm) = tabres(ji,jj,jk,jn) / e3t(ji,jj,jk,Kmm)
223                     END IF
224                  END DO
225               END DO
226            END DO
227         END DO
228         !
229         IF  ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN
230            tr(i1:i2,j1:j2,k1:k2,n1:n2,Kbb)  = tr(i1:i2,j1:j2,k1:k2,n1:n2,Kmm)
231         ENDIF
232         !
233      ENDIF
234      !
235   END SUBROUTINE updateTRC
236#endif
237
238#else
239   !!----------------------------------------------------------------------
240   !!   Empty module                                           no TOP AGRIF
241   !!----------------------------------------------------------------------
242CONTAINS
243   SUBROUTINE agrif_top_update_empty
244      WRITE(*,*)  'agrif_top_update : You should not have seen this print! error?'
245   END SUBROUTINE agrif_top_update_empty
246#endif
247
248   !!======================================================================
249END MODULE agrif_top_update
Note: See TracBrowser for help on using the repository browser.