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

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

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


Revision 1.2 - (hide 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 jmc 1.2 C $Header: /u/gcmpack/MITgcm/model/src/temp_integrate.F,v 1.1 2013/11/19 16:58:38 jmc Exp $
2 jmc 1.1 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 jmc 1.2 #ifdef ALLOW_AUTODIFF
136     rTrans(i,j) = 0. _d 0
137     #endif /* ALLOW_AUTODIFF */
138 jmc 1.1 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