/[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.1 - (show annotations) (download)
Mon Feb 15 23:46:04 2010 UTC (14 years, 4 months ago) by gforget
Branch: MAIN
CVS Tags: checkpoint63p, checkpoint63q, checkpoint63r, checkpoint63l, checkpoint63m, checkpoint63n, checkpoint63o, checkpoint63h, checkpoint63i, checkpoint63j, checkpoint63k, checkpoint63d, checkpoint63e, checkpoint63f, checkpoint63g, checkpoint63a, checkpoint63b, checkpoint63c, checkpoint63, checkpoint62c, checkpoint62g, checkpoint62f, checkpoint62e, checkpoint62d, checkpoint62k, checkpoint62j, checkpoint62i, checkpoint62h, checkpoint62o, checkpoint62n, checkpoint62m, checkpoint62l, checkpoint62s, checkpoint62r, checkpoint62q, checkpoint62p, checkpoint62w, checkpoint62v, checkpoint62u, checkpoint62t, checkpoint62z, checkpoint62y, checkpoint62x
smooth package -- going to the main repository

1 C $Header: /u/gcmpack/MITgcm_contrib/gael/pkg/smooth2/smooth_rhs.F,v 1.2 2010/02/15 19:18:55 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 #include "PARAMS.h"
42 #include "GRID.h"
43 #include "SURFACE.h"
44
45 #ifdef ALLOW_AUTODIFF_TAMC
46 #include "tamc.h"
47 #include "tamc_keys.h"
48 #endif /* ALLOW_AUTODIFF_TAMC */
49
50 #include "SMOOTH.h"
51
52 C !INPUT PARAMETERS: ===================================================
53
54
55 INTEGER bi,bj,iMin,iMax,jMin,jMax
56 _RS xA (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
57 _RS yA (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
58 _RS maskUp(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
59 _RL dTdz (nSx,nSy)
60 _RL dTdx (nSx,nSy)
61 _RL dTdy (nSx,nSy)
62 INTEGER myThid
63 INTEGER i,j,k
64
65 _RL fZon (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nR,nSx,nSy)
66 _RL fMer (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nR,nSx,nSy)
67 _RL fVerT (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nR,nSx,nSy)
68 _RL df (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
69 _RL fld_in(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr,nSx,nSy)
70 _RL gt_in(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr,nSx,nSy)
71
72 DO bj=myByLo(myThid),myByHi(myThid)
73 DO bi=myBxLo(myThid),myBxHi(myThid)
74
75
76 c 1rst k loop: initialization
77 DO k=1,Nr
78 DO j=1-OLy,sNy+OLy
79 DO i=1-OLx,sNx+OLx
80 fZon(i,j,k,bi,bj) = 0. _d 0
81 fMer(i,j,k,bi,bj) = 0. _d 0
82 fVerT(i,j,k,bi,bj) = 0. _d 0
83 gt_in(i,j,k,bi,bj) = 0. _d 0
84 ENDDO
85 ENDDO
86 ENDDO
87
88
89 iMin = 1-OLx+1
90 iMax = sNx+OLx-1
91 jMin = 1-OLy+1
92 jMax = sNy+OLy-1
93
94 c 2nd k loop: flux computation
95 DO k=1,Nr
96
97 DO j=1-OLy,sNy+OLy
98 DO i=1-OLx,sNx+OLx
99 df(i,j,bi,bj) = 0. _d 0
100 xA(i,j,bi,bj) = _dyG(i,j,bi,bj)
101 & *drF(k)*_hFacW(i,j,k,bi,bj)
102 yA(i,j,bi,bj) = _dxG(i,j,bi,bj)
103 & *drF(k)*_hFacS(i,j,k,bi,bj)
104 IF (K .EQ. 1) THEN
105 maskUp(i,j,bi,bj) = 0.
106 ELSE
107 maskUp(i,j,bi,bj) =
108 & maskC(i,j,k-1,bi,bj)*maskC(i,j,k,bi,bj)
109 ENDIF
110 ENDDO
111 ENDDO
112
113
114 c ///gmredi_xtr///
115
116 DO j=jMin,jMax
117 DO i=iMin,iMax
118 df(i,j,bi,bj) = df(i,j,bi,bj)
119 & -xA(i,j,bi,bj)
120 & *smooth3D_Kux(i,j,k,bi,bj)
121 & *recip_dxC(i,j,bi,bj)
122 & *(fld_in(i,j,k,bi,bj)-fld_in(i-1,j,k,bi,bj))
123 ENDDO
124 ENDDO
125
126 DO j=jMin,jMax
127 DO i=iMin,iMax
128 dTdz(bi,bj) = 0.5*(
129 & +0.5*recip_drC(k)*
130 & ( maskC(i-1,j,k,bi,bj)*
131 & (fld_in(i-1,j, MAX(k-1,1) ,bi,bj)-fld_in(i-1,j,k,bi,bj))
132 & +maskC( i ,j,k,bi,bj)*
133 & (fld_in( i ,j, MAX(k-1,1) ,bi,bj)-fld_in( i ,j,k,bi,bj))
134 & )
135 & +0.5*recip_drC(MIN(k+1,Nr))*
136 & ( maskC(i-1,j,MIN(k+1,Nr),bi,bj)*
137 & (fld_in(i-1,j,k,bi,bj)-fld_in(i-1,j,MIN(k+1,Nr),bi,bj))
138 & +maskC( i ,j,MIN(k+1,Nr),bi,bj)*
139 & (fld_in( i ,j,k,bi,bj)-fld_in( i ,j,MIN(k+1,Nr),bi,bj))
140 & ) )
141 df(i,j,bi,bj) = df(i,j,bi,bj)
142 & - xA(i,j,bi,bj)*smooth3D_Kuz(i,j,k,bi,bj)*dTdz(bi,bj)
143 ENDDO
144 ENDDO
145
146 DO j=jMin,jMax
147 DO i=iMin,iMax
148 dTdy(bi,bj) = 0.5*(
149 & +0.5*(maskS(i,j,k,bi,bj)
150 & *recip_dyC(i,j,bi,bj)*
151 & (fld_in(i,j,k,bi,bj)-fld_in(i,j-1,k,bi,bj))
152 & +maskS(i,j+1,k,bi,bj)
153 & *recip_dyC(i,j+1,bi,bj)*
154 & (fld_in(i,j+1,k,bi,bj)-fld_in(i,j,k,bi,bj)))
155 & +0.5*(maskS(i-1,j,k,bi,bj)
156 & *recip_dyC(i,j,bi,bj)*
157 & (fld_in(i-1,j,k,bi,bj)-fld_in(i-1,j-1,k,bi,bj))
158 & +maskS(i-1,j+1,k,bi,bj)
159 & *recip_dyC(i,j+1,bi,bj)*
160 & (fld_in(i-1,j+1,k,bi,bj)-fld_in(i-1,j,k,bi,bj)))
161 & )
162 df(i,j,bi,bj) = df(i,j,bi,bj)
163 & - xA(i,j,bi,bj)*smooth3D_Kuy(i,j,k,bi,bj)*dTdy(bi,bj)
164 ENDDO
165 ENDDO
166
167
168 c /// end for x ///
169
170 DO j=jMin,jMax
171 DO i=iMin,iMax
172 fZon(i,j,k,bi,bj) = fZon(i,j,k,bi,bj) + df(i,j,bi,bj)
173 ENDDO
174 ENDDO
175
176 DO j=jMin,jMax
177 DO i=iMin,iMax
178 df(i,j,bi,bj) = 0.
179 ENDDO
180 ENDDO
181
182 c ///gmredi_ytr///
183
184 DO j=jMin,jMax
185 DO i=iMin,iMax
186 df(i,j,bi,bj) = df(i,j,bi,bj)
187 & -yA(i,j,bi,bj)
188 & *smooth3D_Kvy(i,j,k,bi,bj)
189 & *recip_dyC(i,j,bi,bj)
190 & *(fld_in(i,j,k,bi,bj)-fld_in(i,j-1,k,bi,bj))
191 ENDDO
192 ENDDO
193
194 DO j=jMin,jMax
195 DO i=iMin,iMax
196 dTdz(bi,bj) = 0.5*(
197 & +0.5*recip_drC(k)*
198 & ( maskC(i,j-1,k,bi,bj)*
199 & (fld_in(i,j-1,MAX(k-1,1),bi,bj)-fld_in(i,j-1,k,bi,bj))
200 & +maskC(i, j ,k,bi,bj)*
201 & (fld_in(i, j ,MAX(k-1,1),bi,bj)-fld_in(i, j ,k,bi,bj))
202 & )
203 & +0.5*recip_drC(MIN(k+1,Nr))*
204 & ( maskC(i,j-1,MIN(k+1,Nr),bi,bj)*
205 & (fld_in(i,j-1,k,bi,bj)-fld_in(i,j-1,MIN(k+1,Nr),bi,bj))
206 & +maskC(i, j ,MIN(k+1,Nr),bi,bj)*
207 & (fld_in(i, j ,k,bi,bj)-fld_in(i, j ,MIN(k+1,Nr),bi,bj))
208 & ) )
209 df(i,j,bi,bj) = df(i,j,bi,bj)
210 & - yA(i,j,bi,bj)*smooth3D_Kvz(i,j,k,bi,bj)*dTdz(bi,bj)
211 ENDDO
212 ENDDO
213
214 DO j=jMin,jMax
215 DO i=iMin,iMax
216 dTdx(bi,bj) = 0.5*(
217 & +0.5*(maskW(i+1,j,k,bi,bj)
218 & *recip_dxC(i+1,j,bi,bj)*
219 & (fld_in(i+1,j,k,bi,bj)-fld_in(i,j,k,bi,bj))
220 & +maskW(i,j,k,bi,bj)
221 & *recip_dxC(i,j,bi,bj)*
222 & (fld_in(i,j,k,bi,bj)-fld_in(i-1,j,k,bi,bj)))
223 & +0.5*(maskW(i+1,j-1,k,bi,bj)
224 & *recip_dxC(i+1,j,bi,bj)*
225 & (fld_in(i+1,j-1,k,bi,bj)-fld_in(i,j-1,k,bi,bj))
226 & +maskW(i,j-1,k,bi,bj)
227 & *recip_dxC(i,j,bi,bj)*
228 & (fld_in(i,j-1,k,bi,bj)-fld_in(i-1,j-1,k,bi,bj)))
229 & )
230 df(i,j,bi,bj) = df(i,j,bi,bj)
231 & - yA(i,j,bi,bj)*smooth3D_Kvx(i,j,k,bi,bj)*dTdx(bi,bj)
232 ENDDO
233 ENDDO
234
235 c /// end for y ///
236
237 DO j=jMin,jMax
238 DO i=iMin,iMax
239 fMer(i,j,k,bi,bj) = fMer(i,j,k,bi,bj) + df(i,j,bi,bj)
240 ENDDO
241 ENDDO
242
243 DO j=jMin,jMax
244 DO i=iMin,iMax
245 df(i,j,bi,bj) = 0.
246 ENDDO
247 ENDDO
248
249 c /// GAD_DIFF_R ///
250
251 if (.NOT. smooth3DdoImpldiff ) then
252
253 IF (k.gt.1) then
254 DO j=jMin,jMax
255 DO i=iMin,iMax
256 df(i,j,bi,bj) =
257 & -_rA(i,j,bi,bj)
258 & *smooth3D_kappaR(i,j,k,bi,bj)*recip_drC(k)
259 & *(fld_in(i,j,k,bi,bj)
260 & -fld_in(i,j,k-1,bi,bj))*rkSign
261 ENDDO
262 ENDDO
263 ENDIF
264
265 endif
266
267 c ///gmredi rtrans///
268
269 IF (K.GT.1) THEN
270 DO j=jMin,jMax
271 DO i=iMin,iMax
272 dTdx(bi,bj) = 0.5*(
273 & +0.5*(maskW(i+1,j,k,bi,bj)
274 & *recip_dxC(i+1,j,bi,bj)*
275 & (fld_in(i+1,j,k,bi,bj)-fld_in(i,j,k,bi,bj))
276 & +maskW(i,j,k,bi,bj)
277 & *recip_dxC(i,j,bi,bj)*
278 & (fld_in(i,j,k,bi,bj)-fld_in(i-1,j,k,bi,bj)))
279 & +0.5*(maskW(i+1,j,k-1,bi,bj)
280 & *recip_dxC(i+1,j,bi,bj)*
281 & (fld_in(i+1,j,k-1,bi,bj)-fld_in(i,j,k-1,bi,bj))
282 & +maskW(i,j,k-1,bi,bj)
283 & *recip_dxC(i,j,bi,bj)*
284 & (fld_in(i,j,k-1,bi,bj)-fld_in(i-1,j,k-1,bi,bj)))
285 & )
286
287 dTdy(bi,bj) = 0.5*(
288 & +0.5*(maskS(i,j,k,bi,bj)
289 & *recip_dyC(i,j,bi,bj)*
290 & (fld_in(i,j,k,bi,bj)-fld_in(i,j-1,k,bi,bj))
291 & +maskS(i,j+1,k,bi,bj)
292 & *recip_dyC(i,j+1,bi,bj)*
293 & (fld_in(i,j+1,k,bi,bj)-fld_in(i,j,k,bi,bj)))
294 & +0.5*(maskS(i,j,k-1,bi,bj)
295 & *recip_dyC(i,j,bi,bj)*
296 & (fld_in(i,j,k-1,bi,bj)-fld_in(i,j-1,k-1,bi,bj))
297 & +maskS(i,j+1,k-1,bi,bj)
298 & *recip_dyC(i,j+1,bi,bj)*
299 & (fld_in(i,j+1,k-1,bi,bj)-fld_in(i,j,k-1,bi,bj)))
300 & )
301
302 df(i,j,bi,bj) = df(i,j,bi,bj)
303 & - rA(i,j,bi,bj)
304 & *( smooth3D_Kwx(i,j,k,bi,bj)*dTdx(bi,bj)
305 & +smooth3D_Kwy(i,j,k,bi,bj)*dTdy(bi,bj) )
306
307 ENDDO
308 ENDDO
309
310 ENDIF
311
312
313 c /// end for r ///
314
315 IF (K.GT.1) THEN
316 DO j=jMin,jMax
317 DO i=iMin,iMax
318 fVerT(i,j,k-1,bi,bj) = fVerT(i,j,k-1,bi,bj) +
319 & df(i,j,bi,bj)*maskUp(i,j,bi,bj)
320 ENDDO
321 ENDDO
322 ENDIF
323
324 DO j=jMin,jMax
325 DO i=iMin,iMax
326 df(i,j,bi,bj) = 0.
327 ENDDO
328 ENDDO
329
330 ENDDO
331
332 ENDDO
333 ENDDO
334
335 c these exchanges are crucial:
336 CALL EXCH_UV_XYZ_RL(fZon,fMer,.TRUE.,myThid)
337 _EXCH_XYZ_RL ( fVerT, myThid )
338
339 DO bj=myByLo(myThid),myByHi(myThid)
340 DO bi=myBxLo(myThid),myBxHi(myThid)
341 c 3rd k loop: Divergence of fluxes
342 DO k=1,Nr
343 IF (K.GT.1) THEN
344 DO j=jMin,jMax
345 DO i=iMin,iMax
346 gt_in(i,j,k,bi,bj)=gt_in(i,j,k,bi,bj)
347 & -_recip_hFacC(i,j,k,bi,bj)*recip_drF(k)*recip_rA(i,j,bi,bj)
348 & *( (fZon(i+1,j,k,bi,bj)-fZon(i,j,k,bi,bj))
349 & +(fMer(i,j+1,k,bi,bj)-fMer(i,j,k,bi,bj))
350 & +(fVerT(i,j,k,bi,bj)-fVerT(i,j,k-1,bi,bj))*rkSign
351 & )
352 ENDDO
353 ENDDO
354 ELSE
355 DO j=jMin,jMax
356 DO i=iMin,iMax
357 gt_in(i,j,k,bi,bj)=gt_in(i,j,k,bi,bj)
358 & -_recip_hFacC(i,j,k,bi,bj)*recip_drF(k)*recip_rA(i,j,bi,bj)
359 & *( (fZon(i+1,j,k,bi,bj)-fZon(i,j,k,bi,bj))
360 & +(fMer(i,j+1,k,bi,bj)-fMer(i,j,k,bi,bj))
361 & +(fVerT(i,j,k,bi,bj))*rkSign
362 & )
363 ENDDO
364 ENDDO
365 ENDIF
366 ENDDO
367 ENDDO
368 ENDDO
369
370 _EXCH_XYZ_RL ( gt_in , myThid )
371
372 END
373

  ViewVC Help
Powered by ViewVC 1.1.22