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

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

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

revision 1.11 by cnh, Tue Jun 9 16:48:03 1998 UTC revision 1.27 by cnh, Wed Sep 26 18:09:16 2001 UTC
# Line 1  Line 1 
1  C $Header$  C $Header$
2    C $Name$
3    
4  #include "CPP_OPTIONS.h"  #include "CPP_OPTIONS.h"
5    
6  C     /==========================================================\  CBOP
7  C     | S/R TIMESTEP                                             |  C     !ROUTINE: TIMESTEP
8  C     | o Step model fields forward in time                      |  C     !INTERFACE:
9  C     \==========================================================/        SUBROUTINE TIMESTEP( bi, bj, iMin, iMax, jMin, jMax, K,
10        SUBROUTINE TIMESTEP( bi, bj, iMin, iMax, jMin, jMax,       I                     phiHyd, phiSurfX, phiSurfY,
11       I                     K,       I                     myIter, myThid )
12       I                     myThid )  C     !DESCRIPTION: \bv
13        implicit none  C     *==========================================================*
14  ! Common  C     | S/R TIMESTEP                                              
15    C     | o Step model fields forward in time                      
16    C     *==========================================================*
17    C     \ev
18    
19    C     !USES:
20          IMPLICIT NONE
21    C     == Global variables ==
22  #include "SIZE.h"  #include "SIZE.h"
23  #include "DYNVARS.h"  #include "DYNVARS.h"
24  #include "EEPARAMS.h"  #include "EEPARAMS.h"
25  #include "PARAMS.h"  #include "PARAMS.h"
26  #include "GRID.h"  #include "GRID.h"
27    #include "SURFACE.h"
28    
29    C     !INPUT/OUTPUT PARAMETERS:
30  C     == Routine Arguments ==  C     == Routine Arguments ==
31    C     phiHyd    - Hydrostatic Potential (ocean: pressure/rho)
32    C                                       (atmos: geopotentiel)
33    C     phiSurfX, - gradient of Surface potentiel (Pressure/rho, ocean)
34    C     phiSurfY             or geopotentiel (atmos) in X and Y direction
35        INTEGER bi,bj,iMin,iMax,jMin,jMax        INTEGER bi,bj,iMin,iMax,jMin,jMax
36        INTEGER K        INTEGER K
37        INTEGER myThid        _RL     phiHyd(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
38          _RL     phiSurfX(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
39          _RL     phiSurfY(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
40          INTEGER myIter, myThid
41    
42    C     !LOCAL VARIABLES:
43  C     == Local variables ==  C     == Local variables ==
44        INTEGER i,j        INTEGER i,j
45        _RL ab15,ab05        _RL ab15,ab05
46          _RL phxFac,phyFac, psFac
47          _RL     gUtmp(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
48          _RL     gVtmp(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
49    CEOP
50    
51  C     Adams-Bashforth timestepping weights  C     Adams-Bashforth timestepping weights
52        ab15=1.5+abeps        IF (myIter .EQ. 0) THEN
53        ab05=-0.5-abeps         ab15=1.0
54           ab05=0.0
55          ELSE
56           ab15=1.5+abeps
57           ab05=-0.5-abeps
58          ENDIF
59    
60  C     Step forward zonal velocity (store in Gu)  C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
61    C-    Compute effective gU term (including Adams-Bashforth weights) :
62          DO j=jMin,jMax
63           DO i=iMin,iMax
64            gUtmp(i,j) = ab15*gU(i,j,k,bi,bj)
65         &             + ab05*gUNm1(i,j,k,bi,bj)  
66    #ifdef INCLUDE_CD_CODE
67         &             + guCD(i,j,k,bi,bj)
68    #endif
69           ENDDO
70          ENDDO
71          
72    #ifdef NONLIN_FRSURF
73          IF (.NOT. vectorInvariantMomentum
74         &    .AND. nonlinFreeSurf.GT.1) THEN
75         DO j=jMin,jMax         DO j=jMin,jMax
76          DO i=iMin,iMax          DO i=iMin,iMax
77           gUNm1(i,j,k,bi,bj)=uVel(i,j,k,bi,bj)           IF ( k.EQ.ksurfW(i,j,bi,bj) ) THEN
78       &    +deltaTmom*(ab15*gU(i,j,k,bi,bj)+ab05*gUNm1(i,j,k,bi,bj)             gUtmp(i,j) = gUtmp(i,j)
79  #ifdef ALLOW_CD       &         *hFacW(i,j,k,bi,bj)/hFac_surfW(i,j,bi,bj)
80       &    +guCD(i,j,k,bi,bj)           ENDIF
 #endif  
      &        )*_maskW(i,j,k,bi,bj)  
81          ENDDO          ENDDO
82         ENDDO         ENDDO
83  C     Step forward meridional velocity (store in Gv)        ENDIF
84         DO j=jMin,jMax  #endif
85    
86    C     Step forward zonal velocity (store in Gu)
87          psFac = pfFacMom*(1. _d 0 - implicSurfPress)
88          DO j=jMin,jMax
89          DO i=iMin,iMax          DO i=iMin,iMax
90           gVNm1(i,j,k,bi,bj)=vVel(i,j,k,bi,bj)            gUNm1(i,j,k,bi,bj) = uVel(i,j,k,bi,bj)
91       &    +deltaTmom*(ab15*gV(i,j,k,bi,bj)+ab05*gVNm1(i,j,k,bi,bj)       &     +deltaTmom*(
92  #ifdef ALLOW_CD       &         gUtmp(i,j)
93       &    +gvCD(i,j,k,bi,bj)       &       - psFac*phiSurfX(i,j)
94  #endif       &        )*_maskW(i,j,k,bi,bj)
95       &        )*_maskS(i,j,k,bi,bj)          ENDDO
96          ENDDO
97    
98          IF (staggerTimeStep) THEN
99    C--   -grad Phi_Hyd has not been incorporated to gU and is added here:
100            phxFac = pfFacMom*deltaTmom
101            DO j=jMin,jMax
102              DO i=iMin,iMax
103                gUNm1(i,j,k,bi,bj)=gUNm1(i,j,k,bi,bj)
104         &       - _recip_dxC(i,j,bi,bj)
105         &         *(phiHyd(i,j,k)-phiHyd(i-1,j,k))*phxFac
106         &         *_maskW(i,j,k,bi,bj)
107              ENDDO
108          ENDDO          ENDDO
109          ENDIF
110    
111    C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
112    C-    Compute effective gV term (including Adams-Bashforth weights) :
113          DO j=jMin,jMax
114           DO i=iMin,iMax
115            gVtmp(i,j) = ab15*gV(i,j,k,bi,bj)
116         &             + ab05*gVNm1(i,j,k,bi,bj)  
117    #ifdef INCLUDE_CD_CODE
118         &             + gvCD(i,j,k,bi,bj)
119    #endif
120         ENDDO         ENDDO
121  C     Step forward temperature        ENDDO
122          
123    #ifdef NONLIN_FRSURF
124          IF (.NOT. vectorInvariantMomentum
125         &    .AND. nonlinFreeSurf.GT.1) THEN
126         DO j=jMin,jMax         DO j=jMin,jMax
127          DO i=iMin,iMax          DO i=iMin,iMax
128           theta(i,j,k,bi,bj)=theta(i,j,k,bi,bj)           IF ( k.EQ.ksurfS(i,j,bi,bj) ) THEN
129       &    +deltaTtracer*(ab15*gT(i,j,k,bi,bj)+ab05*gTNm1(i,j,k,bi,bj))             gVtmp(i,j) = gVtmp(i,j)
130           gTNm1(i,j,k,bi,bj)=gT(i,j,k,bi,bj)       &         *hFacS(i,j,k,bi,bj)/hFac_surfS(i,j,bi,bj)
131             ENDIF
132          ENDDO          ENDDO
133         ENDDO         ENDDO
134          ENDIF
135    #endif
136    
137    C     Step forward meridional velocity (store in Gv)
138          psFac = pfFacMom*(1. _d 0 - implicSurfPress)
139          DO j=jMin,jMax
140            DO i=iMin,iMax
141              gVNm1(i,j,k,bi,bj) = vVel(i,j,k,bi,bj)
142         &     +deltaTmom*(
143         &         gVtmp(i,j)
144         &       - psFac*phiSurfY(i,j)
145         &        )*_maskS(i,j,k,bi,bj)
146            ENDDO
147          ENDDO
148    
149          IF (staggerTimeStep) THEN
150    C--   -grad Phi_Hyd has not been incorporated to gV and is added here:
151            phyFac = pfFacMom*deltaTmom
152            DO j=jMin,jMax
153              DO i=iMin,iMax
154                gVNm1(i,j,k,bi,bj)=gVNm1(i,j,k,bi,bj)
155         &       - _recip_dyC(i,j,bi,bj)
156         &         *(phiHyd(i,j,k)-phiHyd(i,j-1,k))*phyFac
157         &         *_maskS(i,j,k,bi,bj)
158              ENDDO
159            ENDDO
160          ENDIF
161    
162        RETURN        RETURN
163        END        END

Legend:
Removed from v.1.11  
changed lines
  Added in v.1.27

  ViewVC Help
Powered by ViewVC 1.1.22