/[MITgcm]/MITgcm/pkg/mom_fluxform/mom_cdscheme.F
ViewVC logotype

Annotation of /MITgcm/pkg/mom_fluxform/mom_cdscheme.F

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

Revision 1.2 - (hide annotations) (download)
Tue May 29 14:01:38 2001 UTC (23 years ago) by adcroft
Branch: MAIN
CVS Tags: checkpoint40pre1, checkpoint40pre9, checkpoint40pre7, checkpoint40pre2, checkpoint40pre5, checkpoint40pre6, checkpoint40pre8, checkpoint40pre4, checkpoint40pre3, checkpoint40
Changes since 1.1: +186 -0 lines 
Merge from branch pre38:
 o essential mods for cubed sphere
 o debugged atmosphere, dynamcis + physics (aim)
 o new packages (mom_vecinv, mom_fluxform, ...)
1 adcroft 1.2 C $Header: /u/gcmpack/models/MITgcmUV/pkg/mom_fluxform/Attic/mom_cdscheme.F,v 1.1.2.1 2001/03/28 19:51:14 adcroft Exp $
2     C $Name: pre38-close $
3    
4     #include "CPP_OPTIONS.h"
5    
6     SUBROUTINE MOM_CDSCHEME(
7     I bi,bj,k,phi_hyd,
8     I myThid)
9    
10     IMPLICIT NONE
11    
12     C == Global variables ==
13     #include "SIZE.h"
14     #include "DYNVARS.h"
15     #include "EEPARAMS.h"
16     #include "PARAMS.h"
17     #include "GRID.h"
18     #include "SURFACE.h"
19    
20     C == Routine arguments ==
21     C myThid - Instance number for this innvocation of CALC_MOM_RHS
22     INTEGER bi,bj,K
23     _RL phi_hyd(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
24     INTEGER myThid
25    
26     #ifdef INCLUDE_CD_CODE
27    
28     C == Local variables ==
29     _RL pF(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
30     _RL vF(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
31     _RL aF(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
32     INTEGER i,j,iMin,iMax,jMin,jMax
33     _RL ab15,ab05
34    
35     C Compute ranges
36     iMin=1-Olx+1
37     iMax=sNx+Olx-1
38     jMin=1-Oly+1
39     jMax=sNy+Oly-1
40    
41     C Adams-Bashforth weighting factors
42     ab15 = 1.5 + abEps
43     ab05 = -0.5 - abEps
44    
45     C Pressure extrapolated forward in time
46     DO j=1-Oly,sNy+Oly
47     DO i=1-Olx,sNx+Olx
48     pf(i,j) =
49     & ab15*( etaN(i,j,bi,bj)*Bo_surf(i,j,bi,bj) )
50     & +ab05*(etaNm1(i,j,bi,bj)*Bo_surf(i,j,bi,bj) )
51     ENDDO
52     ENDDO
53     IF (staggerTimeStep) THEN
54     DO j=jMin,jMax
55     DO i=iMin,iMax
56     pf(i,j) = pf(i,j)+phi_hyd(i,j,k)
57     ENDDO
58     ENDDO
59     ENDIF
60    
61     C-- Zonal velocity coriolis term
62     C Note. As coded here, coriolis will not work with "thin walls"
63     C-- Coriolis with CD scheme allowed
64     C grady(p) + gV
65     DO j=1-Oly+1,sNy+Oly
66     DO i=1-Olx,sNx+Olx
67     af(i,j) = -_maskS(i,j,k,bi,bj)
68     & *_recip_dyC(i,j,bi,bj)
69     & *(pf(i,j)-pf(i,j-1))
70     & +gV(i,j,k,bi,bj)
71     ENDDO
72     ENDDO
73     C Average to Vd point and add coriolis
74     DO j=jMin,jMax
75     DO i=iMin,iMax
76     vf(i,j) =
77     & 0.25*( af(i ,j)+af(i ,j+1)
78     & +af(i-1,j)+af(i-1,j+1)
79     & )*_maskW(i,j,k,bi,bj)
80     & -0.5*(_fCori(i,j,bi,bj)+
81     & _fCori(i-1,j,bi,bj))
82     & *uVel(i,j,k,bi,bj)
83     ENDDO
84     ENDDO
85     C Step forward Vd
86     DO j=jMin,jMax
87     DO i=iMin,iMax
88     vVelD(i,j,k,bi,bj) = vVelD(i,j,k,bi,bj) +
89     & deltaTmom*vf(i,j)
90     ENDDO
91     ENDDO
92     C Relax D grid V to C grid V
93     DO j=jMin,jMax
94     DO i=iMin,iMax
95     vVelD(i,j,k,bi,bj) = rCD*vVelD(i,j,k,bi,bj)
96     & +(1. - rCD)*(
97     & ab15*0.25*(
98     & vVel(i ,j ,k,bi,bj)+vVel(i ,j+1,k,bi,bj)
99     & +vVel(i-1,j ,k,bi,bj)+vVel(i-1,j+1,k,bi,bj)
100     & )*_maskW(i,j,k,bi,bj)
101     & +
102     & ab05*0.25*(
103     & vNM1(i ,j ,k,bi,bj)+vNM1(i ,j+1,k,bi,bj)
104     & +vNM1(i-1,j ,k,bi,bj)+vNM1(i-1,j+1,k,bi,bj)
105     & )*_maskW(i,j,k,bi,bj)
106     & )
107     ENDDO
108     ENDDO
109     C Calculate coriolis force on U
110     DO j=jMin,jMax
111     DO i=iMin,iMax
112     guCD(i,j,k,bi,bj) =
113     & 0.5*( _fCori(i ,j,bi,bj) +
114     & _fCori(i-1,j,bi,bj) )
115     & *vVelD(i,j,k,bi,bj)*cfFacMom
116     ENDDO
117     ENDDO
118    
119     C-- Meridional velocity coriolis term
120     C gradx(p)+gU
121     DO j=1-Oly,sNy+Oly
122     DO i=1-Olx+1,sNx+Olx
123     af(i,j) = -_maskW(i,j,k,bi,bj)
124     & *_recip_dxC(i,j,bi,bj)*
125     & (pf(i,j)-pf(i-1,j))
126     & +gU(i,j,k,bi,bj)
127     ENDDO
128     ENDDO
129     C Average to Ud point and add coriolis
130     DO j=jMin,jMax
131     DO i=iMin,iMax
132     vf(i,j) =
133     & 0.25*( af(i ,j)+af(i ,j-1)
134     & +af(i+1,j)+af(i+1,j-1)
135     & )*_maskS(i,j,k,bi,bj)
136     & +0.5*( _fCori(i,j,bi,bj)
137     & +_fCori(i,j-1,bi,bj))
138     & *vVel(i,j,k,bi,bj)
139     ENDDO
140     ENDDO
141     C Step forward Ud
142     DO j=jMin,jMax
143     DO i=iMin,iMax
144     uVelD(i,j,k,bi,bj) = uVelD(i,j,k,bi,bj) +
145     & deltaTmom*vf(i,j)*_maskS(i,j,k,bi,bj)
146     ENDDO
147     ENDDO
148     C Relax D grid U to C grid U
149     DO j=jMin,jMax
150     DO i=iMin,iMax
151     uVelD(i,j,k,bi,bj) = rCD*uVelD(i,j,k,bi,bj)
152     & +(1. - rCD)*(
153     & ab15*0.25*(
154     & uVel(i,j ,k,bi,bj)+uVel(i+1,j ,k,bi,bj)
155     & +uVel(i,j-1,k,bi,bj)+uVel(i+1,j-1,k,bi,bj)
156     & )*_maskS(i,j,k,bi,bj)
157     & +
158     & ab05*0.25*(
159     & uNM1(i,j ,k,bi,bj)+uNM1(i+1,j ,k,bi,bj)
160     & +uNM1(i,j-1,k,bi,bj)+uNM1(i+1,j-1,k,bi,bj)
161     & )*_maskS(i,j,k,bi,bj)
162     & )
163     ENDDO
164     ENDDO
165     C Calculate coriolis force on V
166     DO j=jMin,jMax
167     DO i=iMin,iMax
168     gvCD(i,j,k,bi,bj) =
169     & -0.5*( _fCori(i ,j,bi,bj)
170     & +_fCori(i,j-1,bi,bj) )
171     & *uVelD(i,j,k,bi,bj)*_maskS(i,j,k,bi,bj)*cfFacMom
172     ENDDO
173     ENDDO
174    
175     C-- Save "previous time level" variables
176     DO j=1-OLy,sNy+OLy
177     DO i=1-OLx,sNx+OLx
178     uNM1(i,j,k,bi,bj) = uVel(i,j,k,bi,bj)
179     vNM1(i,j,k,bi,bj) = vVel(i,j,k,bi,bj)
180     ENDDO
181     ENDDO
182    
183     #endif /* INCLUDE_CD_CODE */
184    
185     RETURN
186     END
  ViewVC Help
Powered by ViewVC 1.1.22