/[MITgcm]/MITgcm/pkg/ptracers/ptracers_integrate.F
ViewVC logotype

Contents of /MITgcm/pkg/ptracers/ptracers_integrate.F

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


Revision 1.49 - (show annotations) (download)
Thu Nov 21 00:17:39 2013 UTC (10 years, 7 months ago) by jmc
Branch: MAIN
Changes since 1.48: +8 -3 lines
add initialisation of sTrans (only for AD)

1 C $Header: /u/gcmpack/MITgcm/pkg/ptracers/ptracers_integrate.F,v 1.48 2013/11/19 17:01:39 jmc Exp $
2 C $Name: $
3
4 #include "PTRACERS_OPTIONS.h"
5
6 CBOP
7 C !ROUTINE: PTRACERS_INTEGRATE
8
9 C !INTERFACE: ==========================================================
10 SUBROUTINE PTRACERS_INTEGRATE(
11 I bi, bj,
12 I uFld, vFld, wFld,
13 U KappaRk,
14 I myTime, myIter, myThid )
15
16 C !DESCRIPTION:
17 C Calculates tendency for passive tracers and integrates forward
18 C in time.
19
20 C !USES: ===============================================================
21 IMPLICIT NONE
22 #include "SIZE.h"
23 #include "EEPARAMS.h"
24 #include "PARAMS.h"
25 #ifdef ALLOW_LONGSTEP
26 #include "LONGSTEP_PARAMS.h"
27 #endif
28 #include "PTRACERS_SIZE.h"
29 #include "PTRACERS_PARAMS.h"
30 #include "PTRACERS_START.h"
31 #include "PTRACERS_FIELDS.h"
32 #include "GAD.h"
33 #ifdef ALLOW_AUTODIFF_TAMC
34 # include "tamc.h"
35 # include "tamc_keys.h"
36 #endif
37
38 C !INPUT PARAMETERS: ===================================================
39 C bi, bj :: tile indices
40 C uFld, vFld, wFld :: Local copy of velocity field (3 components)
41 C KappaRk :: vertical diffusion used for one passive tracer
42 C myTime :: model time
43 C myIter :: time-step number
44 C myThid :: thread number
45 INTEGER bi, bj
46 _RL uFld (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
47 _RL vFld (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
48 _RL wFld (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
49 _RL KappaRk (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
50 _RL myTime
51 INTEGER myIter
52 INTEGER myThid
53
54 C !OUTPUT PARAMETERS: ==================================================
55 C none
56
57 #ifdef ALLOW_PTRACERS
58
59 C !LOCAL VARIABLES: ====================================================
60 C iTracer :: tracer index
61 C iMin,iMax,jMin,jMax :: loop ranges
62 C k :: vertical level number
63 C kUp,kDown :: toggle indices for even/odd level fluxes
64 C kM1 :: =min(1,k-1)
65 C GAD_TR :: passive tracer id (GAD_TR1+iTracer-1)
66 C xA :: face area at U points in level k
67 C yA :: face area at V points in level k
68 C maskUp :: mask for vertical transport
69 C uTrans :: zonal transport in level k
70 C vTrans :: meridional transport in level k
71 C rTrans :: vertical volume transport at interface k
72 C rTransKp :: vertical volume transport at interface k+1
73 C fVerT :: passive tracer vertical flux
74 INTEGER iTracer
75 INTEGER iMin,iMax,jMin,jMax
76 INTEGER i, j, k
77 INTEGER kUp, kDown, kM1
78 INTEGER GAD_TR
79 _RS xA (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
80 _RS yA (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
81 _RS maskUp (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
82 _RL uTrans (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
83 _RL vTrans (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
84 _RL rTrans (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
85 _RL rTransKp(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
86 _RL fVerT (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)
87 _RL gTr_AB (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
88 LOGICAL calcAdvection
89 INTEGER iterNb
90 _RL dummy(Nr)
91 #ifdef ALLOW_DIAGNOSTICS
92 CHARACTER*8 diagName
93 CHARACTER*4 diagSufx
94 C- Functions:
95 CHARACTER*4 GAD_DIAG_SUFX
96 EXTERNAL GAD_DIAG_SUFX
97 #endif /* ALLOW_DIAGNOSTICS */
98 CEOP
99
100 C- Loop ranges for daughter routines
101 iMin = 1-OLx+2
102 iMax = sNx+OLx-1
103 jMin = 1-OLy+2
104 jMax = sNy+OLy-1
105
106 C- Compute iter at beginning of ptracer time step
107 #ifdef ALLOW_LONGSTEP
108 iterNb = myIter - LS_nIter + 1
109 IF (LS_whenToSample.GE.2) iterNb = myIter - LS_nIter
110 #else
111 iterNb = myIter
112 IF (staggerTimeStep) iterNb = myIter - 1
113 #endif
114
115 C-- Loop over tracers
116 DO iTracer=1,PTRACERS_numInUse
117 IF ( PTRACERS_StepFwd(iTracer) ) THEN
118 GAD_TR = GAD_TR1 + iTracer - 1
119 calcAdvection = .NOT.PTRACERS_MultiDimAdv(iTracer)
120
121 #ifdef ALLOW_AUTODIFF_TAMC
122 act0 = iTracer - 1
123 max0 = PTRACERS_num
124 act1 = bi - myBxLo(myThid)
125 max1 = myBxHi(myThid) - myBxLo(myThid) + 1
126 act2 = bj - myByLo(myThid)
127 max2 = myByHi(myThid) - myByLo(myThid) + 1
128 act3 = myThid - 1
129 max3 = nTx*nTy
130 act4 = ikey_dynamics - 1
131 iptrkey = (act0 + 1)
132 & + act1*max0
133 & + act2*max0*max1
134 & + act3*max0*max1*max2
135 & + act4*max0*max1*max2*max3
136 #endif /* ALLOW_AUTODIFF_TAMC */
137
138 DO j=1-OLy,sNy+OLy
139 DO i=1-OLx,sNx+OLx
140 fVerT(i,j,1) = 0. _d 0
141 fVerT(i,j,2) = 0. _d 0
142 ENDDO
143 ENDDO
144 #ifdef ALLOW_AUTODIFF
145 DO j=1-OLy,sNy+OLy
146 DO i=1-OLx,sNx+OLx
147 rTrans(i,j) = 0. _d 0
148 ENDDO
149 ENDDO
150 DO k=1,Nr
151 DO j=1-OLy,sNy+OLy
152 DO i=1-OLx,sNx+OLx
153 kappaRk(i,j,k) = 0. _d 0
154 ENDDO
155 ENDDO
156 ENDDO
157 #endif /* ALLOW_AUTODIFF */
158
159 IF ( .NOT.implicitDiffusion ) THEN
160 CALL CALC_3D_DIFFUSIVITY(
161 I bi, bj, iMin,iMax,jMin,jMax,
162 I GAD_TR,
163 I PTRACERS_useGMRedi(iTracer), PTRACERS_useKPP(iTracer),
164 O kappaRk,
165 I myThid)
166 ENDIF
167
168 DO k=Nr,1,-1
169 #ifdef ALLOW_AUTODIFF_TAMC
170 kkey = (iptrkey-1)*Nr + k
171 #endif /* ALLOW_AUTODIFF_TAMC */
172
173 kM1 = MAX(1,k-1)
174 kUp = 1+MOD(k+1,2)
175 kDown= 1+MOD(k,2)
176
177 #ifdef ALLOW_AUTODIFF_TAMC
178 CADJ STORE rtrans(:,:) = comlev1_bibj_k_ptracers,
179 CADJ & key = kkey, byte = isbyte, kind = isbyte
180 CADJ STORE fVerT(:,:,:) = comlev1_bibj_k_ptracers,
181 CADJ & key = kkey, byte = isbyte, kind = isbyte
182 CADJ STORE gPtr(:,:,k,bi,bj,iTracer) = comlev1_bibj_k_ptracers,
183 CADJ & key = kkey, byte = isbyte, kind = isbyte
184 CADJ STORE gpTrNm1(:,:,k,bi,bj,iTracer) = comlev1_bibj_k_ptracers,
185 CADJ & key = kkey, byte = isbyte, kind = isbyte
186 #endif /* ALLOW_AUTODIFF_TAMC */
187 CALL CALC_ADV_FLOW(
188 I uFld, vFld, wFld,
189 U rTrans,
190 O uTrans, vTrans, rTransKp,
191 O maskUp, xA, yA,
192 I k, bi, bj, myThid )
193
194 C- Calculate active tracer tendencies (gPtr) due to internal processes
195 C (advection, [explicit] diffusion, parameterizations,...)
196 CALL GAD_CALC_RHS(
197 I bi,bj, iMin,iMax,jMin,jMax, k,kM1, kUp,kDown,
198 I xA, yA, maskUp, uFld(1-OLx,1-OLy,k),
199 I vFld(1-OLx,1-OLy,k), wFld(1-OLx,1-OLy,k),
200 I uTrans, vTrans, rTrans, rTransKp,
201 I PTRACERS_diffKh(iTracer),
202 I PTRACERS_diffK4(iTracer),
203 I KappaRk(1-OLx,1-OLy,k), dummy,
204 I gpTrNm1(1-OLx,1-OLy,1,1,1,iTracer),
205 I pTracer(1-OLx,1-OLy,1,1,1,iTracer),
206 I PTRACERS_dTLev, GAD_TR,
207 I PTRACERS_advScheme(iTracer),
208 I PTRACERS_advScheme(iTracer),
209 I calcAdvection, PTRACERS_ImplVertAdv(iTracer),
210 I .FALSE., .FALSE.,
211 I PTRACERS_useGMRedi(iTracer),
212 I PTRACERS_useKPP(iTracer),
213 U fVerT, gPtr(1-OLx,1-OLy,1,1,1,iTracer),
214 I myTime, myIter, myThid )
215
216 C- External forcing term(s)
217 IF ( tracForcingOutAB.NE.1 )
218 & CALL PTRACERS_FORCING(
219 I bi,bj,iMin,iMax,jMin,jMax,k,iTracer,
220 U gPtr(1-OLx,1-OLy,1,1,1,iTracer),
221 I surfaceForcingPTr(1-OLx,1-OLy,1,1,iTracer),
222 I myIter, myTime, myThid )
223
224 C- If using Adams-Bashforth II, then extrapolate tendencies
225 C gPtr is now the tracer tendency for explicit advection/diffusion
226
227 C If matrix is being computed, skip call to S/R ADAMS_BASHFORTH2 to
228 C prevent gPtr from being replaced by the average of gPtr and gpTrNm1.
229 IF ( .NOT.useMATRIX .AND.
230 & PTRACERS_AdamsBashGtr(iTracer) ) THEN
231
232 CALL ADAMS_BASHFORTH2(
233 I bi, bj, k, Nr,
234 U gPtr(1-OLx,1-OLy,1,1,1,iTracer),
235 U gpTrNm1(1-OLx,1-OLy,1,1,1,iTracer),
236 U gTr_AB,
237 I PTRACERS_startAB(iTracer), iterNb, myThid )
238 #ifdef ALLOW_DIAGNOSTICS
239 IF ( useDiagnostics ) THEN
240 diagSufx = GAD_DIAG_SUFX( GAD_TR, myThid )
241 diagName = 'AB_g'//diagSufx
242 CALL DIAGNOSTICS_FILL(gTr_AB,diagName,k,1,2,bi,bj,myThid)
243 ENDIF
244 #endif /* ALLOW_DIAGNOSTICS */
245 ENDIF
246
247 C- External forcing term(s)
248 IF ( tracForcingOutAB.EQ.1 )
249 & CALL PTRACERS_FORCING(
250 I bi,bj,iMin,iMax,jMin,jMax,k,iTracer,
251 U gPtr(1-OLx,1-OLy,1,1,1,iTracer),
252 I surfaceForcingPTr(1-OLx,1-OLy,1,1,iTracer),
253 I myIter, myTime, myThid )
254
255 #ifdef NONLIN_FRSURF
256 C- Account for change in level thickness
257 IF (nonlinFreeSurf.GT.0) THEN
258 CALL FREESURF_RESCALE_G(
259 I bi,bj,k,
260 U gPtr(1-OLx,1-OLy,1,1,1,iTracer),
261 I myThid )
262 IF ( PTRACERS_AdamsBashGtr(iTracer) )
263 & CALL FREESURF_RESCALE_G(
264 I bi,bj,k,
265 U gpTrNm1(1-OLx,1-OLy,1,1,1,iTracer),
266 I myThid )
267 ENDIF
268 #endif /* NONLIN_FRSURF */
269
270 C- Integrate forward in time, storing in gPtr: G=T+dt*G
271 CALL TIMESTEP_TRACER(
272 I bi, bj, k, PTRACERS_dTLev(k),
273 I pTracer(1-OLx,1-OLy,1,1,1,iTracer),
274 U gPtr(1-OLx,1-OLy,1,1,1,iTracer),
275 I myIter, myThid )
276
277 C- end of vertical index (k) loop (Nr:1)
278 ENDDO
279
280 C-- end of tracer loop
281 ENDIF
282 ENDDO
283
284 #endif /* ALLOW_PTRACERS */
285
286 RETURN
287 END

  ViewVC Help
Powered by ViewVC 1.1.22