/[MITgcm]/MITgcm/verification/hs94.cs-32x32x5/code/external_forcing.F
ViewVC logotype

Annotation of /MITgcm/verification/hs94.cs-32x32x5/code/external_forcing.F

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


Revision 1.7 - (hide annotations) (download)
Fri Sep 24 20:43:35 2010 UTC (13 years, 7 months ago) by jmc
Branch: MAIN
CVS Tags: checkpoint65b, checkpoint65a, checkpoint62o, checkpoint62n, checkpoint62m, checkpoint62l, checkpoint62s, checkpoint62r, checkpoint62q, checkpoint62p, checkpoint62w, checkpoint62v, checkpoint62u, checkpoint62t, checkpoint62z, checkpoint62y, checkpoint62x, checkpoint63i, checkpoint64, checkpoint65, checkpoint63, checkpoint63p, checkpoint63q, checkpoint63r, checkpoint63s, checkpoint63l, checkpoint63m, checkpoint63n, checkpoint63o, checkpoint63h, checkpoint63j, checkpoint63k, checkpoint63d, checkpoint63e, checkpoint63f, checkpoint63g, checkpoint63a, checkpoint63b, checkpoint63c, checkpoint64y, checkpoint64x, checkpoint64z, checkpoint64q, checkpoint64p, checkpoint64s, checkpoint64r, checkpoint64u, checkpoint64t, checkpoint64w, checkpoint64v, checkpoint64i, checkpoint64h, checkpoint64k, checkpoint64j, checkpoint64m, checkpoint64l, checkpoint64o, checkpoint64n, checkpoint64a, checkpoint64c, checkpoint64b, checkpoint64e, checkpoint64d, checkpoint64g, checkpoint64f
Changes since 1.6: +63 -15 lines
compatible with sigma (and hybrid-sigma) coordinate

1 jmc 1.7 C $Header: /u/gcmpack/MITgcm/verification/hs94.cs-32x32x5/code/external_forcing.F,v 1.6 2005/07/17 18:40:10 jmc Exp $
2 adcroft 1.3 C $Name: $
3 adcroft 1.2
4     #include "CPP_OPTIONS.h"
5    
6 jmc 1.6 CBOP
7     C !ROUTINE: EXTERNAL_FORCING_U
8     C !INTERFACE:
9 adcroft 1.2 SUBROUTINE EXTERNAL_FORCING_U(
10 jmc 1.6 I iMin,iMax, jMin,jMax, bi,bj, kLev,
11     I myTime, myThid )
12     C !DESCRIPTION: \bv
13     C *==========================================================*
14     C | S/R EXTERNAL_FORCING_U
15     C | o Contains problem specific forcing for zonal velocity.
16     C *==========================================================*
17     C | Adds terms to gU for forcing by external sources
18     C | e.g. wind stress, bottom friction etc ...
19     C *==========================================================*
20     C \ev
21    
22     C !USES:
23 adcroft 1.2 IMPLICIT NONE
24     C == Global data ==
25     #include "SIZE.h"
26     #include "EEPARAMS.h"
27     #include "PARAMS.h"
28     #include "GRID.h"
29 jmc 1.7 #include "SURFACE.h"
30 adcroft 1.2 #include "DYNVARS.h"
31     #include "FFIELDS.h"
32    
33 jmc 1.6 C !INPUT/OUTPUT PARAMETERS:
34 adcroft 1.2 C == Routine arguments ==
35 jmc 1.6 C iMin,iMax :: Working range of x-index for applying forcing.
36     C jMin,jMax :: Working range of y-index for applying forcing.
37     C bi,bj :: Current tile indices
38     C kLev :: Current vertical level index
39     C myTime :: Current time in simulation
40     C myThid :: Thread Id number
41 adcroft 1.2 INTEGER iMin, iMax, jMin, jMax, kLev, bi, bj
42 jmc 1.6 _RL myTime
43 adcroft 1.2 INTEGER myThid
44    
45 jmc 1.6 C !LOCAL VARIABLES:
46 adcroft 1.2 C == Local variables ==
47 jmc 1.6 C i,j :: Loop counters
48     INTEGER i, j
49     CEOP
50 jmc 1.7 _RL recip_P0g, termP, rFullDepth
51     _RL kV, kF, sigma_b
52 adcroft 1.2
53     C-- Forcing term(s)
54 adcroft 1.3 kF=1. _d 0/86400. _d 0
55 jmc 1.5 sigma_b = 0.7 _d 0
56 jmc 1.7 rFullDepth = rF(1)-rF(Nr+1)
57 jmc 1.6 c DO j=1,sNy
58     C-jmc: Without CD-scheme, this is OK ; but with CD-scheme, needs to cover [0:sNy+1]
59     DO j=0,sNy+1
60     DO i=1,sNx+1
61 jmc 1.5 IF ( hFacW(i,j,kLev,bi,bj) .GT. 0. ) THEN
62 jmc 1.7 IF ( selectSigmaCoord.EQ.0 ) THEN
63     recip_P0g=MAX(recip_Rcol(i,j,bi,bj),recip_Rcol(i-1,j,bi,bj))
64     termP=0.5 _d 0*( MIN(rF(kLev)*recip_P0g,1. _d 0)
65     & +rF(kLev+1)*recip_P0g )
66     ELSE
67     C-- Pressure at U.point :
68     c midP = rLowW(i,j,bi,bj) + aHybSigmC(k)*rFullDepth
69     c & + bHybSigmC(k)
70     c & *(etaHw(i,j,bi,bj)+rSurfW(i,j,bi,bj)-rLowW(i,j,bi,bj))
71     C-- Sigma at U.point :
72     c termP = ( midP - rLowW(i,j,bi,bj))
73     c & /(etaHw(i,j,bi,bj)+rSurfW(i,j,bi,bj)-rLowW(i,j,bi,bj))
74     C- which simplifies to:
75     termP = aHybSigmC(kLev)*rFullDepth
76     & /(etaHw(i,j,bi,bj)+rSurfW(i,j,bi,bj)-rLowW(i,j,bi,bj))
77     & + bHybSigmC(kLev)
78     ENDIF
79 jmc 1.5 kV=kF*MAX( 0. _d 0, (termP-sigma_b)/(1. _d 0-sigma_b) )
80 adcroft 1.2 gU(i,j,kLev,bi,bj)=gU(i,j,kLev,bi,bj)
81 jmc 1.6 & -kV*uVel(i,j,kLev,bi,bj)
82 adcroft 1.2 ENDIF
83     ENDDO
84     ENDDO
85    
86     RETURN
87     END
88 jmc 1.6
89     C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
90     CBOP
91     C !ROUTINE: EXTERNAL_FORCING_V
92     C !INTERFACE:
93 adcroft 1.2 SUBROUTINE EXTERNAL_FORCING_V(
94 jmc 1.6 I iMin,iMax, jMin,jMax, bi,bj, kLev,
95     I myTime, myThid )
96     C !DESCRIPTION: \bv
97     C *==========================================================*
98     C | S/R EXTERNAL_FORCING_V
99     C | o Contains problem specific forcing for merid velocity.
100     C *==========================================================*
101     C | Adds terms to gV for forcing by external sources
102     C | e.g. wind stress, bottom friction etc ...
103     C *==========================================================*
104     C \ev
105    
106     C !USES:
107 adcroft 1.2 IMPLICIT NONE
108     C == Global data ==
109     #include "SIZE.h"
110     #include "EEPARAMS.h"
111     #include "PARAMS.h"
112     #include "GRID.h"
113 jmc 1.7 #include "SURFACE.h"
114 adcroft 1.2 #include "DYNVARS.h"
115     #include "FFIELDS.h"
116    
117 jmc 1.6 C !INPUT/OUTPUT PARAMETERS:
118 adcroft 1.2 C == Routine arguments ==
119 jmc 1.6 C iMin,iMax :: Working range of x-index for applying forcing.
120     C jMin,jMax :: Working range of y-index for applying forcing.
121     C bi,bj :: Current tile indices
122     C kLev :: Current vertical level index
123     C myTime :: Current time in simulation
124     C myThid :: Thread Id number
125 adcroft 1.2 INTEGER iMin, iMax, jMin, jMax, kLev, bi, bj
126 jmc 1.6 _RL myTime
127 adcroft 1.2 INTEGER myThid
128 jmc 1.6
129     C !LOCAL VARIABLES:
130 adcroft 1.2 C == Local variables ==
131 jmc 1.6 C i,j :: Loop counters
132     INTEGER i, j
133     CEOP
134 jmc 1.7 _RL recip_P0g, termP, rFullDepth
135     _RL kV, kF, sigma_b
136 adcroft 1.2
137     C-- Forcing term(s)
138 adcroft 1.3 kF=1. _d 0/86400. _d 0
139 jmc 1.5 sigma_b = 0.7 _d 0
140 jmc 1.7 rFullDepth = rF(1)-rF(Nr+1)
141 jmc 1.6 DO j=1,sNy+1
142     c DO i=1,sNx
143     C-jmc: Without CD-scheme, this is OK ; but with CD-scheme, needs to cover [0:sNx+1]
144     DO i=0,sNx+1
145 jmc 1.5 IF ( hFacS(i,j,kLev,bi,bj) .GT. 0. ) THEN
146 jmc 1.7 IF ( selectSigmaCoord.EQ.0 ) THEN
147     recip_P0g=MAX(recip_Rcol(i,j,bi,bj),recip_Rcol(i,j-1,bi,bj))
148     termP=0.5 _d 0*( MIN(rF(kLev)*recip_P0g,1. _d 0)
149     & +rF(kLev+1)*recip_P0g )
150     ELSE
151     C-- Pressure at V.point :
152     c midP = rLowS(i,j,bi,bj) + aHybSigmC(k)*rFullDepth
153     c & + bHybSigmC(k)
154     c & *(etaHs(i,j,bi,bj)+rSurfS(i,j,bi,bj)-rLowS(i,j,bi,bj))
155     C-- Sigma at V.point :
156     c termP = ( midP - rLowS(i,j,bi,bj))
157     c & /(etaHs(i,j,bi,bj)+rSurfS(i,j,bi,bj)-rLowS(i,j,bi,bj))
158     C- which simplifies to:
159     termP = aHybSigmC(kLev)*rFullDepth
160     & /(etaHs(i,j,bi,bj)+rSurfS(i,j,bi,bj)-rLowS(i,j,bi,bj))
161     & + bHybSigmC(kLev)
162     ENDIF
163 jmc 1.5 kV=kF*MAX( 0. _d 0, (termP-sigma_b)/(1. _d 0-sigma_b) )
164 adcroft 1.2 gV(i,j,kLev,bi,bj)=gV(i,j,kLev,bi,bj)
165     & -kV*vVel(i,j,kLev,bi,bj)
166     ENDIF
167     ENDDO
168     ENDDO
169    
170     RETURN
171     END
172 jmc 1.6
173     C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
174     CBOP
175     C !ROUTINE: EXTERNAL_FORCING_T
176     C !INTERFACE:
177 adcroft 1.2 SUBROUTINE EXTERNAL_FORCING_T(
178 jmc 1.6 I iMin,iMax, jMin,jMax, bi,bj, kLev,
179     I myTime, myThid )
180     C !DESCRIPTION: \bv
181     C *==========================================================*
182     C | S/R EXTERNAL_FORCING_T
183     C | o Contains problem specific forcing for temperature.
184     C *==========================================================*
185     C | Adds terms to gT for forcing by external sources
186     C | e.g. heat flux, climatalogical relaxation, etc ...
187     C *==========================================================*
188     C \ev
189    
190     C !USES:
191 adcroft 1.2 IMPLICIT NONE
192     C == Global data ==
193     #include "SIZE.h"
194     #include "EEPARAMS.h"
195     #include "PARAMS.h"
196     #include "GRID.h"
197     #include "DYNVARS.h"
198     #include "FFIELDS.h"
199    
200 jmc 1.6 C !INPUT/OUTPUT PARAMETERS:
201 adcroft 1.2 C == Routine arguments ==
202 jmc 1.6 C iMin,iMax :: Working range of x-index for applying forcing.
203     C jMin,jMax :: Working range of y-index for applying forcing.
204     C bi,bj :: Current tile indices
205     C kLev :: Current vertical level index
206     C myTime :: Current time in simulation
207     C myThid :: Thread Id number
208 adcroft 1.2 INTEGER iMin, iMax, jMin, jMax, kLev, bi, bj
209 jmc 1.6 _RL myTime
210 adcroft 1.2 INTEGER myThid
211    
212 jmc 1.6 C !LOCAL VARIABLES:
213 adcroft 1.2 C == Local variables ==
214 jmc 1.6 C i,j :: Loop counters
215     INTEGER i, j
216     CEOP
217 jmc 1.7 _RL thetaLim,kT,ka,ks,sigma_b,term1,term2,thetaEq
218     _RL termP, rFullDepth
219 adcroft 1.2
220     C-- Forcing term(s)
221 adcroft 1.3 ka=1. _d 0/(40. _d 0*86400. _d 0)
222     ks=1. _d 0/(4. _d 0 *86400. _d 0)
223 jmc 1.5 sigma_b = 0.7 _d 0
224 jmc 1.7 rFullDepth = rF(1)-rF(Nr+1)
225 jmc 1.6 DO j=1,sNy
226     DO i=1,sNx
227     term1=60. _d 0*(SIN(yC(i,j,bi,bj)*deg2rad)**2)
228 jmc 1.5 termP=0.5 _d 0*( rF(kLev) + rF(kLev+1) )
229 jmc 1.6 term2=10. _d 0*LOG(termP/atm_po)
230     & *(COS(yC(i,j,bi,bj)*deg2rad)**2)
231 jmc 1.5 thetaLim = 200. _d 0/ ((termP/atm_po)**atm_kappa)
232     thetaEq=315. _d 0-term1-term2
233     thetaEq=MAX(thetaLim,thetaEq)
234 jmc 1.7 IF ( selectSigmaCoord.EQ.0 ) THEN
235     termP=0.5 _d 0*( MIN(rF(kLev),Ro_surf(i,j,bi,bj))+rF(kLev+1) )
236     & *recip_Rcol(i,j,bi,bj)
237     ELSE
238     C-- Pressure at T.point :
239     c midP = R_low(i,j,bi,bj) + aHybSigmC(k)*rFullDepth
240     c & + bHybSigmC(k)
241     c & *(etaH(i,j,bi,bj)+Ro_surf(i,j,bi,bj)-R_low(i,j,bi,bj))
242     C-- Sigma at T.point :
243     c termP = ( midP - R_low(i,j,bi,bj))
244     c & /(etaH(i,j,bi,bj)+Ro_surf(i,j,bi,bj)-R_low(i,j,bi,bj))
245     C- which simplifies to:
246     termP = aHybSigmC(kLev)*rFullDepth
247     & /(etaH(i,j,bi,bj)+Ro_surf(i,j,bi,bj)-R_low(i,j,bi,bj))
248     & + bHybSigmC(kLev)
249     ENDIF
250 jmc 1.5 kT=ka+(ks-ka)
251 jmc 1.7 & *MAX(0. _d 0, (termP-sigma_b)/(1. _d 0-sigma_b) )
252 jmc 1.6 & *COS((yC(i,j,bi,bj)*deg2rad))**4
253 adcroft 1.2 gT(i,j,kLev,bi,bj)=gT(i,j,kLev,bi,bj)
254 jmc 1.6 & - kT*( theta(i,j,kLev,bi,bj)-thetaEq )
255 adcroft 1.2 & *maskC(i,j,kLev,bi,bj)
256     ENDDO
257     ENDDO
258    
259     RETURN
260     END
261 jmc 1.6
262     C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
263     CBOP
264     C !ROUTINE: EXTERNAL_FORCING_S
265     C !INTERFACE:
266 adcroft 1.2 SUBROUTINE EXTERNAL_FORCING_S(
267 jmc 1.6 I iMin,iMax, jMin,jMax, bi,bj, kLev,
268     I myTime, myThid )
269    
270     C !DESCRIPTION: \bv
271     C *==========================================================*
272     C | S/R EXTERNAL_FORCING_S
273     C | o Contains problem specific forcing for merid velocity.
274     C *==========================================================*
275     C | Adds terms to gS for forcing by external sources
276     C | e.g. fresh-water flux, climatalogical relaxation, etc ...
277     C *==========================================================*
278     C \ev
279    
280     C !USES:
281 adcroft 1.2 IMPLICIT NONE
282     C == Global data ==
283     #include "SIZE.h"
284     #include "EEPARAMS.h"
285     #include "PARAMS.h"
286     #include "GRID.h"
287     #include "DYNVARS.h"
288     #include "FFIELDS.h"
289    
290 jmc 1.6 C !INPUT/OUTPUT PARAMETERS:
291 adcroft 1.2 C == Routine arguments ==
292 jmc 1.6 C iMin,iMax :: Working range of x-index for applying forcing.
293     C jMin,jMax :: Working range of y-index for applying forcing.
294     C bi,bj :: Current tile indices
295     C kLev :: Current vertical level index
296     C myTime :: Current time in simulation
297     C myThid :: Thread Id number
298 adcroft 1.2 INTEGER iMin, iMax, jMin, jMax, kLev, bi, bj
299 jmc 1.6 _RL myTime
300 adcroft 1.2 INTEGER myThid
301    
302 jmc 1.6 C !LOCAL VARIABLES:
303 adcroft 1.2 C == Local variables ==
304 jmc 1.6 C i,j :: Loop counters
305     c INTEGER i, j
306     CEOP
307 adcroft 1.2
308     C-- Forcing term(s)
309    
310     RETURN
311     END

  ViewVC Help
Powered by ViewVC 1.1.22