/[MITgcm]/MITgcm/model/src/temp_integrate.F
ViewVC logotype

Contents of /MITgcm/model/src/temp_integrate.F

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


Revision 1.2 - (show annotations) (download)
Thu Nov 21 00:16:55 2013 UTC (10 years, 6 months ago) by jmc
Branch: MAIN
Changes since 1.1: +4 -1 lines
add initialisation of rTrans (only for AD)

1 C $Header: /u/gcmpack/MITgcm/model/src/temp_integrate.F,v 1.1 2013/11/19 16:58:38 jmc Exp $
2 C $Name: $
3
4 #include "PACKAGES_CONFIG.h"
5 #include "CPP_OPTIONS.h"
6
7 CBOP
8 C !ROUTINE: TEMP_INTEGRATE
9 C !INTERFACE:
10 SUBROUTINE TEMP_INTEGRATE(
11 I bi, bj, iMin, iMax, jMin, jMax,
12 I uFld, vFld, wFld, KappaRk,
13 I myTime, myIter, myThid )
14 C !DESCRIPTION: \bv
15 C *==========================================================*
16 C | SUBROUTINE TEMP_INTEGRATE
17 C | o Calculate tendency for temperature
18 C | and integrates forward in time.
19 C *==========================================================*
20 C | A procedure called EXTERNAL_FORCING_T is called from
21 C | here. These procedures can be used to add per problem
22 C | heat flux source terms.
23 C | Note: Although it is slightly counter-intuitive the
24 C | EXTERNAL_FORCING routine is not the place to put
25 C | file I/O. Instead files that are required to
26 C | calculate the external source terms are generally
27 C | read during the model main loop. This makes the
28 C | logistics of multi-processing simpler and also
29 C | makes the adjoint generation simpler. It also
30 C | allows for I/O to overlap computation where that
31 C | is supported by hardware.
32 C | Aside from the problem specific term the code here
33 C | forms the tendency terms due to advection and mixing
34 C | The baseline implementation here uses a centered
35 C | difference form for the advection term and a tensorial
36 C | divergence of a flux form for the diffusive term. The
37 C | diffusive term is formulated so that isopycnal mixing and
38 C | GM-style subgrid-scale terms can be incorporated b simply
39 C | setting the diffusion tensor terms appropriately.
40 C *==========================================================*
41 C \ev
42
43 C !USES:
44 IMPLICIT NONE
45 C == GLobal variables ==
46 #include "SIZE.h"
47 #include "DYNVARS.h"
48 #include "EEPARAMS.h"
49 #include "PARAMS.h"
50 #include "RESTART.h"
51 #ifdef ALLOW_GENERIC_ADVDIFF
52 #include "GAD.h"
53 #endif
54 #ifdef ALLOW_AUTODIFF_TAMC
55 # include "tamc.h"
56 # include "tamc_keys.h"
57 #endif
58
59 C !INPUT/OUTPUT PARAMETERS:
60 C == Routine arguments ==
61 C bi, bj, :: tile indices
62 C iMin,iMax :: loop range for called routines
63 C jMin,jMax :: loop range for called routines
64 C uFld,vFld :: Local copy of horizontal velocity field
65 C wFld :: Local copy of vertical velocity field
66 C KappaRk :: Vertical diffusion for Tempertature
67 C myTime :: current time
68 C myIter :: current iteration number
69 C myThid :: my Thread Id. number
70 INTEGER bi, bj, iMin, iMax, jMin, jMax
71 _RL uFld (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
72 _RL vFld (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
73 _RL wFld (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
74 _RL KappaRk(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
75 _RL myTime
76 INTEGER myIter
77 INTEGER myThid
78 CEOP
79
80 #ifdef ALLOW_GENERIC_ADVDIFF
81 C !LOCAL VARIABLES:
82 C k :: vertical index
83 C kM1 :: =k-1 for k>1, =1 for k=1
84 C kUp :: index into 2 1/2D array, toggles between 1|2
85 C kDown :: index into 2 1/2D array, toggles between 2|1
86 C xA :: Tracer cell face area normal to X
87 C yA :: Tracer cell face area normal to X
88 C maskUp :: Land/water mask for Wvel points (interface k)
89 C uTrans :: Zonal volume transport through cell face
90 C vTrans :: Meridional volume transport through cell face
91 C rTrans :: Vertical volume transport at interface k
92 C rTransKp :: Vertical volume transport at inteface k+1
93 C fVerT :: Flux of temperature (T) in the vertical direction
94 C at the upper(U) and lower(D) faces of a cell.
95 INTEGER i, j, k
96 INTEGER kUp, kDown, kM1
97 _RS xA (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
98 _RS yA (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
99 _RS maskUp (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
100 _RL uTrans (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
101 _RL vTrans (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
102 _RL rTrans (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
103 _RL rTransKp(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
104 _RL fVerT (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)
105 _RL gt_AB (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
106 LOGICAL calcAdvection
107 INTEGER iterNb
108 #ifdef ALLOW_ADAMSBASHFORTH_3
109 INTEGER m1, m2
110 #endif
111
112 #ifdef ALLOW_AUTODIFF_TAMC
113 act1 = bi - myBxLo(myThid)
114 max1 = myBxHi(myThid) - myBxLo(myThid) + 1
115 act2 = bj - myByLo(myThid)
116 max2 = myByHi(myThid) - myByLo(myThid) + 1
117 act3 = myThid - 1
118 max3 = nTx*nTy
119 act4 = ikey_dynamics - 1
120 itdkey = (act1 + 1) + act2*max1
121 & + act3*max1*max2
122 & + act4*max1*max2*max3
123 #endif /* ALLOW_AUTODIFF_TAMC */
124
125 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
126
127 calcAdvection = tempAdvection .AND. .NOT.tempMultiDimAdvec
128 iterNb = myIter
129 IF (staggerTimeStep) iterNb = myIter -1
130
131 DO j=1-OLy,sNy+OLy
132 DO i=1-OLx,sNx+OLx
133 fVerT(i,j,1) = 0. _d 0
134 fVerT(i,j,2) = 0. _d 0
135 #ifdef ALLOW_AUTODIFF
136 rTrans(i,j) = 0. _d 0
137 #endif /* ALLOW_AUTODIFF */
138 ENDDO
139 ENDDO
140
141 DO k=Nr,1,-1
142 #ifdef ALLOW_AUTODIFF_TAMC
143 kkey = (itdkey-1)*Nr + k
144 #endif /* ALLOW_AUTODIFF_TAMC */
145 kM1 = MAX(1,k-1)
146 kUp = 1+MOD(k+1,2)
147 kDown= 1+MOD(k,2)
148
149 #ifdef ALLOW_AUTODIFF_TAMC
150 CADJ STORE rtrans(:,:) = comlev1_bibj_k, key=kkey,
151 CADJ & byte=isbyte, kind = isbyte
152 CADJ STORE fVerT(:,:,:) = comlev1_bibj_k, key=kkey,
153 CADJ & byte=isbyte, kind = isbyte
154 CADJ STORE gT(:,:,k,bi,bj) = comlev1_bibj_k, key=kkey,
155 CADJ & byte=isbyte, kind = isbyte
156 # ifdef ALLOW_ADAMSBASHFORTH_3
157 CADJ STORE gtNm(:,:,k,bi,bj,1) = comlev1_bibj_k, key=kkey,
158 CADJ & byte=isbyte, kind = isbyte
159 CADJ STORE gtNm(:,:,k,bi,bj,2) = comlev1_bibj_k, key=kkey,
160 CADJ & byte=isbyte, kind = isbyte
161 # else
162 CADJ STORE gtNm1(:,:,k,bi,bj) = comlev1_bibj_k, key=kkey,
163 CADJ & byte=isbyte, kind = isbyte
164 # endif
165 #endif /* ALLOW_AUTODIFF_TAMC */
166 CALL CALC_ADV_FLOW(
167 I uFld, vFld, wFld,
168 U rTrans,
169 O uTrans, vTrans, rTransKp,
170 O maskUp, xA, yA,
171 I k, bi, bj, myThid )
172
173 #ifdef ALLOW_ADAMSBASHFORTH_3
174 m1 = 1 + MOD(iterNb+1,2)
175 m2 = 1 + MOD( iterNb ,2)
176 CALL GAD_CALC_RHS(
177 I bi, bj, iMin,iMax,jMin,jMax, k, kM1, kUp, kDown,
178 I xA, yA, maskUp, uFld(1-OLx,1-OLy,k),
179 I vFld(1-OLx,1-OLy,k), wFld(1-OLx,1-OLy,k),
180 I uTrans, vTrans, rTrans, rTransKp,
181 I diffKhT, diffK4T, KappaRk(1-OLx,1-OLy,k), diffKr4T,
182 I gtNm(1-OLx,1-OLy,1,1,1,m2), theta, dTtracerLev,
183 I GAD_TEMPERATURE, tempAdvScheme, tempVertAdvScheme,
184 I calcAdvection, tempImplVertAdv, AdamsBashforth_T,
185 I tempVertDiff4, useGMRedi, useKPP,
186 U fVerT, gT,
187 I myTime, myIter, myThid )
188 #else /* ALLOW_ADAMSBASHFORTH_3 */
189 CALL GAD_CALC_RHS(
190 I bi, bj, iMin,iMax,jMin,jMax, k, kM1, kUp, kDown,
191 I xA, yA, maskUp, uFld(1-OLx,1-OLy,k),
192 I vFld(1-OLx,1-OLy,k), wFld(1-OLx,1-OLy,k),
193 I uTrans, vTrans, rTrans, rTransKp,
194 I diffKhT, diffK4T, KappaRk(1-OLx,1-OLy,k), diffKr4T,
195 I gtNm1, theta, dTtracerLev,
196 I GAD_TEMPERATURE, tempAdvScheme, tempVertAdvScheme,
197 I calcAdvection, tempImplVertAdv, AdamsBashforth_T,
198 I tempVertDiff4, useGMRedi, useKPP,
199 U fVerT, gT,
200 I myTime, myIter, myThid )
201 #endif
202
203 C-- External thermal forcing term(s) inside Adams-Bashforth:
204 IF ( tempForcing .AND. tracForcingOutAB.NE.1 )
205 & CALL EXTERNAL_FORCING_T(
206 I iMin, iMax, jMin, jMax, bi, bj, k,
207 I myTime, myThid )
208
209 IF ( AdamsBashforthGt ) THEN
210 #ifdef ALLOW_ADAMSBASHFORTH_3
211 CALL ADAMS_BASHFORTH3(
212 I bi, bj, k, Nr,
213 U gT, gtNm, gt_AB,
214 I tempStartAB, iterNb, myThid )
215 #else
216 CALL ADAMS_BASHFORTH2(
217 I bi, bj, k, Nr,
218 U gT, gtNm1, gt_AB,
219 I tempStartAB, iterNb, myThid )
220 #endif
221 #ifdef ALLOW_DIAGNOSTICS
222 IF ( useDiagnostics ) THEN
223 CALL DIAGNOSTICS_FILL(gt_AB,'AB_gT ',k,1,2,bi,bj,myThid)
224 ENDIF
225 #endif /* ALLOW_DIAGNOSTICS */
226 ENDIF
227
228 C-- External thermal forcing term(s) outside Adams-Bashforth:
229 IF ( tempForcing .AND. tracForcingOutAB.EQ.1 )
230 & CALL EXTERNAL_FORCING_T(
231 I iMin, iMax, jMin, jMax, bi, bj, k,
232 I myTime, myThid )
233
234 #ifdef NONLIN_FRSURF
235 IF (nonlinFreeSurf.GT.0) THEN
236 CALL FREESURF_RESCALE_G(
237 I bi, bj, k,
238 U gT,
239 I myThid )
240 IF ( AdamsBashforthGt ) THEN
241 #ifdef ALLOW_ADAMSBASHFORTH_3
242 # ifdef ALLOW_AUTODIFF_TAMC
243 CADJ STORE gtNm(:,:,k,bi,bj,1) = comlev1_bibj_k, key=kkey,
244 CADJ & byte=isbyte, kind = isbyte
245 CADJ STORE gtNm(:,:,k,bi,bj,2) = comlev1_bibj_k, key=kkey,
246 CADJ & byte=isbyte, kind = isbyte
247 # endif
248 CALL FREESURF_RESCALE_G(
249 I bi, bj, k,
250 U gtNm(1-OLx,1-OLy,1,1,1,1),
251 I myThid )
252 CALL FREESURF_RESCALE_G(
253 I bi, bj, k,
254 U gtNm(1-OLx,1-OLy,1,1,1,2),
255 I myThid )
256 #else
257 CALL FREESURF_RESCALE_G(
258 I bi, bj, k,
259 U gtNm1,
260 I myThid )
261 #endif
262 ENDIF
263 ENDIF
264 #endif /* NONLIN_FRSURF */
265
266 #ifdef ALLOW_ADAMSBASHFORTH_3
267 IF ( AdamsBashforth_T ) THEN
268 CALL TIMESTEP_TRACER(
269 I bi, bj, k, dTtracerLev(k),
270 I gtNm(1-OLx,1-OLy,1,1,1,m2),
271 U gT,
272 I myIter, myThid )
273 ELSE
274 #endif
275 CALL TIMESTEP_TRACER(
276 I bi, bj, k, dTtracerLev(k),
277 I theta,
278 U gT,
279 I myIter, myThid )
280 #ifdef ALLOW_ADAMSBASHFORTH_3
281 ENDIF
282 #endif
283
284 C- end of vertical index (k) loop (Nr:1)
285 ENDDO
286
287 #endif /* ALLOW_GENERIC_ADVDIFF */
288
289 RETURN
290 END

  ViewVC Help
Powered by ViewVC 1.1.22