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

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

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

Revision 1.36 - (hide annotations) (download)
Wed Nov 10 03:02:00 2004 UTC (19 years, 8 months ago) by jmc
Branch: MAIN
CVS Tags: checkpoint57g_pre, checkpoint57b_post, checkpoint56b_post, checkpoint57d_post, checkpoint57, checkpoint56, checkpoint57f_post, checkpoint57a_post, checkpoint57c_post, checkpoint57c_pre, checkpoint55j_post, checkpoint57e_post, eckpoint57e_pre, checkpoint56a_post, checkpoint57f_pre, checkpoint56c_post, checkpoint57a_pre
Changes since 1.35: +55 -10 lines 
isolate dissipation tendency (allow to keep it out off AB)
 note: only implemented in vector-invariant form.
1 jmc 1.36 C $Header: /u/gcmpack/MITgcm/model/src/timestep.F,v 1.35 2003/10/30 18:44:26 heimbach Exp $
2 jmc 1.21 C $Name: $
3 cnh 1.1
4 edhill 1.34 #include "PACKAGES_CONFIG.h"
5 adcroft 1.10 #include "CPP_OPTIONS.h"
6 cnh 1.1
7 cnh 1.27 CBOP
8     C !ROUTINE: TIMESTEP
9     C !INTERFACE:
10 jmc 1.32 SUBROUTINE TIMESTEP( bi, bj, iMin, iMax, jMin, jMax, k,
11 jmc 1.31 I dPhiHydX,dPhiHydY, phiSurfX, phiSurfY,
12 jmc 1.36 I guDissip, gvDissip,
13 jmc 1.32 I myTime, myIter, myThid )
14 cnh 1.27 C !DESCRIPTION: \bv
15     C *==========================================================*
16     C | S/R TIMESTEP
17     C | o Step model fields forward in time
18     C *==========================================================*
19     C \ev
20    
21     C !USES:
22 jmc 1.21 IMPLICIT NONE
23 heimbach 1.18 C == Global variables ==
24 cnh 1.1 #include "SIZE.h"
25     #include "DYNVARS.h"
26 cnh 1.11 #include "EEPARAMS.h"
27 cnh 1.1 #include "PARAMS.h"
28     #include "GRID.h"
29 jmc 1.26 #include "SURFACE.h"
30 heimbach 1.18
31 cnh 1.27 C !INPUT/OUTPUT PARAMETERS:
32 cnh 1.1 C == Routine Arguments ==
33 jmc 1.29 C dPhiHydX,Y :: Gradient (X & Y directions) of Hydrostatic Potential
34     C phiSurfX :: gradient of Surface potential (Pressure/rho, ocean)
35     C phiSurfY :: or geopotential (atmos) in X and Y direction
36 jmc 1.36 C guDissip :: dissipation tendency (all explicit terms), u component
37     C gvDissip :: dissipation tendency (all explicit terms), v component
38    
39 cnh 1.1 INTEGER bi,bj,iMin,iMax,jMin,jMax
40 jmc 1.32 INTEGER k
41 jmc 1.29 _RL dPhiHydX(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
42     _RL dPhiHydY(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
43 jmc 1.21 _RL phiSurfX(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
44     _RL phiSurfY(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
45 jmc 1.36 _RL guDissip(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
46     _RL gvDissip(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
47 jmc 1.32 _RL myTime
48 adcroft 1.17 INTEGER myIter, myThid
49 heimbach 1.18
50 cnh 1.27 C !LOCAL VARIABLES:
51 cnh 1.1 C == Local variables ==
52 jmc 1.32 LOGICAL momForcing_In_AB
53 jmc 1.36 LOGICAL momDissip_In_AB
54 adcroft 1.4 INTEGER i,j
55 adcroft 1.7 _RL ab15,ab05
56 jmc 1.21 _RL phxFac,phyFac, psFac
57 jmc 1.26 _RL gUtmp(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
58     _RL gVtmp(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
59 edhill 1.34 #ifdef ALLOW_CD_CODE
60 jmc 1.32 _RL guCor(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
61     _RL gvCor(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
62     #endif
63 cnh 1.27 CEOP
64 cnh 1.1
65     C Adams-Bashforth timestepping weights
66 adcroft 1.25 IF (myIter .EQ. 0) THEN
67     ab15=1.0
68     ab05=0.0
69     ELSE
70 adcroft 1.17 ab15=1.5+abeps
71     ab05=-0.5-abeps
72 adcroft 1.25 ENDIF
73 cnh 1.1
74 jmc 1.29 C-- explicit part of the surface potential gradient is added in this S/R
75     psFac = pfFacMom*(1. _d 0 - implicSurfPress)
76    
77 jmc 1.32 C-- including or excluding momentum forcing from Adams-Bashforth:
78     momForcing_In_AB = forcing_In_AB
79     momForcing_In_AB = .TRUE.
80 jmc 1.36 momDissip_In_AB = .TRUE.
81 jmc 1.32
82 jmc 1.26 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
83 jmc 1.32
84     C- Initialize local arrays (not really necessary but safer)
85     DO j=1-Oly,sNy+Oly
86     DO i=1-Olx,sNx+Olx
87     gUtmp(i,j) = 0. _d 0
88     gVtmp(i,j) = 0. _d 0
89 edhill 1.34 #ifdef ALLOW_CD_CODE
90 jmc 1.32 guCor(i,j) = 0. _d 0
91     gvCor(i,j) = 0. _d 0
92     #endif
93     ENDDO
94     ENDDO
95    
96 jmc 1.36 C-- Stagger time step: grad Phi_Hyp will be added later
97     IF (staggerTimeStep) THEN
98     phxFac = pfFacMom
99     phyFac = pfFacMom
100     ELSE
101     C-- Synchronous time step: add grad Phi_Hyp to gU,gV before doing Adams-Bashforth
102     C note: already done in S/R mom_vecinv and mom_fluxform but would be better
103     C to add it to gU,gV here.
104     c DO j=jMin,jMax
105     c DO i=iMin,iMax
106     c gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj) - pfFacMom*dPhiHydX(i,j)
107     c gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj) - pfFacMom*dPhiHydY(i,j)
108     c ENDDO
109     c ENDDO
110     phxFac = 0.
111     phyFac = 0.
112     ENDIF
113    
114     #ifdef ALLOW_MOM_VECINV
115     C-- Dissipation term inside the Adams-Bashforth:
116     C note: already in gU,gV if using fluxform
117     IF ( momViscosity .AND. momDissip_In_AB
118     & .AND. vectorInvariantMomentum ) THEN
119     DO j=jMin,jMax
120     DO i=iMin,iMax
121     gU(i,j,k,bi,bj) = gU(i,j,k,bi,bj) + guDissip(i,j)
122     gV(i,j,k,bi,bj) = gV(i,j,k,bi,bj) + gvDissip(i,j)
123     ENDDO
124     ENDDO
125     ENDIF
126     #endif
127    
128 jmc 1.32 C-- Forcing term inside the Adams-Bashforth:
129     IF (momForcing .AND. momForcing_In_AB) THEN
130     CALL EXTERNAL_FORCING_U(
131     I iMin,iMax,jMin,jMax,bi,bj,k,
132     I myTime,myThid)
133     CALL EXTERNAL_FORCING_V(
134     I iMin,iMax,jMin,jMax,bi,bj,k,
135     I myTime,myThid)
136     ENDIF
137    
138     C- Compute effective gU,gV_[n+1/2] terms (including Adams-Bashforth weights)
139     C and save gU,gV_[n] into guNm1,gvNm1 for the next time step.
140 jmc 1.26 DO j=jMin,jMax
141     DO i=iMin,iMax
142     gUtmp(i,j) = ab15*gU(i,j,k,bi,bj)
143 jmc 1.32 & + ab05*guNm1(i,j,k,bi,bj)
144     gVtmp(i,j) = ab15*gV(i,j,k,bi,bj)
145     & + ab05*gvNm1(i,j,k,bi,bj)
146     guNm1(i,j,k,bi,bj)= gU(i,j,k,bi,bj)
147     gvNm1(i,j,k,bi,bj)= gV(i,j,k,bi,bj)
148     gU(i,j,k,bi,bj) = gUtmp(i,j)
149     gV(i,j,k,bi,bj) = gVtmp(i,j)
150 jmc 1.26 ENDDO
151     ENDDO
152    
153 jmc 1.32 C-- Forcing term outside the Adams-Bashforth:
154     C (not recommanded with CD-scheme ON)
155     IF (momForcing .AND. .NOT.momForcing_In_AB) THEN
156     CALL EXTERNAL_FORCING_U(
157     I iMin,iMax,jMin,jMax,bi,bj,k,
158     I myTime,myThid)
159     CALL EXTERNAL_FORCING_V(
160     I iMin,iMax,jMin,jMax,bi,bj,k,
161     I myTime,myThid)
162     IF (useCDscheme) THEN
163     C- for CD-scheme, compute gU,Vtmp = gU,V^n + forcing
164 jmc 1.28 DO j=jMin,jMax
165     DO i=iMin,iMax
166 jmc 1.32 gUtmp(i,j) = gU(i,j,k,bi,bj)-gUtmp(i,j)
167     & + guNm1(i,j,k,bi,bj)
168     gVtmp(i,j) = gV(i,j,k,bi,bj)-gVtmp(i,j)
169     & + gvNm1(i,j,k,bi,bj)
170 jmc 1.28 ENDDO
171     ENDDO
172     ELSE
173     DO j=jMin,jMax
174     DO i=iMin,iMax
175 jmc 1.32 gUtmp(i,j) = gU(i,j,k,bi,bj)
176     gVtmp(i,j) = gV(i,j,k,bi,bj)
177 jmc 1.28 ENDDO
178 jmc 1.26 ENDDO
179 jmc 1.28 ENDIF
180 jmc 1.32 ELSEIF ( useCDscheme) THEN
181     DO j=jMin,jMax
182     DO i=iMin,iMax
183     gUtmp(i,j) = guNm1(i,j,k,bi,bj)
184     gVtmp(i,j) = gvNm1(i,j,k,bi,bj)
185     ENDDO
186     ENDDO
187 jmc 1.26 ENDIF
188    
189 edhill 1.34 #ifdef ALLOW_CD_CODE
190 jmc 1.32 IF (useCDscheme) THEN
191     C- Step forward D-grid velocity using C-grid gU,Vtmp = gU,V^n + forcing
192     C and return coriolis terms on C-grid (guCor,gvCor)
193 heimbach 1.35 CALL CD_CODE_SCHEME(
194 edhill 1.34 I bi,bj,k, dPhiHydX,dPhiHydY, gUtmp,gVtmp,
195     O guCor,gvCor,
196     I myTime, myIter, myThid)
197 jmc 1.32 DO j=jMin,jMax
198     DO i=iMin,iMax
199     gUtmp(i,j) = gU(i,j,k,bi,bj)
200     & + guCor(i,j)
201     gVtmp(i,j) = gV(i,j,k,bi,bj)
202     & + gvCor(i,j)
203     ENDDO
204 cnh 1.1 ENDDO
205 jmc 1.32 ENDIF
206 edhill 1.34 #endif /* ALLOW_CD_CODE */
207 adcroft 1.19
208 jmc 1.26 #ifdef NONLIN_FRSURF
209     IF (.NOT. vectorInvariantMomentum
210     & .AND. nonlinFreeSurf.GT.1) THEN
211 jmc 1.28 IF (select_rStar.GT.0) THEN
212     DO j=jMin,jMax
213     DO i=iMin,iMax
214 jmc 1.32 gUtmp(i,j) = gUtmp(i,j)/rStarExpW(i,j,bi,bj)
215 jmc 1.28 gVtmp(i,j) = gVtmp(i,j)/rStarExpS(i,j,bi,bj)
216     ENDDO
217     ENDDO
218     ELSE
219     DO j=jMin,jMax
220     DO i=iMin,iMax
221 jmc 1.32 IF ( k.EQ.ksurfW(i,j,bi,bj) ) THEN
222     gUtmp(i,j) = gUtmp(i,j)
223     & *hFacW(i,j,k,bi,bj)/hFac_surfW(i,j,bi,bj)
224     ENDIF
225 jmc 1.28 IF ( k.EQ.ksurfS(i,j,bi,bj) ) THEN
226 jmc 1.26 gVtmp(i,j) = gVtmp(i,j)
227     & *hFacS(i,j,k,bi,bj)/hFac_surfS(i,j,bi,bj)
228 jmc 1.28 ENDIF
229     ENDDO
230 jmc 1.26 ENDDO
231 jmc 1.28 ENDIF
232 jmc 1.26 ENDIF
233 jmc 1.32 #endif /* NONLIN_FRSURF */
234    
235     C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
236    
237 jmc 1.36 #ifdef ALLOW_MOM_VECINV
238     C-- Dissipation term outside the Adams-Bashforth:
239     C note: only implemented with vecinv formulation
240     IF ( momViscosity .AND. .NOT.momDissip_In_AB
241     & .AND. vectorInvariantMomentum ) THEN
242     DO j=jMin,jMax
243     DO i=iMin,iMax
244     gUtmp(i,j) = gUtmp(i,j) + guDissip(i,j)
245     gVtmp(i,j) = gVtmp(i,j) + gvDissip(i,j)
246     ENDDO
247     ENDDO
248     ENDIF
249     #endif
250    
251 jmc 1.32 C Step forward zonal velocity (store in Gu)
252     DO j=jMin,jMax
253     DO i=iMin,iMax
254     gU(i,j,k,bi,bj) = uVel(i,j,k,bi,bj)
255     & +deltaTmom*(
256     & gUtmp(i,j)
257     & - psFac*phiSurfX(i,j)
258     & - phxFac*dPhiHydX(i,j)
259     & )*_maskW(i,j,k,bi,bj)
260     ENDDO
261     ENDDO
262 jmc 1.26
263 adcroft 1.7 C Step forward meridional velocity (store in Gv)
264 adcroft 1.19 DO j=jMin,jMax
265 cnh 1.1 DO i=iMin,iMax
266 jmc 1.32 gV(i,j,k,bi,bj) = vVel(i,j,k,bi,bj)
267 jmc 1.26 & +deltaTmom*(
268     & gVtmp(i,j)
269     & - psFac*phiSurfY(i,j)
270 jmc 1.29 & - phyFac*dPhiHydY(i,j)
271 adcroft 1.7 & )*_maskS(i,j,k,bi,bj)
272 cnh 1.1 ENDDO
273 adcroft 1.19 ENDDO
274 cnh 1.1
275     RETURN
276     END
  ViewVC Help
Powered by ViewVC 1.1.22