/[MITgcm]/MITgcm_contrib/gael/pkg/smooth/smooth_rhs.F
ViewVC logotype

Annotation of /MITgcm_contrib/gael/pkg/smooth/smooth_rhs.F

Parent Directory Parent Directory | Revision Log Revision Log | View Revision Graph Revision Graph


Revision 1.2 - (hide annotations) (download)
Fri Oct 16 03:36:34 2009 UTC (15 years, 9 months ago) by gforget
Branch: MAIN
CVS Tags: HEAD
Changes since 1.1: +5 -5 lines
bring pkg/smooth up to date

1 gforget 1.1 #include "CPP_OPTIONS.h"
2     #include "GAD_OPTIONS.h"
3    
4     C !INTERFACE: ==========================================================
5     SUBROUTINE smooth_rhs(fld_in,gt_in,myThid)
6    
7     C !DESCRIPTION:
8     C Calculates the tendency of a tracer due to advection and diffusion.
9     C It calculates the fluxes in each direction indepentently and then
10     C sets the tendency to the divergence of these fluxes. The advective
11     C fluxes are only calculated here when using the linear advection schemes
12     C otherwise only the diffusive and parameterized fluxes are calculated.
13     C
14     C Contributions to the flux are calculated and added:
15     C \begin{equation*}
16     C {\bf F} = {\bf F}_{adv} + {\bf F}_{diff} +{\bf F}_{GM} + {\bf F}_{KPP}
17     C \end{equation*}
18     C
19     C The tendency is the divergence of the fluxes:
20     C \begin{equation*}
21     C G_\theta = G_\theta + \nabla \cdot {\bf F}
22     C \end{equation*}
23     C
24     C The tendency is assumed to contain data on entry.
25    
26     C !USES: ===============================================================
27     IMPLICIT NONE
28     #include "SIZE.h"
29     #include "EEPARAMS.h"
30     #include "PARAMS.h"
31     #include "GRID.h"
32     #include "SURFACE.h"
33     #include "GAD.h"
34     #include "GMREDI.h"
35    
36     #ifdef ALLOW_AUTODIFF_TAMC
37     #include "tamc.h"
38     #include "tamc_keys.h"
39     #endif /* ALLOW_AUTODIFF_TAMC */
40    
41     #include "smooth.h"
42    
43     C !INPUT PARAMETERS: ===================================================
44    
45     c ... xA, yA, maskUp => in commons ?determined at init?
46     c dTdy...: peuvent devenir (bi,bj)
47    
48     INTEGER bi,bj,iMin,iMax,jMin,jMax
49     _RS xA (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
50     _RS yA (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
51     _RS maskUp(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
52     _RL dTdz (nSx,nSy)
53     _RL dTdx (nSx,nSy)
54     _RL dTdy (nSx,nSy)
55     INTEGER myThid
56     INTEGER i,j,k
57    
58     _RL fZon (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nR,nSx,nSy)
59     _RL fMer (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nR,nSx,nSy)
60     _RL fVerT (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nR,nSx,nSy)
61     _RL df (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
62     _RL fld_in(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr,nSx,nSy)
63     _RL gt_in(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr,nSx,nSy)
64     integer wc01_di, wc01_dj, wc01_dk, ii, jj, kk
65     _RL wc01_tmp(2)
66    
67     DO bj=myByLo(myThid),myByHi(myThid)
68     DO bi=myBxLo(myThid),myBxHi(myThid)
69    
70    
71     c 1rst k loop: initialization
72     DO k=1,Nr
73     DO j=1-OLy,sNy+OLy
74     DO i=1-OLx,sNx+OLx
75     fZon(i,j,k,bi,bj) = 0. _d 0
76     fMer(i,j,k,bi,bj) = 0. _d 0
77     fVerT(i,j,k,bi,bj) = 0. _d 0
78     gt_in(i,j,k,bi,bj) = 0. _d 0
79     ENDDO
80     ENDDO
81     ENDDO
82    
83 gforget 1.2 cnot needed et !!!!!!!!bibj _EXCH_XYZ_RL( gt_in , myThid )
84 gforget 1.1
85     iMin = 1-OLx+1
86     iMax = sNx+OLx-1
87     jMin = 1-OLy+1
88     jMax = sNy+OLy-1
89    
90     c 2nd k loop: flux computation
91     DO k=1,Nr
92    
93     DO j=1-OLy,sNy+OLy
94     DO i=1-OLx,sNx+OLx
95     df(i,j,bi,bj) = 0. _d 0
96     xA(i,j,bi,bj) = _dyG(i,j,bi,bj)
97     & *drF(k)*_hFacW(i,j,k,bi,bj)
98     yA(i,j,bi,bj) = _dxG(i,j,bi,bj)
99     & *drF(k)*_hFacS(i,j,k,bi,bj)
100     IF (K .EQ. 1) THEN
101     maskUp(i,j,bi,bj) = 0.
102     ELSE
103     maskUp(i,j,bi,bj) =
104     & maskC(i,j,k-1,bi,bj)*maskC(i,j,k,bi,bj)
105     ENDIF
106     ENDDO
107     ENDDO
108    
109    
110     c ///gmredi_xtr///
111    
112     DO j=jMin,jMax
113     DO i=iMin,iMax
114     df(i,j,bi,bj) = df(i,j,bi,bj)
115     & -xA(i,j,bi,bj)
116     & *wc01_Kux(i,j,k,bi,bj)
117     & *recip_dxC(i,j,bi,bj)
118     & *(fld_in(i,j,k,bi,bj)-fld_in(i-1,j,k,bi,bj))
119     ENDDO
120     ENDDO
121    
122     DO j=jMin,jMax
123     DO i=iMin,iMax
124     dTdz(bi,bj) = op5*(
125     & +op5*recip_drC(k)*
126     & ( maskC(i-1,j,k,bi,bj)*
127     & (fld_in(i-1,j, MAX(k-1,1) ,bi,bj)-fld_in(i-1,j,k,bi,bj))
128     & +maskC( i ,j,k,bi,bj)*
129     & (fld_in( i ,j, MAX(k-1,1) ,bi,bj)-fld_in( i ,j,k,bi,bj))
130     & )
131     & +op5*recip_drC(MIN(k+1,Nr))*
132     & ( maskC(i-1,j,MIN(k+1,Nr),bi,bj)*
133     & (fld_in(i-1,j,k,bi,bj)-fld_in(i-1,j,MIN(k+1,Nr),bi,bj))
134     & +maskC( i ,j,MIN(k+1,Nr),bi,bj)*
135     & (fld_in( i ,j,k,bi,bj)-fld_in( i ,j,MIN(k+1,Nr),bi,bj))
136     & ) )
137     df(i,j,bi,bj) = df(i,j,bi,bj)
138     & - xA(i,j,bi,bj)*wc01_Kuz(i,j,k,bi,bj)*dTdz(bi,bj)
139     ENDDO
140     ENDDO
141    
142     DO j=jMin,jMax
143     DO i=iMin,iMax
144     dTdy(bi,bj) = op5*(
145     & +op5*(maskS(i,j,k,bi,bj)
146     & *recip_dyC(i,j,bi,bj)*
147     & (fld_in(i,j,k,bi,bj)-fld_in(i,j-1,k,bi,bj))
148     & +maskS(i,j+1,k,bi,bj)
149     & *recip_dyC(i,j+1,bi,bj)*
150     & (fld_in(i,j+1,k,bi,bj)-fld_in(i,j,k,bi,bj)))
151     & +op5*(maskS(i-1,j,k,bi,bj)
152     & *recip_dyC(i,j,bi,bj)*
153     & (fld_in(i-1,j,k,bi,bj)-fld_in(i-1,j-1,k,bi,bj))
154     & +maskS(i-1,j+1,k,bi,bj)
155     & *recip_dyC(i,j+1,bi,bj)*
156     & (fld_in(i-1,j+1,k,bi,bj)-fld_in(i-1,j,k,bi,bj)))
157     & )
158     df(i,j,bi,bj) = df(i,j,bi,bj)
159     & - xA(i,j,bi,bj)*wc01_Kuy(i,j,k,bi,bj)*dTdy(bi,bj)
160     ENDDO
161     ENDDO
162    
163    
164     c /// end for x ///
165    
166     DO j=jMin,jMax
167     DO i=iMin,iMax
168     fZon(i,j,k,bi,bj) = fZon(i,j,k,bi,bj) + df(i,j,bi,bj)
169     ENDDO
170     ENDDO
171    
172     DO j=jMin,jMax
173     DO i=iMin,iMax
174     df(i,j,bi,bj) = 0.
175     ENDDO
176     ENDDO
177    
178     c ///gmredi_ytr///
179    
180     DO j=jMin,jMax
181     DO i=iMin,iMax
182     df(i,j,bi,bj) = df(i,j,bi,bj)
183     & -yA(i,j,bi,bj)
184     & *wc01_Kvy(i,j,k,bi,bj)
185     & *recip_dyC(i,j,bi,bj)
186     & *(fld_in(i,j,k,bi,bj)-fld_in(i,j-1,k,bi,bj))
187     ENDDO
188     ENDDO
189    
190     DO j=jMin,jMax
191     DO i=iMin,iMax
192     dTdz(bi,bj) = op5*(
193     & +op5*recip_drC(k)*
194     & ( maskC(i,j-1,k,bi,bj)*
195     & (fld_in(i,j-1,MAX(k-1,1),bi,bj)-fld_in(i,j-1,k,bi,bj))
196     & +maskC(i, j ,k,bi,bj)*
197     & (fld_in(i, j ,MAX(k-1,1),bi,bj)-fld_in(i, j ,k,bi,bj))
198     & )
199     & +op5*recip_drC(MIN(k+1,Nr))*
200     & ( maskC(i,j-1,MIN(k+1,Nr),bi,bj)*
201     & (fld_in(i,j-1,k,bi,bj)-fld_in(i,j-1,MIN(k+1,Nr),bi,bj))
202     & +maskC(i, j ,MIN(k+1,Nr),bi,bj)*
203     & (fld_in(i, j ,k,bi,bj)-fld_in(i, j ,MIN(k+1,Nr),bi,bj))
204     & ) )
205     df(i,j,bi,bj) = df(i,j,bi,bj)
206     & - yA(i,j,bi,bj)*wc01_Kvz(i,j,k,bi,bj)*dTdz(bi,bj)
207     ENDDO
208     ENDDO
209    
210     DO j=jMin,jMax
211     DO i=iMin,iMax
212     dTdx(bi,bj) = op5*(
213     & +op5*(maskW(i+1,j,k,bi,bj)
214     & *recip_dxC(i+1,j,bi,bj)*
215     & (fld_in(i+1,j,k,bi,bj)-fld_in(i,j,k,bi,bj))
216     & +maskW(i,j,k,bi,bj)
217     & *recip_dxC(i,j,bi,bj)*
218     & (fld_in(i,j,k,bi,bj)-fld_in(i-1,j,k,bi,bj)))
219     & +op5*(maskW(i+1,j-1,k,bi,bj)
220     & *recip_dxC(i+1,j,bi,bj)*
221     & (fld_in(i+1,j-1,k,bi,bj)-fld_in(i,j-1,k,bi,bj))
222     & +maskW(i,j-1,k,bi,bj)
223     & *recip_dxC(i,j,bi,bj)*
224     & (fld_in(i,j-1,k,bi,bj)-fld_in(i-1,j-1,k,bi,bj)))
225     & )
226     df(i,j,bi,bj) = df(i,j,bi,bj)
227     & - yA(i,j,bi,bj)*wc01_Kvx(i,j,k,bi,bj)*dTdx(bi,bj)
228     ENDDO
229     ENDDO
230    
231     c /// end for y ///
232    
233     DO j=jMin,jMax
234     DO i=iMin,iMax
235     fMer(i,j,k,bi,bj) = fMer(i,j,k,bi,bj) + df(i,j,bi,bj)
236     ENDDO
237     ENDDO
238    
239     DO j=jMin,jMax
240     DO i=iMin,iMax
241     df(i,j,bi,bj) = 0.
242     ENDDO
243     ENDDO
244    
245     c /// GAD_DIFF_R ///
246    
247     if (.NOT. wc01_implicitDiffusion ) then
248    
249     IF (k.gt.1) then
250     DO j=jMin,jMax
251     DO i=iMin,iMax
252     df(i,j,bi,bj) =
253     & -_rA(i,j,bi,bj)
254     & *KappaRwc01(i,j,k,bi,bj)*recip_drC(k)
255     & *(fld_in(i,j,k,bi,bj)
256     & -fld_in(i,j,k-1,bi,bj))*rkSign
257     ENDDO
258     ENDDO
259     ENDIF
260    
261     endif
262    
263     c ///gmredi rtrans///
264    
265     IF (K.GT.1) THEN
266     DO j=jMin,jMax
267     DO i=iMin,iMax
268     dTdx(bi,bj) = op5*(
269     & +op5*(maskW(i+1,j,k,bi,bj)
270     & *recip_dxC(i+1,j,bi,bj)*
271     & (fld_in(i+1,j,k,bi,bj)-fld_in(i,j,k,bi,bj))
272     & +maskW(i,j,k,bi,bj)
273     & *recip_dxC(i,j,bi,bj)*
274     & (fld_in(i,j,k,bi,bj)-fld_in(i-1,j,k,bi,bj)))
275     & +op5*(maskW(i+1,j,k-1,bi,bj)
276     & *recip_dxC(i+1,j,bi,bj)*
277     & (fld_in(i+1,j,k-1,bi,bj)-fld_in(i,j,k-1,bi,bj))
278     & +maskW(i,j,k-1,bi,bj)
279     & *recip_dxC(i,j,bi,bj)*
280     & (fld_in(i,j,k-1,bi,bj)-fld_in(i-1,j,k-1,bi,bj)))
281     & )
282    
283     dTdy(bi,bj) = op5*(
284     & +op5*(maskS(i,j,k,bi,bj)
285     & *recip_dyC(i,j,bi,bj)*
286     & (fld_in(i,j,k,bi,bj)-fld_in(i,j-1,k,bi,bj))
287     & +maskS(i,j+1,k,bi,bj)
288     & *recip_dyC(i,j+1,bi,bj)*
289     & (fld_in(i,j+1,k,bi,bj)-fld_in(i,j,k,bi,bj)))
290     & +op5*(maskS(i,j,k-1,bi,bj)
291     & *recip_dyC(i,j,bi,bj)*
292     & (fld_in(i,j,k-1,bi,bj)-fld_in(i,j-1,k-1,bi,bj))
293     & +maskS(i,j+1,k-1,bi,bj)
294     & *recip_dyC(i,j+1,bi,bj)*
295     & (fld_in(i,j+1,k-1,bi,bj)-fld_in(i,j,k-1,bi,bj)))
296     & )
297    
298     df(i,j,bi,bj) = df(i,j,bi,bj)
299     & - rA(i,j,bi,bj)
300     & *( wc01_Kwx(i,j,k,bi,bj)*dTdx(bi,bj)
301     & +wc01_Kwy(i,j,k,bi,bj)*dTdy(bi,bj) )
302    
303     ENDDO
304     ENDDO
305    
306     ENDIF
307    
308    
309     c /// end for r ///
310    
311     IF (K.GT.1) THEN
312     DO j=jMin,jMax
313     DO i=iMin,iMax
314     fVerT(i,j,k-1,bi,bj) = fVerT(i,j,k-1,bi,bj) +
315     & df(i,j,bi,bj)*maskUp(i,j,bi,bj)
316     ENDDO
317     ENDDO
318     ENDIF
319    
320     DO j=jMin,jMax
321     DO i=iMin,iMax
322     df(i,j,bi,bj) = 0.
323     ENDDO
324     ENDDO
325    
326     ENDDO
327    
328     ENDDO
329     ENDDO
330    
331     c these exchanges are crucial... two bibj loops
332 gforget 1.2 _EXCH_XYZ_RL ( fZon , myThid )
333     _EXCH_XYZ_RL ( fMer , myThid )
334     _EXCH_XYZ_RL ( fVerT, myThid )
335 gforget 1.1
336     DO bj=myByLo(myThid),myByHi(myThid)
337     DO bi=myBxLo(myThid),myBxHi(myThid)
338     c 3rd k loop: Divergence of fluxes
339     DO k=1,Nr
340     IF (K.GT.1) THEN
341     DO j=jMin,jMax
342     DO i=iMin,iMax
343     gt_in(i,j,k,bi,bj)=gt_in(i,j,k,bi,bj)
344     & -_recip_hFacC(i,j,k,bi,bj)*recip_drF(k)*recip_rA(i,j,bi,bj)
345     & *( (fZon(i+1,j,k,bi,bj)-fZon(i,j,k,bi,bj))
346     & +(fMer(i,j+1,k,bi,bj)-fMer(i,j,k,bi,bj))
347     & +(fVerT(i,j,k,bi,bj)-fVerT(i,j,k-1,bi,bj))*rkSign
348     & )
349     ENDDO
350     ENDDO
351     ELSE
352     DO j=jMin,jMax
353     DO i=iMin,iMax
354     gt_in(i,j,k,bi,bj)=gt_in(i,j,k,bi,bj)
355     & -_recip_hFacC(i,j,k,bi,bj)*recip_drF(k)*recip_rA(i,j,bi,bj)
356     & *( (fZon(i+1,j,k,bi,bj)-fZon(i,j,k,bi,bj))
357     & +(fMer(i,j+1,k,bi,bj)-fMer(i,j,k,bi,bj))
358     & +(fVerT(i,j,k,bi,bj))*rkSign
359     & )
360     ENDDO
361     ENDDO
362     ENDIF
363     ENDDO
364     ENDDO
365     ENDDO
366    
367 gforget 1.2 _EXCH_XYZ_RL ( gt_in , myThid )
368 gforget 1.1
369     END
370    

  ViewVC Help
Powered by ViewVC 1.1.22