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

Contents of /MITgcm/pkg/smooth/smooth_rhs.F

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


Revision 1.3 - (show annotations) (download)
Tue Sep 4 14:37:18 2012 UTC (11 years, 9 months ago) by gforget
Branch: MAIN
CVS Tags: checkpoint64y, checkpoint64x, checkpoint64z, checkpoint64q, checkpoint64p, checkpoint64s, checkpoint64r, checkpoint64u, checkpoint64t, checkpoint64w, checkpoint64v, checkpoint64i, checkpoint64h, checkpoint64k, checkpoint64j, checkpoint64m, checkpoint64l, checkpoint64o, checkpoint64n, checkpoint64a, checkpoint64c, checkpoint64b, checkpoint64e, checkpoint64d, checkpoint64g, checkpoint64f, checkpoint63s, checkpoint64, checkpoint65, checkpoint65h, checkpoint65i, checkpoint65b, checkpoint65c, checkpoint65a, checkpoint65f, checkpoint65g, checkpoint65d, checkpoint65e
Changes since 1.2: +204 -210 lines
- remove in-necessary includes.
- remove ALLOW_SMOOTH*D and ALLOW_SMOOTH_CORREL*D brakets.
  Those CPP options were never defined, and not necessary.

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

  ViewVC Help
Powered by ViewVC 1.1.22