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

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

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


Revision 1.10 - (hide annotations) (download)
Mon Jun 22 15:26:25 1998 UTC (25 years, 11 months ago) by adcroft
Branch: MAIN
CVS Tags: checkpoint11, checkpoint10, checkpoint9, checkpoint8, checkpoint12
Changes since 1.9: +5 -3 lines
Various changes including time-dependant forcing:
 o logic for controlling external forcing fields now allows
   for time-dependant forcing: load_external_fields.F
 o genmake.dec needed a special line for the above file.
 o theta* and salt* time-stepping algorithm were re-implemented.
The 4x4 global configuration has been "double-checked" against
CNH's version. However, we do not assume any responsibility for
the correctness of this code ...  8-)

1 adcroft 1.10 C $Header: /u/gcmpack/models/MITgcmUV/model/src/calc_gs.F,v 1.9 1998/06/17 21:07:01 adcroft Exp $
2 cnh 1.1
3     #include "CPP_EEOPTIONS.h"
4    
5     CStartOfInterFace
6     SUBROUTINE CALC_GS(
7     I bi,bj,iMin,iMax,jMin,jMax,k,kM1,kUp,kDown,
8 adcroft 1.10 I xA,yA,uTrans,vTrans,wTrans,maskup,maskC,
9 adcroft 1.7 I K13,K23,KappaZS,KapGM,
10     U af,df,fZon,fMer,fVerS,
11 cnh 1.1 I myThid )
12     C /==========================================================\
13     C | SUBROUTINE CALC_GS |
14 adcroft 1.7 C | o Calculate the salt tendency terms. |
15 cnh 1.1 C |==========================================================|
16     C | A procedure called EXTERNAL_FORCING_S is called from |
17     C | here. These procedures can be used to add per problem |
18 adcroft 1.7 C | E-P flux source terms. |
19 cnh 1.1 C | Note: Although it is slightly counter-intuitive the |
20     C | EXTERNAL_FORCING routine is not the place to put |
21     C | file I/O. Instead files that are required to |
22     C | calculate the external source terms are generally |
23     C | read during the model main loop. This makes the |
24     C | logisitics of multi-processing simpler and also |
25     C | makes the adjoint generation simpler. It also |
26     C | allows for I/O to overlap computation where that |
27     C | is supported by hardware. |
28     C | Aside from the problem specific term the code here |
29     C | forms the tendency terms due to advection and mixing |
30     C | The baseline implementation here uses a centered |
31     C | difference form for the advection term and a tensorial |
32     C | divergence of a flux form for the diffusive term. The |
33     C | diffusive term is formulated so that isopycnal mixing and|
34     C | GM-style subgrid-scale terms can be incorporated b simply|
35     C | setting the diffusion tensor terms appropriately. |
36     C \==========================================================/
37     IMPLICIT NONE
38    
39     C == GLobal variables ==
40     #include "SIZE.h"
41     #include "DYNVARS.h"
42     #include "EEPARAMS.h"
43     #include "PARAMS.h"
44     #include "GRID.h"
45 cnh 1.8 #include "FFIELDS.h"
46 cnh 1.1
47     C == Routine arguments ==
48     C fZon - Work array for flux of temperature in the east-west
49     C direction at the west face of a cell.
50     C fMer - Work array for flux of temperature in the north-south
51     C direction at the south face of a cell.
52 adcroft 1.7 C fVerS - Flux of salt (S) in the vertical
53 cnh 1.1 C direction at the upper(U) and lower(D) faces of a cell.
54     C maskUp - Land mask used to denote base of the domain.
55 adcroft 1.10 C maskC - Land mask for salt cells (used in TOP_LAYER only)
56 cnh 1.1 C xA - Tracer cell face area normal to X
57     C yA - Tracer cell face area normal to X
58     C uTrans - Zonal volume transport through cell face
59     C vTrans - Meridional volume transport through cell face
60     C wTrans - Vertical volume transport through cell face
61     C af - Advective flux component work array
62     C df - Diffusive flux component work array
63     C bi, bj, iMin, iMax, jMin, jMax - Range of points for which calculation
64     C results will be set.
65 adcroft 1.7 C myThid - Instance number for this innvocation of CALC_GT
66 cnh 1.1 _RL fZon (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
67     _RL fMer (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
68     _RL fVerS (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)
69     _RS xA (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
70     _RS yA (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
71     _RL uTrans(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
72     _RL vTrans(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
73     _RL wTrans(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
74     _RS maskUp(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
75 adcroft 1.10 _RS maskC (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
76 adcroft 1.7 _RL K13 (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nz)
77     _RL K23 (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nz)
78     _RL KappaZS(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nz)
79     _RL KapGM (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
80 cnh 1.1 _RL af (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
81     _RL df (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
82 adcroft 1.7 INTEGER k,kUp,kDown,kM1
83 cnh 1.1 INTEGER bi,bj,iMin,iMax,jMin,jMax
84     INTEGER myThid
85     CEndOfInterface
86    
87     C == Local variables ==
88     C I, J, K - Loop counters
89 adcroft 1.7 INTEGER i,j
90     LOGICAL TOP_LAYER
91     _RL afFacS, dfFacS
92     _RL dSdx(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
93     _RL dSdy(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
94 cnh 1.1
95     afFacS = 1. _d 0
96     dfFacS = 1. _d 0
97 adcroft 1.7 TOP_LAYER = K .EQ. 1
98 cnh 1.1
99     C--- Calculate advective and diffusive fluxes between cells.
100    
101     C-- Zonal flux (fZon is at west face of "salt" cell)
102     C Advective component of zonal flux
103     DO j=jMin,jMax
104     DO i=iMin,iMax
105     af(i,j) =
106     & uTrans(i,j)*(salt(i,j,k,bi,bj)+salt(i-1,j,k,bi,bj))*0.5 _d 0
107     ENDDO
108     ENDDO
109 adcroft 1.7 C Zonal tracer gradient
110     DO j=jMin,jMax
111     DO i=iMin,iMax
112     dSdx(i,j) = _rdxC(i,j,bi,bj)*
113     & (salt(i,j,k,bi,bj)-salt(i-1,j,k,bi,bj))
114     ENDDO
115     ENDDO
116 cnh 1.1 C Diffusive component of zonal flux
117     DO j=jMin,jMax
118     DO i=iMin,iMax
119 adcroft 1.7 df(i,j) = -(diffKhS+0.5*(KapGM(i,j)+KapGM(i-1,j)))*
120     & xA(i,j)*dSdx(i,j)
121 cnh 1.1 ENDDO
122     ENDDO
123     C Net zonal flux
124     DO j=jMin,jMax
125     DO i=iMin,iMax
126     fZon(i,j) = afFacS*af(i,j) + dfFacS*df(i,j)
127     ENDDO
128     ENDDO
129    
130     C-- Meridional flux (fMer is at south face of "salt" cell)
131     C Advective component of meridional flux
132     DO j=jMin,jMax
133     DO i=iMin,iMax
134     C Advective component of meridional flux
135     af(i,j) =
136     & vTrans(i,j)*(salt(i,j,k,bi,bj)+salt(i,j-1,k,bi,bj))*0.5 _d 0
137     ENDDO
138     ENDDO
139 adcroft 1.7 C Zonal tracer gradient
140     DO j=jMin,jMax
141     DO i=iMin,iMax
142     dSdy(i,j) = _rdyC(i,j,bi,bj)*
143     & (salt(i,j,k,bi,bj)-salt(i,j-1,k,bi,bj))
144     ENDDO
145     ENDDO
146 cnh 1.1 C Diffusive component of meridional flux
147     DO j=jMin,jMax
148     DO i=iMin,iMax
149 adcroft 1.7 df(i,j) = -(diffKhS+0.5*(KapGM(i,j)+KapGM(i,j-1)))*
150     & yA(i,j)*dSdy(i,j)
151 cnh 1.1 ENDDO
152     ENDDO
153     C Net meridional flux
154     DO j=jMin,jMax
155     DO i=iMin,iMax
156     fMer(i,j) = afFacS*af(i,j) + dfFacS*df(i,j)
157     ENDDO
158     ENDDO
159    
160 adcroft 1.7 C-- Interpolate terms for Redi/GM scheme
161     DO j=jMin,jMax
162     DO i=iMin,iMax
163     dSdx(i,j) = 0.5*(
164     & +0.5*(_maskW(i+1,j,k,bi,bj)*_rdxC(i+1,j,bi,bj)*
165     & (salt(i+1,j,k,bi,bj)-salt(i,j,k,bi,bj))
166     & +_maskW(i,j,k,bi,bj)*_rdxC(i,j,bi,bj)*
167     & (salt(i,j,k,bi,bj)-salt(i-1,j,k,bi,bj)))
168     & +0.5*(_maskW(i+1,j,km1,bi,bj)*_rdxC(i+1,j,bi,bj)*
169     & (salt(i+1,j,km1,bi,bj)-salt(i,j,km1,bi,bj))
170     & +_maskW(i,j,km1,bi,bj)*_rdxC(i,j,bi,bj)*
171     & (salt(i,j,km1,bi,bj)-salt(i-1,j,km1,bi,bj)))
172     & )
173     ENDDO
174     ENDDO
175     DO j=jMin,jMax
176     DO i=iMin,iMax
177     dSdy(i,j) = 0.5*(
178     & +0.5*(_maskS(i,j,k,bi,bj)*_rdyC(i,j,bi,bj)*
179     & (salt(i,j,k,bi,bj)-salt(i,j-1,k,bi,bj))
180     & +_maskS(i,j+1,k,bi,bj)*_rdyC(i,j+1,bi,bj)*
181     & (salt(i,j+1,k,bi,bj)-salt(i,j,k,bi,bj)))
182     & +0.5*(_maskS(i,j,km1,bi,bj)*_rdyC(i,j,bi,bj)*
183     & (salt(i,j,km1,bi,bj)-salt(i,j-1,km1,bi,bj))
184     & +_maskS(i,j+1,km1,bi,bj)*_rdyC(i,j+1,bi,bj)*
185     & (salt(i,j+1,km1,bi,bj)-salt(i,j,km1,bi,bj)))
186     & )
187     ENDDO
188     ENDDO
189    
190 cnh 1.1 C-- Vertical flux (fVerS) above
191     C Advective component of vertical flux
192 adcroft 1.7 C Note: For K=1 then KM1=1 this gives a barZ(T) = T
193     C (this plays the role of the free-surface correction)
194 cnh 1.1 DO j=jMin,jMax
195     DO i=iMin,iMax
196     af(i,j) =
197     & wTrans(i,j)*(salt(i,j,k,bi,bj)+salt(i,j,kM1,bi,bj))*0.5 _d 0
198     ENDDO
199     ENDDO
200     C Diffusive component of vertical flux
201 adcroft 1.7 C Note: For K=1 then KM1=1 this gives a dS/dz = 0 upper
202     C boundary condition.
203 cnh 1.1 DO j=jMin,jMax
204     DO i=iMin,iMax
205 adcroft 1.7 df(i,j) = _zA(i,j,bi,bj)*(
206     & -KapGM(i,j)*K13(i,j,k)*dSdx(i,j)
207     & -KapGM(i,j)*K23(i,j,k)*dSdy(i,j)
208     & )
209 cnh 1.1 ENDDO
210     ENDDO
211 adcroft 1.7 IF (.NOT.implicitDiffusion) THEN
212     DO j=jMin,jMax
213     DO i=iMin,iMax
214     df(i,j) = df(i,j) + _zA(i,j,bi,bj)*(
215     & -KappaZS(i,j,k)*rdzC(k)
216     & *(salt(i,j,kM1,bi,bj)-salt(i,j,k,bi,bj))
217     & )
218     ENDDO
219     ENDDO
220     ENDIF
221 cnh 1.1 C Net vertical flux
222     DO j=jMin,jMax
223     DO i=iMin,iMax
224 adcroft 1.7 fVerS(i,j,kUp) = ( afFacS*af(i,j)+ dfFacS*df(i,j) )*maskUp(i,j)
225 cnh 1.1 ENDDO
226     ENDDO
227 adcroft 1.7 IF ( TOP_LAYER ) THEN
228     DO j=jMin,jMax
229     DO i=iMin,iMax
230     fVerS(i,j,kUp) = afFacS*af(i,j)*freeSurfFac
231     ENDDO
232     ENDDO
233     ENDIF
234 cnh 1.1
235     C-- Tendency is minus divergence of the fluxes.
236     C Note. Tendency terms will only be correct for range
237     C i=iMin+1:iMax-1, j=jMin+1:jMax-1. Edge points
238     C will contain valid floating point numbers but
239     C they are not algorithmically correct. These points
240     C are not used.
241     DO j=jMin,jMax
242     DO i=iMin,iMax
243 adcroft 1.7 C & -_rhFacC(i,j,k,bi,bj)*rdzF(k)*_rdxF(i,j,bi,bj)*_rdyF(i,j,bi,bj)
244     C & -_rhFacC(i,j,k,bi,bj)*rdzF(k)/_zA(i,j,bi,bj)
245     C #define _rVolS(i,j,k,bi,bj) _rhFacC(i,j,k,bi,bj)*rdzF(k)*_rdxF(i,j,bi,bj)*_rdyF(i,j,bi,bj)
246     #define _rVolS(i,j,k,bi,bj) _rhFacC(i,j,k,bi,bj)*rdzF(k)/_zA(i,j,bi,bj)
247 cnh 1.1 gS(i,j,k,bi,bj)=
248 adcroft 1.7 & -_rVolS(i,j,k,bi,bj)
249 cnh 1.1 & *(
250     & +( fZon(i+1,j)-fZon(i,j) )
251     & +( fMer(i,j+1)-fMer(i,j) )
252     & +( fVerS(i,j,kUp)-fVerS(i,j,kDown) )
253     & )
254     ENDDO
255     ENDDO
256    
257 adcroft 1.7 C-- External P-E forcing term(s)
258 cnh 1.8 C o Surface relaxation term
259     IF ( TOP_LAYER ) THEN
260     DO j=jMin,jMax
261     DO i=iMin,iMax
262     gS(i,j,k,bi,bj)=gS(i,j,k,bi,bj)
263 adcroft 1.10 & +maskC(i,j)*(
264 adcroft 1.9 & -lambdaSaltClimRelax*(salt(i,j,k,bi,bj)-SSS(i,j,bi,bj))
265     & -EmPpR(i,j,bi,bj) )
266 cnh 1.8 ENDDO
267     ENDDO
268     ENDIF
269    
270 cnh 1.1
271     RETURN
272     END

  ViewVC Help
Powered by ViewVC 1.1.22