/[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.34 - (hide annotations) (download)
Tue Oct 28 22:57:59 2003 UTC (20 years, 8 months ago) by edhill
Branch: MAIN
CVS Tags: checkpoint51p_post
Changes since 1.33: +10 -11 lines 
 o add a "cd_code" package and update all the verification tests
   so that they use the new package instead of "INCLUDE_CD_CODE"
1 edhill 1.34 C $Header: /u/u3/gcmpack/MITgcm/model/src/timestep.F,v 1.33 2003/10/09 04:19:18 edhill 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.32 I myTime, myIter, myThid )
13 cnh 1.27 C !DESCRIPTION: \bv
14     C *==========================================================*
15     C | S/R TIMESTEP
16     C | o Step model fields forward in time
17     C *==========================================================*
18     C \ev
19    
20     C !USES:
21 jmc 1.21 IMPLICIT NONE
22 heimbach 1.18 C == Global variables ==
23 cnh 1.1 #include "SIZE.h"
24     #include "DYNVARS.h"
25 cnh 1.11 #include "EEPARAMS.h"
26 cnh 1.1 #include "PARAMS.h"
27     #include "GRID.h"
28 jmc 1.26 #include "SURFACE.h"
29 heimbach 1.18
30 cnh 1.27 C !INPUT/OUTPUT PARAMETERS:
31 cnh 1.1 C == Routine Arguments ==
32 jmc 1.29 C dPhiHydX,Y :: Gradient (X & Y directions) of Hydrostatic Potential
33     C phiSurfX :: gradient of Surface potential (Pressure/rho, ocean)
34     C phiSurfY :: or geopotential (atmos) in X and Y direction
35 cnh 1.1 INTEGER bi,bj,iMin,iMax,jMin,jMax
36 jmc 1.32 INTEGER k
37 jmc 1.29 _RL dPhiHydX(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
38     _RL dPhiHydY(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
39 jmc 1.21 _RL phiSurfX(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
40     _RL phiSurfY(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
41 jmc 1.32 _RL myTime
42 adcroft 1.17 INTEGER myIter, myThid
43 heimbach 1.18
44 cnh 1.27 C !LOCAL VARIABLES:
45 cnh 1.1 C == Local variables ==
46 jmc 1.32 LOGICAL momForcing_In_AB
47 adcroft 1.4 INTEGER i,j
48 adcroft 1.7 _RL ab15,ab05
49 jmc 1.21 _RL phxFac,phyFac, psFac
50 jmc 1.26 _RL gUtmp(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
51     _RL gVtmp(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
52 edhill 1.34 #ifdef ALLOW_CD_CODE
53 jmc 1.32 _RL guCor(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
54     _RL gvCor(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
55     #endif
56 cnh 1.27 CEOP
57 cnh 1.1
58     C Adams-Bashforth timestepping weights
59 adcroft 1.25 IF (myIter .EQ. 0) THEN
60     ab15=1.0
61     ab05=0.0
62     ELSE
63 adcroft 1.17 ab15=1.5+abeps
64     ab05=-0.5-abeps
65 adcroft 1.25 ENDIF
66 cnh 1.1
67 jmc 1.29 C-- stagger time step: grad Phi_Hyp is not in gU,gV => add it in this S/R
68 jmc 1.30 IF (staggerTimeStep) THEN
69 jmc 1.29 phxFac = pfFacMom
70     phyFac = pfFacMom
71     ELSE
72     phxFac = 0.
73     phyFac = 0.
74     ENDIF
75    
76     C-- explicit part of the surface potential gradient is added in this S/R
77     psFac = pfFacMom*(1. _d 0 - implicSurfPress)
78    
79 jmc 1.32 C-- including or excluding momentum forcing from Adams-Bashforth:
80     momForcing_In_AB = forcing_In_AB
81     momForcing_In_AB = .TRUE.
82    
83 jmc 1.26 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
84 jmc 1.32
85     C- Initialize local arrays (not really necessary but safer)
86     DO j=1-Oly,sNy+Oly
87     DO i=1-Olx,sNx+Olx
88     gUtmp(i,j) = 0. _d 0
89     gVtmp(i,j) = 0. _d 0
90 edhill 1.34 #ifdef ALLOW_CD_CODE
91 jmc 1.32 guCor(i,j) = 0. _d 0
92     gvCor(i,j) = 0. _d 0
93     #endif
94     ENDDO
95     ENDDO
96    
97     C-- Forcing term inside the Adams-Bashforth:
98     IF (momForcing .AND. momForcing_In_AB) THEN
99     CALL EXTERNAL_FORCING_U(
100     I iMin,iMax,jMin,jMax,bi,bj,k,
101     I myTime,myThid)
102     CALL EXTERNAL_FORCING_V(
103     I iMin,iMax,jMin,jMax,bi,bj,k,
104     I myTime,myThid)
105     ENDIF
106    
107     C- Compute effective gU,gV_[n+1/2] terms (including Adams-Bashforth weights)
108     C and save gU,gV_[n] into guNm1,gvNm1 for the next time step.
109 jmc 1.26 DO j=jMin,jMax
110     DO i=iMin,iMax
111     gUtmp(i,j) = ab15*gU(i,j,k,bi,bj)
112 jmc 1.32 & + ab05*guNm1(i,j,k,bi,bj)
113     gVtmp(i,j) = ab15*gV(i,j,k,bi,bj)
114     & + ab05*gvNm1(i,j,k,bi,bj)
115     guNm1(i,j,k,bi,bj)= gU(i,j,k,bi,bj)
116     gvNm1(i,j,k,bi,bj)= gV(i,j,k,bi,bj)
117     gU(i,j,k,bi,bj) = gUtmp(i,j)
118     gV(i,j,k,bi,bj) = gVtmp(i,j)
119 jmc 1.26 ENDDO
120     ENDDO
121    
122 jmc 1.32 C-- Forcing term outside the Adams-Bashforth:
123     C (not recommanded with CD-scheme ON)
124     IF (momForcing .AND. .NOT.momForcing_In_AB) THEN
125     CALL EXTERNAL_FORCING_U(
126     I iMin,iMax,jMin,jMax,bi,bj,k,
127     I myTime,myThid)
128     CALL EXTERNAL_FORCING_V(
129     I iMin,iMax,jMin,jMax,bi,bj,k,
130     I myTime,myThid)
131     IF (useCDscheme) THEN
132     C- for CD-scheme, compute gU,Vtmp = gU,V^n + forcing
133 jmc 1.28 DO j=jMin,jMax
134     DO i=iMin,iMax
135 jmc 1.32 gUtmp(i,j) = gU(i,j,k,bi,bj)-gUtmp(i,j)
136     & + guNm1(i,j,k,bi,bj)
137     gVtmp(i,j) = gV(i,j,k,bi,bj)-gVtmp(i,j)
138     & + gvNm1(i,j,k,bi,bj)
139 jmc 1.28 ENDDO
140     ENDDO
141     ELSE
142     DO j=jMin,jMax
143     DO i=iMin,iMax
144 jmc 1.32 gUtmp(i,j) = gU(i,j,k,bi,bj)
145     gVtmp(i,j) = gV(i,j,k,bi,bj)
146 jmc 1.28 ENDDO
147 jmc 1.26 ENDDO
148 jmc 1.28 ENDIF
149 jmc 1.32 ELSEIF ( useCDscheme) THEN
150     DO j=jMin,jMax
151     DO i=iMin,iMax
152     gUtmp(i,j) = guNm1(i,j,k,bi,bj)
153     gVtmp(i,j) = gvNm1(i,j,k,bi,bj)
154     ENDDO
155     ENDDO
156 jmc 1.26 ENDIF
157    
158 edhill 1.34 #ifdef ALLOW_CD_CODE
159 jmc 1.32 IF (useCDscheme) THEN
160     C- Step forward D-grid velocity using C-grid gU,Vtmp = gU,V^n + forcing
161     C and return coriolis terms on C-grid (guCor,gvCor)
162 edhill 1.34 CALL CD_SCHEME(
163     I bi,bj,k, dPhiHydX,dPhiHydY, gUtmp,gVtmp,
164     O guCor,gvCor,
165     I myTime, myIter, myThid)
166 jmc 1.32 DO j=jMin,jMax
167     DO i=iMin,iMax
168     gUtmp(i,j) = gU(i,j,k,bi,bj)
169     & + guCor(i,j)
170     gVtmp(i,j) = gV(i,j,k,bi,bj)
171     & + gvCor(i,j)
172     ENDDO
173 cnh 1.1 ENDDO
174 jmc 1.32 ENDIF
175 edhill 1.34 #endif /* ALLOW_CD_CODE */
176 adcroft 1.19
177 jmc 1.26 #ifdef NONLIN_FRSURF
178     IF (.NOT. vectorInvariantMomentum
179     & .AND. nonlinFreeSurf.GT.1) THEN
180 jmc 1.28 IF (select_rStar.GT.0) THEN
181     DO j=jMin,jMax
182     DO i=iMin,iMax
183 jmc 1.32 gUtmp(i,j) = gUtmp(i,j)/rStarExpW(i,j,bi,bj)
184 jmc 1.28 gVtmp(i,j) = gVtmp(i,j)/rStarExpS(i,j,bi,bj)
185     ENDDO
186     ENDDO
187     ELSE
188     DO j=jMin,jMax
189     DO i=iMin,iMax
190 jmc 1.32 IF ( k.EQ.ksurfW(i,j,bi,bj) ) THEN
191     gUtmp(i,j) = gUtmp(i,j)
192     & *hFacW(i,j,k,bi,bj)/hFac_surfW(i,j,bi,bj)
193     ENDIF
194 jmc 1.28 IF ( k.EQ.ksurfS(i,j,bi,bj) ) THEN
195 jmc 1.26 gVtmp(i,j) = gVtmp(i,j)
196     & *hFacS(i,j,k,bi,bj)/hFac_surfS(i,j,bi,bj)
197 jmc 1.28 ENDIF
198     ENDDO
199 jmc 1.26 ENDDO
200 jmc 1.28 ENDIF
201 jmc 1.26 ENDIF
202 jmc 1.32 #endif /* NONLIN_FRSURF */
203    
204     C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
205    
206     C Step forward zonal velocity (store in Gu)
207     DO j=jMin,jMax
208     DO i=iMin,iMax
209     gU(i,j,k,bi,bj) = uVel(i,j,k,bi,bj)
210     & +deltaTmom*(
211     & gUtmp(i,j)
212     & - psFac*phiSurfX(i,j)
213     & - phxFac*dPhiHydX(i,j)
214     & )*_maskW(i,j,k,bi,bj)
215     ENDDO
216     ENDDO
217 jmc 1.26
218 adcroft 1.7 C Step forward meridional velocity (store in Gv)
219 adcroft 1.19 DO j=jMin,jMax
220 cnh 1.1 DO i=iMin,iMax
221 jmc 1.32 gV(i,j,k,bi,bj) = vVel(i,j,k,bi,bj)
222 jmc 1.26 & +deltaTmom*(
223     & gVtmp(i,j)
224     & - psFac*phiSurfY(i,j)
225 jmc 1.29 & - phyFac*dPhiHydY(i,j)
226 adcroft 1.7 & )*_maskS(i,j,k,bi,bj)
227 cnh 1.1 ENDDO
228 adcroft 1.19 ENDDO
229 cnh 1.1
230     RETURN
231     END
  ViewVC Help
Powered by ViewVC 1.1.22