1 |
#include "ctrparam.h" |
2 |
|
3 |
! ========================================================== |
4 |
! |
5 |
! EDDYPA_NEW.F: THIS ROUTINE DEALS WITH THE PARAMETERIZATION OF |
6 |
! EDDY TRANSPORTS AND VARIANCES |
7 |
! Modified version, includes real averaging |
8 |
! |
9 |
! |
10 |
SUBROUTINE EDDYPA (UT,VT,TT,PT,QT,U,V,T,P,Q,DTE,EDDFLX) 2181. |
11 |
C 2181.5 |
12 |
C THIS ROUTINE DEALS WITH THE PARAMETERIZATION OF 2182. |
13 |
C EDDY TRANSPORTS AND VARIANCES 2182.5 |
14 |
C 2183. |
15 |
#if ( defined CPL_CHEM ) |
16 |
! |
17 |
#include "chem_para" |
18 |
#include "chem_com" |
19 |
! |
20 |
#endif |
21 |
|
22 |
#include "BD2G04.COM" |
23 |
|
24 |
COMMON/WORK3/PHI(IM0,JM0,LM0) 2184. |
25 |
COMMON/WORK4/FD(IM0,JM0) 2184.5 |
26 |
DIMENSION UT(IM0,JM0,LM0),VT(IM0,JM0,LM0),TT(IM0,JM0,LM0), |
27 |
& PT(IM0,JM0), |
28 |
* QT(IM0,JM0,LM0) 2185.5 |
29 |
DIMENSION SINL(JM0),COSL(JM0), 2186. |
30 |
* DXU(JM0),DYU(JM0),DXYU(JM0),PHIS(IO0,JM0),UX(IO0,JM0,2*LM0) 2186.5 |
31 |
EQUIVALENCE (SINL(1),SINP(1)),(COSL(1),COSP(1)), 2187. |
32 |
* (DXU(1),DXV(1)),(DYU(1),DYV(1)),(DXYU(1),DXYV(1)), 2187.5 |
33 |
* (PHIS,FDATA) 2188. |
34 |
c REAL*4 KAPAP1 2188.5 |
35 |
REAL KAPAP1,MOMCOEF |
36 |
COMMON/EPARA/VTH(JM0,LM0),WTH(JM0,LM0),VU(JM0,LM0),VV(JM0,LM0) |
37 |
& ,DQSDT(JM0,LM0) 2189. |
38 |
* ,DWV(JM0),PHIT(JM0,LM0),TPRIM2(JM0,LM0),WU(JM0,LM0),CKS,CKN 2189.5 |
39 |
* ,WQ(JM0,LM0),VQ(JM0,LM0) 2190. |
40 |
DIMENSION PHITAV(LM0),TPRIM0(JM0) |
41 |
COMMON/SPEC2/KM,KINC,COEK,C3LAND(IO0,JM0),C3OICE(IO0,JM0) 2190.5 |
42 |
* ,C3LICE(IO0,JM0),WMGE(IO0,JM0),TSSFC(1,JM0,4) 2191. |
43 |
DIMENSION TSZM(JM0),THSZM(JM0),FO(JM0),WUDATA(JM0,LM0 ) |
44 |
& ,VUDATA(JM0,LM0), 2191.5 |
45 |
& VSH(JM0,LM0),WSH(JM0,LM0),TEMZ(JM0,LM0),COEKY(JM0) |
46 |
& ,DTHDY(JM0,LM0),TEMJL(JM0,LM0+1) 2192. |
47 |
* ,DSHDY(JM0,LM0),DSHDP(JM0,LM0+1),DTHDZ(JM0,LM0+1) |
48 |
& ,DSHDZ(JM0,LM0+1) 2192.5 |
49 |
* ,SHSAT(JM0,LM0),DSHDYS(JM0,LM0) 2193. |
50 |
DIMENSION DQDYJ(JM0),DQDY(JM0,LM0),UZM(JM0),PSROK(JM0,LM0) |
51 |
& ,DQSDT1(JM0,LM0), 2193.5 |
52 |
* DUDP(JM0,LM0),TZA(JM0,LM0),DTHDP(JM0,LM0+1),PRESR(JM0,LM0+1) |
53 |
& ,PSRK(JM0,LM0+1) 2194. |
54 |
* ,PSI(JM0,LM0),PVOR(JM0,LM0),VQZ(JM0),UNEW(JM0,LM0),UOLD(JM0,LM0) |
55 |
& ,FMIX(JM0,LM0) 2194.5 |
56 |
* ,SCLH(JM0,LM0),GAMA(JM0),SHJ(JM0),HTJ(JM0),ZSTAR(JM0,LM0) |
57 |
& ,BVF(JM0) 2195. |
58 |
* ,PRSRU(JM0,LM0),DTDYJ(JM0),DTDZJ(JM0),DWVJL(JM0,LM0) |
59 |
& ,BVFJL(JM0,LM0), 2195.5 |
60 |
* DUDZ(JM0,LM0),DUDZJ(JM0),COEJL(JM0,LM0),DQDY1(JM0,LM0) |
61 |
& ,COEJL1(JM0,LM0) 2196. |
62 |
* ,UVAR(JM0),TEMJ(JM0),WTDATA(JM0,LM0),WQDATA(JM0,LM0) |
63 |
& ,VTDATA(JM0,LM0), 2196.5 |
64 |
* VQDATA(JM0,LM0),VVDATA(JM0,LM0),TTDATA(JM0,LM0),GAMAJ(JM0) |
65 |
& ,DTHDZJ(JM0) 2197. |
66 |
* ,BVFM(JM0,LM0),BVFMJ(JM0),DQDZ(JM0,LM0),VLAT(JM0),VVN(9) |
67 |
& ,VVS(9) |
68 |
& ,VVNEWA(LM0),FOV(JM0),BETAV(JM0),BETAP(JM0),TEMKGP(JM0) |
69 |
& ,TEMKGV(JM0),RDWVERR(JM0),ADWVERR(JM0) |
70 |
& ,TEMKV(JM0),TEMKP(JM0),GAMAP(JM0),DWVP(JM0),HTJP(JM0) |
71 |
& ,DWVA(JM0),DWVPA(JM0) |
72 |
dimension TZN(JM0,LM0),QZN(JM0,LM0),UZN(JM0,LM0),VZN(JM0,LM0) |
73 |
common/tprmtg/tprmg(JM0),ntprmg(JM0) |
74 |
common/Dscale/DWAV0(JM0) |
75 |
logical first,uvtr,wvtr,uttr,wttr,uqtr,wqtr,EDDFLX |
76 |
&,OLDUVVERT,HPRNT |
77 |
SAVE |
78 |
data first/.true./ |
79 |
DATA VVN/37.6,39.5,57.8,94.8,150.6,236.3,329.5,114.8,43.8/ 2261.5 |
80 |
DATA VVS/32.7,28.3,33.7,44.9,79.4,151.5,193.6,90.9,15.7/ 2262. |
81 |
C 2262.5 |
82 |
QSAT(TM,PR,QL)=3.797915*EXP(QL*(7.93252E-6-2.166847E-3/TM))/PR 2263. |
83 |
C 2263.5 |
84 |
c HPRNT=TAU.gt.50.0.and.TAU.lt.75.0 |
85 |
HPRNT=.FALSE. |
86 |
c HPRNT=TAU.ge.17520.0.and.TAU.lt.17521.0 |
87 |
c HPRNT=.TRUE. |
88 |
JPR=24 |
89 |
uvtr=EDDFLX |
90 |
wvtr=EDDFLX |
91 |
uttr=EDDFLX |
92 |
wttr=EDDFLX |
93 |
uqtr=EDDFLX |
94 |
wqtr=EDDFLX |
95 |
if(first)then |
96 |
LZ=LM 2264. |
97 |
LZM1=LZ-1 2264.5 |
98 |
LZM2=LZ-2 2265. |
99 |
JHALF=JM/2 2265.5 |
100 |
JHAFM1=JHALF-1 2266. |
101 |
JHAFP2=JHALF+2 2266.5 |
102 |
JMM2=JM-2 2267. |
103 |
IX=IM 2267.5 |
104 |
IS=IM 2268. |
105 |
FIS=IS 2268.5 |
106 |
CKS=1. 2268.6 |
107 |
CKN=1. 2268.7 |
108 |
CP=RGAS/KAPA 2269.5 |
109 |
SHA=RGAS/KAPA 2270. |
110 |
CONSTA=2.0 |
111 |
CONSTA=1.5 |
112 |
CONSTA=1.0 |
113 |
RCZZ=1.35 |
114 |
c RCZZ=1.83 |
115 |
MOMCOEF=1. |
116 |
c MOMCOEF=4.*0.5*(0.57+0.71)/CONSTA |
117 |
KAPAP1=KAPA+1. 2271. |
118 |
HLAT=LHE 2271.5 |
119 |
CLH=HLAT/CP 2272. |
120 |
DQDTX=.622*HLAT/RGAS 2272.5 |
121 |
P1000=1000.**KAPA 2273. |
122 |
DTDYMN=1.E-6/P1000 2273.5 |
123 |
DQDYMN=1.E-09 2274. |
124 |
DSCALE=626.*1.e3 |
125 |
c DSCALE=323.*1.e3 |
126 |
c DSCALE=1000.*1.e3 |
127 |
FDR=exp(-DLAT*RADIUS/DSCALE) |
128 |
OLDUVVERT=.false. |
129 |
DWVCRI=1000. |
130 |
DWVCRI=100. |
131 |
JINT=6 |
132 |
IF(JM.EQ.46)JINT=12 |
133 |
c JINT=JM/2-2 |
134 |
JEQ=JM/2+1 |
135 |
JEQP1=JEQ+1 |
136 |
JEQM1=JEQ-1 |
137 |
JEQM2=JEQ-2 |
138 |
LBGN=1 |
139 |
ALPHA=1. |
140 |
NVAIT=2 |
141 |
NVAIT=3 |
142 |
WMGMIN=3.0 |
143 |
WMGMIN=0.0 |
144 |
print *,' Mom. flux param. for dry atmosph.' |
145 |
print *,' New vertical avereging' |
146 |
print *,' Change in YCUT calculations' |
147 |
print *,' with tg difference' |
148 |
print *,' CONSTA=',CONSTA |
149 |
print *,' RCZZ=',RCZZ |
150 |
print *,' DWVCRI=',DWVCRI |
151 |
print *,' MRCH=',MRCH |
152 |
print *,' FDR=',FDR,' L=',DSCALE |
153 |
print *,' JINT=',JINT |
154 |
print *,' NVAIT=',NVAIT |
155 |
print *,' PSF=',PSF |
156 |
print *,' WMGMIN=',WMGMIN |
157 |
print *,' With A for momentum=',CONSTA*MOMCOEF |
158 |
print *,' UVTR=',uvtr,' WVTR=',wvtr |
159 |
print *,' UTTR=',uttr,' WTTR=',wttr |
160 |
print *,' UQTR=',uqtr,' WQTR=',wqtr |
161 |
if(OLDUVVERT)then |
162 |
print *,' Vertical structure of UV flux from GISS |
163 |
& 2D model' |
164 |
endif |
165 |
DO 3 J=1,JM 2274.5 |
166 |
3 FO(J)=F(J)/DXYP(J) 2275. |
167 |
DO 4 J=2,JM 2275.5 |
168 |
VLAT(J)=(LAT(J)+LAT(J-1))*.5 2276. |
169 |
FOV(J)=0.5*(FO(J)+FO(J-1)) |
170 |
BETAV(J)=(FO(J)-FO(J-1))/DYU(J) |
171 |
if(j.NE.JM)BETAP(J)=.5*(FO(J+1)-FO(J-1))/DYU(J) |
172 |
TEMY=RADIUS*(.25*TWOPI-ABS(VLAT(J))) |
173 |
Xk=1./(RADIUS*COSV(j)) |
174 |
Yk=TWOPI/(4.*TEMY) |
175 |
TEMKGV(j)=SQRT(Xk**2+Yk**2) |
176 |
4 CONTINUE |
177 |
TPR00=1.0 |
178 |
TPR90=2.5 |
179 |
TPR11=0.5*(TPR90+TPR00) |
180 |
TPR12=0.5*(TPR90-TPR00) |
181 |
print *,'TPR00=',TPR00,' TPR90=',TPR90 |
182 |
do j=1,jm |
183 |
TEMY=RADIUS*(.25*TWOPI-ABS(LAT(J))) |
184 |
Xk=1./(RADIUS*COSP(j)) |
185 |
Yk=TWOPI/(4.*TEMY) |
186 |
TEMKGP(j)=SQRT(Xk**2+Yk**2) |
187 |
tprimrad = 3.14159*(float(j-1)/22.5) |
188 |
TPRIM0(j)=5.5+(4.5*cos(tprimrad)) ! 2323 |
189 |
TPRIM0(j)=0.0 ! 2324 |
190 |
! TPRIM0(j)=1.5+(1.0*cos(tprimrad)) ! 2327 2332 |
191 |
TPRIM0(j)=TPR11+(TPR12*cos(tprimrad)) ! 2327 2332 |
192 |
! TPRIM0(j)=2.0 ! 2331 |
193 |
enddo |
194 |
c print *,'LAT' |
195 |
c print *,LAT |
196 |
c print *,'VLAT' |
197 |
c print *,VLAT |
198 |
c print *,'COSV' |
199 |
c print *,COSV |
200 |
print *,'TPRIM0' |
201 |
print *,TPRIM0 |
202 |
first=.false. |
203 |
endif |
204 |
C 2276.5 |
205 |
#if ( defined CPL_CHEM ) |
206 |
! |
207 |
do 1 j=1,nlat |
208 |
beta1(j)=0.0 |
209 |
beta2(j)=0.0 |
210 |
do 2 k=1,nlev |
211 |
beta3(j,k)=0.0 |
212 |
beta4(j,k)=0.0 |
213 |
fkt (j,k)=0.0 |
214 |
2 continue |
215 |
1 continue |
216 |
! |
217 |
#endif |
218 |
|
219 |
C reassign T,Q,U and V to internal variables |
220 |
|
221 |
do l=1,lm |
222 |
do j=1,jm |
223 |
TZN(j,l)=T(1,j,l) |
224 |
QZN(j,l)=Q(1,j,l) |
225 |
UZN(j,l)=U(1,j,l) |
226 |
VZN(j,l)=V(1,j,l) |
227 |
enddo |
228 |
enddo |
229 |
if(HPRNT)then |
230 |
print *,' Eddypa UZN before filter' |
231 |
print *,'LZM2=',LZM2 |
232 |
print *,(UZN(13,L),L=1,LM) |
233 |
print *,(UZN(14,L),L=1,LM) |
234 |
endif |
235 |
if(JM.eq.46)then |
236 |
CALL FLTR(TZN,JM,LM,1) |
237 |
CALL FLTR(QZN,JM,LM,1) |
238 |
CALL FLTR(UZN,JM,LM,2) |
239 |
endif |
240 |
|
241 |
DO 6 J=1,JM 2277. |
242 |
DWV(J)=0. 2277.5 |
243 |
DO 6 I=1,IM 2278. |
244 |
6 WMGE(I,J)=WMGMIN 2278.5 |
245 |
DO 10 L=1,LM 2279. |
246 |
DO 10 J=1,JM 2279.5 |
247 |
VU(J,L)=0. 2280. |
248 |
VQ(J,L)=0. 2280.5 |
249 |
VTH(J,L)=0. 2281. |
250 |
WTH(J,L)=0. 2281.5 |
251 |
WU(J,L)=0. 2282. |
252 |
WQ(J,L)=0. 2282.5 |
253 |
VV(J,L)=0. 2283. |
254 |
! |
255 |
!CAS Applying latitudinal depedence of convection params |
256 |
! |
257 |
! tprimrad = 3.14159*(float(j-1)/22.5) |
258 |
! TPRIMCO=6.0+(4.0*cos(tprimrad)) |
259 |
TPRIM2(J,L)=TPRIM0(j)/(EXPBYK(SIG(L)*P(1,J)+PTOP))**2. |
260 |
!CAS TPRIM2(J,L)=2./(EXPBYK(SIG(L)*P(1,J)+PTOP))**2. 2283.5 |
261 |
! AJL(J,L,59)=AJL(J,L,59)+2.*PSF |
262 |
AJL(J,L,59)=AJL(J,L,59)+TPRIM0(j)*PSF |
263 |
c TPRIM2(J,L)=0 |
264 |
PHIT(J,L)=PHI(1,J,L) 2284. |
265 |
10 CONTINUE 2284.5 |
266 |
DO L=1,LM |
267 |
PHITAV(L)=0. |
268 |
CCC=0. |
269 |
DO J=1,JM |
270 |
PHITAV(L)=PHITAV(L)+PHIT(J,L)*COSP(J) |
271 |
CCC=CCC+COSP(J) |
272 |
ENDDO |
273 |
PHITAV(L)=PHITAV(L)/CCC |
274 |
ENDDO |
275 |
C 2285. |
276 |
C COMPUTE PRESSURE, P**KAPA, STATIC STABILITY AND SCALE HEIGHT 2285.5 |
277 |
C 2286. |
278 |
DO 24 L=1,LZ 2286.5 |
279 |
DO 24 J=1,JM 2287. |
280 |
PRESR(J,L+1)= P(1,J)*SIGE(L+1)+PTOP 2287.5 |
281 |
PSRK(J,L+1)=PRESR(J,L+1)**KAPA 2288. |
282 |
IF(J.GE.2) PRSRU(J,L)=.5*(P(1,J)+P(1,J-1))*SIG(L)+PTOP 2288.5 |
283 |
IF(J.GE.2) PSROK(J,L)=PRSRU(J,L)**KAPA 2289. |
284 |
TEM=( P(1,J)*SIG(L)+PTOP)**KAPA 2289.5 |
285 |
c24 TZA(J,L)= T(1,J,L)*TEM 2290. |
286 |
24 TZA(J,L)= TZN(J,L)*TEM |
287 |
DO 23 J=1,JM 2290.5 |
288 |
DO 22 L=1,LZM1 2291. |
289 |
c TEM= T(1,J,L+1)- T(1,J,L) 2291.5 |
290 |
TEM= TZN(J,L+1)- TZN(J,L) |
291 |
IF(TEM.LT..1) TEM=.1 2292. |
292 |
DTHDP(J,L)=-TEM/(P(1,J)*DSIGO(L)) 2292.5 |
293 |
c DSHDP(J,L)=(Q(1,J,L)-Q(1,J,L+1))/(P(1,J)*DSIGO(L)) 2293. |
294 |
DSHDP(J,L)=(QZN(J,L)-QZN(J,L+1))/(P(1,J)*DSIGO(L)) |
295 |
TEMJL(J,L+1)=-TEM/(SIG(L+1)-SIG(L))*RGAS* 2293.5 |
296 |
* PSRK(J,L+1)/PRESR(J,L+1)/ P(1,J) 2294. |
297 |
TEMX=GRAV*PRESR(J,L+1)/(.5*RGAS*(TZA(J,L)+TZA(J,L+1))) 2294.5 |
298 |
|
299 |
#if ( defined CPL_CHEM ) |
300 |
! |
301 |
deltap(j,l)=1./(p(1,j)*dsigo(l)) |
302 |
dp2dz (j,l)=-temx |
303 |
! |
304 |
#endif |
305 |
|
306 |
|
307 |
DTHDZ(J,L+1)=-DTHDP(J,L)*TEMX 2295. |
308 |
DSHDZ(J,L+1)=-DSHDP(J,L)*TEMX 2295.5 |
309 |
22 CONTINUE 2296. |
310 |
TEMJL(J,1)=(TEMJL(J,2)*(DSIG(1)+DSIG(2))-TEMJL(J,3)*DSIG(1))/ 2296.5 |
311 |
* DSIG(2) 2297. |
312 |
IF(TEMJL(J,1).LE..1*TEMJL(J,2)) TEMJL(J,1)=.1*TEMJL(J,2) 2297.5 |
313 |
TEMJL(J,LZ+1)=TEMJL(J,LZ) 2298. |
314 |
PRESR(J,1)= P(1,J)+PTOP 2298.5 |
315 |
PSRK(J,1)=PRESR(J,1)**KAPA 2299. |
316 |
DTHDZ(J,1)=DTHDZ(J,2) 2299.5 |
317 |
DSHDZ(J,1)=DSHDZ(J,2) 2300. |
318 |
DTHDZ(J,LZ+1)=DTHDZ(J,LZ) 2300.5 |
319 |
DSHDZ(J,LZ+1)=DSHDZ(J,LZ) 2301. |
320 |
23 CONTINUE 2301.5 |
321 |
DO 27 L=1,LZ 2302. |
322 |
DO 26 J=2,JMM1 2302.5 |
323 |
26 TEMZ(J,L)=0.25*(TEMJL(J,L)*2.+TEMJL(J+1,L)+TEMJL(J-1,L)) 2303. |
324 |
TEMJL(1,L)=(2.*TEMJL(1,L)+TEMJL(2,L))/3. 2303.5 |
325 |
TEMJL(JM,L)=(2.*TEMJL(JM,L)+TEMJL(JMM1,L))/3. 2304. |
326 |
DO 27 J=2,JMM1 2304.5 |
327 |
TEMJL(J,L)=TEMZ(J,L) 2305. |
328 |
27 CONTINUE 2305.5 |
329 |
DO 21 L=1,LZM1 2306. |
330 |
DO 21 J=1,JM 2306.5 |
331 |
TEM=2.*PRESR(J,L+1)*GRAV/RGAS/(TZA(J,L)+TZA(J,L+1)) 2307. |
332 |
BVFJL(J,L)=TEM*SQRT(TEMJL(J,L+1)) 2307.5 |
333 |
C DRY |
334 |
c BVFM(J,L)=BVFJL(J,L)*BVFJL(J,L)+2.*DSHDZ(J,L+1)*CLH*GRAV/ 2308. |
335 |
c * (TZA(J,L)+TZA(J,L+1)) 2308.5 |
336 |
BVFM(J,L)=BVFJL(J,L)*BVFJL(J,L) |
337 |
c RVFM is now (Nd)**2 not (Nm)**2 |
338 |
C DRY |
339 |
SCLH(J,L)=RGAS*.5*(TZA(J,L)+TZA(J,L+1))/GRAV 2309. |
340 |
21 CONTINUE 2309.5 |
341 |
C 2310. |
342 |
2310.5 |
343 |
C COMPUTE DEL QS/DEL T 2311. |
344 |
C 2311.5 |
345 |
DO 30 L=1,LZM1 2312. |
346 |
DO 30 J=1,JM 2312.5 |
347 |
TEMZ(J,L)=0. 2313. |
348 |
DQSDT1(J,L)=0. 2313.5 |
349 |
DO 25 I=1,IS 2314. |
350 |
PR=(P(I,J)*SIG(L)+PTOP) 2314.5 |
351 |
TM=TZA(J,L) 2315. |
352 |
TEM=QSAT(TM,PR,HLAT) 2315.5 |
353 |
SHSAT(J,L)=TEM 2316. |
354 |
c TEM=DQDTX/TM/TM*Q(1,J,L) 2316.5 |
355 |
c TEM1=DQDTX*Q(1,J,L)/TM/TM 2317. |
356 |
TEM=DQDTX/TM/TM*QZN(J,L) |
357 |
TEM1=DQDTX*QZN(J,L)/TM/TM |
358 |
TEMZ(J,L)=TEMZ(J,L)+TEM 2317.5 |
359 |
DQSDT1(J,L)=DQSDT1(J,L)+TEM1 2318. |
360 |
25 CONTINUE 2318.5 |
361 |
DQSDT1(J,L)=DQSDT1(J,L)*CLH/FIS 2319. |
362 |
30 DQSDT(J,L)=TEMZ(J,L)/FIS 2319.5 |
363 |
C 2320. |
364 |
C COMPUTE MERIDIONAL TEMPERATURE AND MOISTURE GRADIENT 2320.5 |
365 |
C 2321. |
366 |
DO 3140 L=1,LZ 2321.5 |
367 |
DO 3140 J=2,JM 2322. |
368 |
c DTHDY(J,L)=( T(1,J,L)- T(1,J-1,L))/DYU(J) 2322.5 |
369 |
DTHDY(J,L)=( TZN(J,L)- TZN(J-1,L))/DYU(J) |
370 |
3140 CONTINUE 2323. |
371 |
DO 3141 L=1,LZ 2323.5 |
372 |
DO 3142 J=3,JMM1 2324. |
373 |
3142 TEMZ(J,L)=.25*(2.*DTHDY(J,L)+DTHDY(J+1,L)+DTHDY(J-1,L)) 2324.5 |
374 |
DTHDY(2,L)=(2.*DTHDY(2,L)+DTHDY(3,L))/3. 2325. |
375 |
DTHDY(JM,L)=(2.*DTHDY(JM,L)+DTHDY(JMM1,L))/3. 2325.5 |
376 |
DO 3141 J=3,JMM1 2326. |
377 |
3141 DTHDY(J,L)=TEMZ(J,L) 2326.5 |
378 |
DO 3144 L=1,LZ 2327. |
379 |
DO 3144 J=2,JM 2327.5 |
380 |
c DSHDY(J,L)=( Q(1,J,L)- Q(1,J-1,L))/DYU(J) 2328. |
381 |
DSHDY(J,L)=( QZN(J,L)- QZN(J-1,L))/DYU(J) |
382 |
3144 CONTINUE 2328.5 |
383 |
DO 3148 L=1,LZ 2329. |
384 |
DO 3148 J=2,JM 2329.5 |
385 |
DSHDYS(J,L)=(SHSAT(J,L)-SHSAT(J-1,L))/DYU(J) 2330. |
386 |
3148 CONTINUE 2330.5 |
387 |
C DO 3147 L=1,LZ 2331. |
388 |
C DO 3146 J=3,JMM1 2331.5 |
389 |
C3146 TEMZ(J,L)=.25*(2.*DSHDY(J,L)+DSHDY(J+1,L)+DSHDY(J-1,L)) 2332. |
390 |
C DSHDY(2,L)=(2.*DSHDY(2,L)+DSHDY(3,L))/3. 2332.5 |
391 |
C DSHDY(JM,L)=(2.*DSHDY(JM,L)+DSHDY(JMM1,L))/3. 2333. |
392 |
C DO 3147 J=3,JMM1 2333.5 |
393 |
C3147 DSHDY(J,L)=TEMZ(J,L) 2334. |
394 |
C 2334.5 |
395 |
C COMPUTE WAVE DEPTH 2335. |
396 |
C 2335.5 |
397 |
DO 3145 J=2,JM 2336.5 |
398 |
B45=BETAV(J) |
399 |
FCONT=FOV(J) |
400 |
DO 3145 L=1,LZM2 2338. |
401 |
TEM1=.5*(BVFJL(J,L)+BVFJL(J-1,L)) 2338.5 |
402 |
TEM=.5*(SCLH(J,L)+SCLH(J-1,L))*TEM1*TEM1*B45*.25* 2339. |
403 |
c * (T(1,J,L)+T(1,J,L+1)+T(1,J-1,L+1)+ 2339.5 |
404 |
* (TZN(J,L)+TZN(J,L+1)+TZN(J-1,L+1)+ |
405 |
c * T(1,J-1,L))/(FCONT*GRAV*.5*(DTHDY(J,L)+DTHDY(J,L+1))+1.E-20) 2340. |
406 |
* TZN(J-1,L))/(FCONT*GRAV*.5*(DTHDY(J,L)+DTHDY(J,L+1))+1.E-20) |
407 |
c DWVJL(J,L)= .5*(SCLH(J,L)+SCLH(J-1,L))/ 2340.5 |
408 |
c * (1.48*ABS(TEM)+.48) 2341. |
409 |
DWVJL(J,L)= .5*(SCLH(J,L)+SCLH(J-1,L))/ |
410 |
& (2./RCZZ*(1.+ABS(TEM))-1.) |
411 |
IF(DWVJL(J,L).LT.DWVCRI) DWVJL(J,L)=DWVCRI 2341.5 |
412 |
3145 CONTINUE 2342. |
413 |
C 2342.5 |
414 |
C DETERMINE THE LAT.FOR CUTTING OFF THE WAVE PROPAGATION 2343. |
415 |
C 2343.5 |
416 |
DO 3150 J=2,JM 2345. |
417 |
USUM=0. 2345.5 |
418 |
DO 3151 L=1,LM 2346. |
419 |
c3151 USUM=USUM+U(1,J,L)*DSIG(L) 2346.5 |
420 |
3151 USUM=USUM+UZN(J,L)*DSIG(L) |
421 |
3150 UZM(J)=USUM 2347. |
422 |
JB=JEQ+1 |
423 |
JE=JB+JINT 2348. |
424 |
YCUTN=0. 2348.5 |
425 |
DO 3155 J=JB,JE 2349. |
426 |
JR=JB+JE-J 2349.5 |
427 |
3155 IF(UZM(JR).LE.0.) GO TO 3156 2350. |
428 |
GO TO 3157 2350.5 |
429 |
3156 YCUTN=(UZM(JR)*VLAT(JR+1)-UZM(JR+1)*VLAT(JR))/ |
430 |
* (UZM(JR)-UZM(JR+1)-1.E-20) |
431 |
3157 CONTINUE |
432 |
JE=JEQ-1 |
433 |
JB=JE-JINT 2352.5 |
434 |
YCUTS=0. 2353. |
435 |
DO 3158 JR=JB,JE 2353.5 |
436 |
3158 IF(UZM(JR).LE.0.) GO TO 3159 2354. |
437 |
GO TO 3160 2354.5 |
438 |
3159 YCUTS=(UZM(JR)*VLAT(JR-1)-UZM(JR-1)*VLAT(JR))/ |
439 |
* (UZM(JR)-UZM(JR-1)-1.E-20) |
440 |
3160 CONTINUE 2356. |
441 |
c print *,'YCUTN=',180./TWOPI*YCUTN,' YCUTS=',180./TWOPI*YCUTS |
442 |
C 2356.5 |
443 |
C 2356.5 |
444 |
C NEW WEIGHTING VARIABLES BY WAVE DEPTH |
445 |
C |
446 |
c print *,' DWAV0' |
447 |
c print *,DWAV0 |
448 |
DO J=1,JM |
449 |
DWV(J)=DWAV0(J) |
450 |
ENDDO |
451 |
DO II=1,NVAIT |
452 |
DO J=1,JM |
453 |
DO L=1,LZM1 |
454 |
TEMZ(J,L)=(TZA(J,L+1)-TZA(J,L))*GRAV/(PHI(1,J,L+1)-PHI(1,J,L)) |
455 |
ENDDO ! L |
456 |
ENDDO ! J |
457 |
CALL VWEIGHAV(GRAV,TEMZ,DTDZJ,PHI,PHIS,DWV,1,LZM2,1,JM ,PTOP, |
458 |
* SCLH,P,SIGE,DSIGO,CKS,CKN,IO,IM,JM,LM) |
459 |
DO J=2,JM |
460 |
DO L=1,LZM2 |
461 |
TEMZ(J,L)=.5*(DTHDY(J,L)+DTHDY(J,L+1)) |
462 |
ENDDO ! L |
463 |
ENDDO ! J |
464 |
CALL VWEI1AV(GRAV,TEMZ ,DTDYJ,PHI,PHIS,DWV,1,LZM2,2,JM,PTOP, |
465 |
* SCLH,P,SIGE,DSIGO,CKS,CKN,IO,IM,JM,LM) |
466 |
DO L=1,LZM1 |
467 |
DO J=1,JM |
468 |
TEMZ(J,L)=DTHDZ(J,L+1) |
469 |
ENDDO ! L |
470 |
ENDDO ! J |
471 |
CALL VWEIGHAV(GRAV,TEMZ,DTHDZJ,PHI,PHIS,DWV,1,LZM2,1,JM ,PTOP, |
472 |
* SCLH,P,SIGE,DSIGO,CKS,CKN,IO,IM,JM,LM) |
473 |
c Vertical averaging of N**2 |
474 |
CALL VWEIGHAV(GRAV,BVFM,BVFMJ,PHI,PHIS,DWV,1,LZM2,1,JM ,PTOP, |
475 |
* SCLH,P,SIGE,DSIGO,CKS,CKN,IO,IM,JM,LM) |
476 |
c Nw=sqrt((N**2)w) |
477 |
DO j=1,JM |
478 |
BVF(j)=sqrt(BVFMJ(j)) |
479 |
ENDDO |
480 |
c Nw=sqrt((N**2)w) |
481 |
CALL VWEIGHAV(GRAV,SCLH,SHJ, PHI,PHIS,DWV,1,LZM2,1,JM ,PTOP, |
482 |
* SCLH,P,SIGE,DSIGO,CKS,CKN,IO,IM,JM,LM) |
483 |
DO L=1,LZM1 |
484 |
DO J=1,JM |
485 |
c TEMZ(J,L)=.5*(T(1,J,L)+T(1,J,L+1)) |
486 |
TEMZ(J,L)=.5*(TZN(J,L)+TZN(J,L+1)) |
487 |
ENDDO ! L |
488 |
ENDDO ! J |
489 |
CALL VWEIGHAV(GRAV,TEMZ,THSZM, PHI,PHIS,DWV,1,LZM2,1,JM ,PTOP, |
490 |
* SCLH,P,SIGE,DSIGO,CKS,CKN,IO,IM,JM,LM) |
491 |
DO L=1,LZM2 |
492 |
DO J=2,JM |
493 |
c TEMZ(J,L)=.5*(U(1,J,L)+U(1,J,L+1)) |
494 |
TEMZ(J,L)=.5*(UZN(J,L)+UZN(J,L+1)) |
495 |
ENDDO ! L |
496 |
ENDDO ! J |
497 |
if(HPRNT)then |
498 |
J=JM-1 |
499 |
print *,'LZ=',LZ |
500 |
print *,' Eddypa UZN before VWEI1AV' |
501 |
print *,(UZN(13,L),L=1,LM) |
502 |
print *,(UZN(14,L),L=1,LM) |
503 |
print *,DWV(J+1),DWV(J-1),DWV(J) |
504 |
endif |
505 |
CALL VWEI1AV(GRAV,TEMZ,UZM,PHI,PHIS,DWV,1,LZM2,2,JM, |
506 |
* PTOP,SCLH,P,SIGE,DSIGO,CKS,CKN,IO,IM,JM,LM) |
507 |
if(HPRNT)then |
508 |
print *,' Eddypa UZM after VWEI1AV' |
509 |
J=14 |
510 |
print *,UZM(J+1),UZM(J-1),UZM(J) |
511 |
endif |
512 |
GAMAD=GRAV/CP |
513 |
DO J=1,JM |
514 |
DO L=1,LZM1 |
515 |
TEDGE=.5*(TZA(J,L+1)+TZA(J,L)) |
516 |
QSE=QSAT(TEDGE,PRESR(J,L+1),LHE) |
517 |
TEMZ(J,L)=GAMAD*(1.+LHE*QSE/(RGAS*TEDGE))/ |
518 |
* (1.+.622*LHE*LHE*QSE/(CP*RGAS*TEDGE*TEDGE)) |
519 |
ENDDO ! L |
520 |
ENDDO ! J |
521 |
CALL VWEIGHAV(GRAV,TEMZ,GAMAJ,PHI,PHIS,DWV,1,LZM2,1,JM ,PTOP, |
522 |
* SCLH,P,SIGE,DSIGO,CKS,CKN,IO,IM,JM,LM) |
523 |
DO J=2,JM |
524 |
B45=BETAV(J) |
525 |
FCONT=FOV(J) |
526 |
c |
527 |
c BVFMJ=(N**2)w |
528 |
TEM1=.5*(BVFMJ(J)+BVFMJ(J-1)) |
529 |
TEM=.5*(SHJ(J)+SHJ(J-1))*TEM1*B45* .5*(THSZM(J)+THSZM(J-1))/ |
530 |
* (FCONT*GRAV*DTDYJ(J)+1.E-20) |
531 |
c |
532 |
GAMA(J)=ABS(TEM) |
533 |
HTJ(J)=.5*(SHJ(J)+SHJ(J-1))/(1.+GAMA(J)) |
534 |
IF(HTJ(J).LT.DWVCRI) HTJ(J)=DWVCRI |
535 |
c DWV(J)=.5*(SHJ(J)+SHJ(J-1))/(1.48*GAMA(J)+.48) |
536 |
DWV(J)=.5*(SHJ(J)+SHJ(J-1))/(2./RCZZ*(1.+GAMA(J))-1.) |
537 |
IF(DWV(J).LT.DWVCRI) DWV(J)=DWVCRI |
538 |
IF(BVFMJ(J).LT.1.E-5) BVFMJ(J)=1.E-5 |
539 |
C |
540 |
DWV(J)=0.5*(DWV(j)+DWAV0(J)) |
541 |
C |
542 |
ENDDO ! J |
543 |
IF(BVFMJ(1).LT.1.E-5) BVFMJ(1)=1.E-5 |
544 |
c print *,' NEW II=',ii |
545 |
c print *,' DWV' |
546 |
c print *,DWV |
547 |
AERR=0.0 |
548 |
RERR=0.0 |
549 |
AERRMAX=0.0 |
550 |
RERRMAX=0.0 |
551 |
DO J=2,JM |
552 |
ADWVERR(J)=abs(DWV(J)-DWAV0(J)) |
553 |
RDWVERR(J)=ADWVERR(J)/DWAV0(J) |
554 |
AERR=AERR+ADWVERR(J) |
555 |
RERR=RERR+RDWVERR(J) |
556 |
DWAV0(J)=DWV(J) |
557 |
if(ADWVERR(J).gt.AERRMAX)then |
558 |
AERRMAX=ADWVERR(J) |
559 |
JAMAX=j |
560 |
endif |
561 |
if(RDWVERR(J).gt.RERRMAX)then |
562 |
RERRMAX=RDWVERR(J) |
563 |
JRMAX=j |
564 |
endif |
565 |
ENDDO |
566 |
c print *,'Eddy iteration ii=',ii |
567 |
c print *,'AERR=',AERR/(jm-1.0),' RERR=',RERR/(jm-1.0) |
568 |
c print * ,'MAXAERR=',AERRMAX,' JAMAX=',JAMAX |
569 |
c print * ,'MAXRERR=',RERRMAX,' JRMAX=',JRMAX |
570 |
c print *,' DWV D' |
571 |
c print *,DWV |
572 |
c print *,' HTJ d' |
573 |
c print *,HTJ |
574 |
c print *,' SHJ H' |
575 |
c print *,SHJ |
576 |
if(RERRMAX.lt.1.e-3) go to 7788 |
577 |
ENDDO ! iterations |
578 |
7788 continue |
579 |
c if(ii.gt.20)print *,' NIT=',II |
580 |
DO J=1,JM |
581 |
DWAV0(J)=DWV(J) |
582 |
ENDDO |
583 |
c print *,' after DAAVO calculation' |
584 |
c print *,' GAMA ' |
585 |
c print *,GAMA |
586 |
c print *,' DWV' |
587 |
c print *,DWV |
588 |
c print *,' HTJ' |
589 |
c print *,HTJ |
590 |
C |
591 |
C DWV (D) and CAMA on P-grid |
592 |
do j=2,jmm1 |
593 |
TEM=SHJ(J)*BVFMJ(J)*BETAP(J)*THSZM(J)/ |
594 |
& (FO(J)*GRAV*0.5*(DTDYJ(J)+DTDYJ(J+1))) |
595 |
GAMAP(j)=abs(TEM) |
596 |
DWVP(J)=SHJ(J)/(2./RCZZ*(1.+GAMAP(J))-1.) |
597 |
IF(DWVP(J).LT.DWVCRI) DWVP(J)=DWVCRI |
598 |
HTJP(J)=SHJ(J)/(1.+GAMAP(J)) |
599 |
IF(HTJP(J).LT.DWVCRI) HTJP(J)=DWVCRI |
600 |
enddo |
601 |
c print *,' GAMAP ' |
602 |
c print *,(GAMAP(j),j=2,jm-1) |
603 |
c print *,' ' |
604 |
c print *,(0.5*(GAMA(j)+GAMA(J+1)),j=2,jm-1) |
605 |
c print *,' DWVP' |
606 |
c print *,(DWVP(j),j=2,jm-1) |
607 |
c print *,' ' |
608 |
c print *,(0.5*(DWV(j)+DWV(J+1)),j=2,jm-1) |
609 |
C NEW |
610 |
C 2387. |
611 |
c Adjust wave length near pole |
612 |
DO j=2,JM |
613 |
if(j.ne.JEQ)then |
614 |
TEM1=.5*(BVF(J)+BVF(J-1)) |
615 |
TEM=.5*(SHJ(J)+SHJ(J-1))*TEM1/ABS(FOV(J)) |
616 |
TEMH=TEMKGV(j)*TEM |
617 |
TEMKV(J)=sqrt(((1.+GAMA(J))/RCZZ)**2-0.25) |
618 |
C |
619 |
if(TEMKV(J).lt.TEMH)then |
620 |
DWVA(J)=.5*(SHJ(J)+SHJ(J-1))/(sqrt(4.*TEMH**2+1.)-1.) |
621 |
else |
622 |
DWVA(J)=DWV(J) |
623 |
endif |
624 |
endif |
625 |
ENDDO ! J |
626 |
do j=2,jmm1 |
627 |
TEM=BVF(J)*SHJ(J)/abs(FO(J)) |
628 |
TEMH=TEMKGP(j)*TEM |
629 |
TEMKP(J)=sqrt(((1.+GAMAP(J))/RCZZ)**2-0.25) |
630 |
C |
631 |
if(TEMKP(J).lt.TEMH)then |
632 |
DWVPA(J)=SHJ(J)/(sqrt(4.*TEMH**2+1.)-1.) |
633 |
else |
634 |
DWVPA(J)=DWVP(J) |
635 |
endif |
636 |
enddo |
637 |
c print *,' After adjustment' |
638 |
c print *,' DWV' |
639 |
c print *,DWV |
640 |
c print *,' DWVA' |
641 |
c print *,DWVA |
642 |
c print *,' DWVP' |
643 |
c print *,DWVP |
644 |
c print *,' DWVPA' |
645 |
c print *,DWVPA |
646 |
c print *,' HTJ' |
647 |
c print *,HTJ |
648 |
c stop |
649 |
C COMPUTE TEMPERATURE VARIANCES 2387.5 |
650 |
C 2388. |
651 |
DO 4013 J=2,JMM1 2388.5 |
652 |
COEKC=CKS 2389. |
653 |
IF(J.GT.JHALF) COEKC=CKN 2389.5 |
654 |
c TEM=BVF(J)*.25*(HTJ(J)+HTJ(J+1))*(DTDYJ(J)+DTDYJ(J+1))/ 2390. |
655 |
TEM=BVF(J)*.5*HTJP(J)*(DTDYJ(J)+DTDYJ(J+1))/ |
656 |
* (FO(J)+1.E-20) 2390.5 |
657 |
IF ((DTDYJ(J)+DTDYJ(J+1))*FO(J).GT.0.) TEM=0. 2390.6 |
658 |
c IF(FO(J).EQ.0.) TEM=0. 2391. |
659 |
if(abs(FO(J)).lt.1.e-15) TEM=0. |
660 |
PLAND=FDATA(1,J,2) 2395. |
661 |
PWATER=1.-PLAND 2395.5 |
662 |
PLICE=FDATA(1,J,3)*PLAND 2396. |
663 |
PEARTH=PLAND-PLICE 2396.5 |
664 |
POICE=ODATA(1,J,2)*PWATER 2397. |
665 |
POCEAN=PWATER-POICE 2397.5 |
666 |
if(POCEAN.LE.1.E-5)then |
667 |
POCEAN=0. |
668 |
POICE=PWATER |
669 |
endif |
670 |
tland=GDATA(1,J,4) |
671 |
toice=GDATA(1,J,3) |
672 |
tlice=GDATA(1,J,13) |
673 |
tocean=ODATA(1,J,1) |
674 |
TSS=POCEAN*tocean+POICE*toice+PLICE*tlice+PEARTH*tland |
675 |
TPR1=POCEAN*(tocean-TSS)**2 |
676 |
TPR2=POICE*(toice-TSS)**2 |
677 |
TPR3=PLICE*(tlice-TSS)**2 |
678 |
TPR4=PEARTH*(tland-TSS)**2 |
679 |
tprmg(J)=tprmg(J)+(TPR1+TPR2+TPR3+TPR4) |
680 |
ntprmg(J)=ntprmg(J)+1 |
681 |
c TEMH=.5*(DWV(J)+DWV(J+1)) 2400.5 |
682 |
TEMH=DWVPA(J) |
683 |
DO 4013 L=1,LZ 2400. |
684 |
ZPHI=(PHI(1,J,L)-PHIS(1,J))/(GRAV*TEMH+1.E-20) 2401.5 |
685 |
C TEM1=TEM*(1.-.5*ZPHI) 2402. |
686 |
TEM1=TEM*(1.-.5*ZPHI)*(1.+DQSDT1(J,L)) 2402.5 |
687 |
C for analisys |
688 |
! TEMAN=TEM*(1.-.5*ZPHI) |
689 |
C for analisys |
690 |
TPRIM2(J,L)=TPRIM2(J,L)+1.*TEM1*TEM1*EXP(-ZPHI/COEKC)/COEKC 2403. |
691 |
! 2319 |
692 |
! TPRIM2(J,L)=max(TPRIM2(J,L),1.*TEM1*TEM1*EXP(-ZPHI/COEKC)/COEKC) |
693 |
PK2=(EXPBYK(SIG(L)*P(1,J)+PTOP))**2. 2403.5 |
694 |
IF(L.EQ.1) THEN |
695 |
TPRIM2(J,L)=TPRIM2(J,L)+.6065* 2404. |
696 |
& (TPR1+TPR2+TPR3+TPR4)/PK2 2404.5 |
697 |
AJL(J,L,59)=AJL(J,L,59)+.6065* |
698 |
& (TPR1+TPR2+TPR3+TPR4)*P(1,J) |
699 |
ELSEIF(L.EQ.2) THEN |
700 |
TPRIM2(J,L)=TPRIM2(J,L)+.1653* 2405. |
701 |
& (TPR1+TPR2+TPR3+TPR4)/PK2 2405.5 |
702 |
AJL(J,L,59)=AJL(J,L,59)+.1653* |
703 |
& (TPR1+TPR2+TPR3+TPR4)*P(1,J) |
704 |
ENDIF |
705 |
AJL(J,L,59)=AJL(J,L,59)+TPRIM2(J,L) |
706 |
4013 CONTINUE 2406. |
707 |
tprmg(JM)=0. |
708 |
tprmg(1)=0. |
709 |
ntprmg(JM)=1. |
710 |
ntprmg(1)=1. |
711 |
CCC ************ |
712 |
C 2406.5 |
713 |
C COMPUTE WIND VARIANCE AND EDDY TRANSFER COEFFICIENTS 2407. |
714 |
C 2407.5 |
715 |
DO 4015 I=1,IM 2408. |
716 |
DO 4015 J=2,JM 2408.5 |
717 |
COEKD=CKS 2409. |
718 |
IF(J.GT.JHALF) COEKD=CKN 2409.5 |
719 |
TEM2=FOV(J) |
720 |
TEM3=.5*(THSZM(J)+THSZM(J-1)) 2410.5 |
721 |
TEM=HTJ(J)*DTDYJ(J)*GRAV/(TEM2*TEM3+1.E-20) 2411. |
722 |
IF (DTDYJ(J)*TEM2.GT.0.) TEM=0. 2411.1 |
723 |
c IF(TEM2.EQ.0.) TEM=0. 2411.5 |
724 |
if(abs(TEM2).lt.1.e-15) TEM=0. |
725 |
C TEML=3.*1.4142*.5*(BVF(J)+BVF(J-1))*DWV(J)/(TEM2+1.E-20) 2412. |
726 |
C TEML=ABS(TEML) 2412.5 |
727 |
C TEMY=RADIUS*(.25*TWOPI-.5*ABS(LAT(J)+LAT(J-1))) 2413. |
728 |
C RATIO=TEMY/(TEML+1.E-20) 2413.5 |
729 |
C IF(TEML.GT.TEMY) TEM=TEM*RATIO 2414. |
730 |
UVAR(J)=.5*TEM*TEM/COEKD 2414.5 |
731 |
WMGETEMP=.33*TEM*TEM/COEKD |
732 |
c WMGE(I,J)=max(WMGETEMP,WMGMIN) 2415. |
733 |
WMGE(I,J)=WMGETEMP+WMGMIN 2415.1 |
734 |
4015 continue |
735 |
IF(DMOD(TAU,25.).EQ.0..AND.IABS(MRCH).EQ.1) 2415.5 |
736 |
* WRITE(6,905)(WMGE(1,J),J=2,JM) 2416. |
737 |
905 FORMAT(' WMGE=',24F5.1) 2416.5 |
738 |
C 2417. |
739 |
DO 2230 J=2,JM 2417.5 |
740 |
COEKC=CKS 2418. |
741 |
IF(J.GT.JHALF) COEKC=CKN 2418.5 |
742 |
TEM=CONSTA*(BVF(J)+BVF(J-1))*ABS(DTDYJ(J))*HTJ(J)*HTJ(J) 2419. |
743 |
* *GRAV/(THSZM(J)+THSZM(J-1)) 2419.5 |
744 |
COEKY(J)=0. 2420. |
745 |
FVTH=FOV(J) |
746 |
c IF(FVTH.NE.0.) COEKY(J)=TEM/(FVTH*FVTH)/COEKC 2421. |
747 |
if(abs(FVTH).gt.1.e-15)COEKY(J)=TEM/(FVTH*FVTH)/COEKC |
748 |
IF(DTDYJ(J)*FVTH.GT.0.) COEKY(J)=0. 2421.1 |
749 |
c TEMH=DWV(J) 2425.5 |
750 |
TEMH=DWVA(J) 2425.5 |
751 |
c |
752 |
C Old restriction |
753 |
c |
754 |
c TEML=3.*1.4142*(BVF(J)+BVF(J-1))*DWV(J)/((FO(J)+FO(J-1))+1.E-20) 2421.5 |
755 |
c TEML=ABS(TEML) 2422. |
756 |
c TEMY=RADIUS*(.25*TWOPI-.5*ABS(LAT(J)+LAT(J-1))) 2422.5 |
757 |
c TEMK=.7071*.5*TWOPI*.25*(BVF(J)+BVF(J-1))* 2423. |
758 |
c * (SHJ(J)+SHJ(J-1))/(TEMY*ABS(FVTH)+1.E-20) 2423.5 |
759 |
c if(TEML.GT.TEMY) then |
760 |
c if(TEMK.lt.1.e19)then |
761 |
c TEMH=.5*(SHJ(J)+SHJ(J-1))/(SQRT(4.*TEMK*TEMK+1.)-1.) |
762 |
c else |
763 |
c TEMH=.5*(SHJ(J)+SHJ(J-1))/(2.*TEMK-1.) |
764 |
c end if |
765 |
c end if |
766 |
c |
767 |
C Old restriction |
768 |
c |
769 |
|
770 |
DO 2230 L=1,LZM1 2425. |
771 |
ZSTAR(J,L)=.5*(PHI(1,J,L)+PHI(1,J-1,L)-PHIS(1,J)-PHIS(1,J-1)) 2427. |
772 |
* /(GRAV*TEMH+1.E-20) 2427.5 |
773 |
EXPFAC=EXP(-ZSTAR(J,L)/COEKC) 2428. |
774 |
C TEMV=5.+45.*(5.*SIG(L)*EXP(1.-5.*SIG(L)))**4. 2428.5 |
775 |
C * +TEMV 2429.5 |
776 |
COEJL(J,L)=COEKY(J)*2. 2430. |
777 |
COEJL1(J,L)=COEKY(J)*EXPFAC 2430.5 |
778 |
2230 CONTINUE 2431. |
779 |
C 2431.5 |
780 |
C PARA. MERIDIONAL EDDY MOMENTUM TRANSPORT 2432. |
781 |
C 2432.5 |
782 |
DO 35 J=2,JM 2436. |
783 |
BETA=BETAV(J) |
784 |
FAVE=FOV(J) |
785 |
C DRY |
786 |
TEM1=.5*(BVFMJ(J)+BVFMJ(J-1)) 2437.5 |
787 |
c Nd instead of Nm Eq 23 from Part II. |
788 |
c See change in calculation of BVFM above |
789 |
c BVFMJ=(Nd)**2, while BVF=Nd |
790 |
C DRY |
791 |
TEM=FAVE*GRAV*DTDYJ(J)*2./(1. *TEM1*(THSZM(J)+THSZM(J-1)) 2438. |
792 |
* *DWV(J)+1.E-20) 2438.5 |
793 |
if(HPRNT)then |
794 |
if(J.eq.16)then |
795 |
print *,' Eddypa 2350' |
796 |
print *,FAVE,DTDYJ(J),TEM1 |
797 |
print *,THSZM(J),THSZM(J-1),DWV(J) |
798 |
print *,j,TEM |
799 |
endif |
800 |
endif |
801 |
IF (J.GT.2.AND.J.LT.JM) 2439. |
802 |
* DUDY2=(UZM(J+1)+UZM(J-1)-2.*UZM(J))/DYP(J)/DYP(J) 2439.5 |
803 |
* -(SINL(J)+SINL(J-1))/(COSL(J)+COSL(J-1))*.5* 2440. |
804 |
* (UZM(J+1)-UZM(J-1))/(RADIUS*DYU(J)) 2440.5 |
805 |
* -UZM(J)/(RADIUS*.5*(COSL(J)+COSL(J-1)))**2. 2441. |
806 |
IF (J.EQ.2) DUDY2=(UZM(3)-2.*UZM(2))/DYP(2)**2. 2441.5 |
807 |
* -(SINL(2)+SINL(1))/(COSL(2)+COSL(1))* 2442. |
808 |
* (UZM(3)-UZM(2))/(RADIUS*DYP(2)) 2442.5 |
809 |
* -UZM(2)/(RADIUS*.5*(COSL(2)+COSL(1)))**2. 2443. |
810 |
IF (J.EQ.JM) DUDY2=(UZM(JMM1)-2.*UZM(JM))/DYP(JM)**2. 2443.5 |
811 |
* -(SINL(JMM1)+SINL(JM))/(COSL(JMM1)+COSL(JM))* 2444. |
812 |
* (UZM(JM)-UZM(JMM1))/(RADIUS*DYP(JMM1)) 2444.5 |
813 |
* -UZM(JM)/(RADIUS*.5*(COSL(JMM1)+COSL(JM)))**2. 2445. |
814 |
if(HPRNT)then |
815 |
if(J.eq.14)then |
816 |
print *,' Eddypa 35 J=',J |
817 |
print *,BETAV(J),DUDY2 |
818 |
print *,UZM(J+1),UZM(J-1),UZM(J) |
819 |
endif |
820 |
endif |
821 |
DO 35 L=1,LZM1 2445.5 |
822 |
DQDY(J,L)= BETAV(J)-DUDY2 |
823 |
C DRY |
824 |
c |
825 |
c35 DQDY1(J,L)=TEM*(1.+.5*(DQSDT1(J,L)+DQSDT1(J-1,L)))+DQDY(J,L) 2446.5 |
826 |
c |
827 |
c M=0. Eq 23 from Part II. |
828 |
c |
829 |
DQDY1(J,L)=TEM+DQDY(J,L) |
830 |
CCCCC check + or * |
831 |
c |
832 |
C DRY |
833 |
35 CONTINUE |
834 |
DO 2350 J=2,JM 2447. |
835 |
TEM=0.25*(COSL(J)+COSL(J-1))**2. 2447.5 |
836 |
DO 2380 L=1,LZM1 2448. |
837 |
TEMZ(J,L)=DQDY(J,L)*TEM 2448.5 |
838 |
DQDY1(J,L)=COEJL1(J,L)*DQDY1(J,L)*TEM 2449. |
839 |
if(HPRNT)then |
840 |
if(J.eq.14)then |
841 |
print *,' Eddypa 2350' |
842 |
print *,L,TEM |
843 |
print *,DQDY(J,L),COEJL1(J,L),DQDY1(J,L) |
844 |
endif |
845 |
endif |
846 |
2380 CONTINUE 2449.5 |
847 |
2350 CONTINUE 2450. |
848 |
C**** TAKING VERTICAL AVERAGES 2450.5 |
849 |
DO 2365 J=2,JM 2451. |
850 |
SUMSIG=0. 2451.5 |
851 |
TEMJ(J)=0. 2452. |
852 |
TSZM(J)=0. 2452.5 |
853 |
DO 2366 L=LBGN,LZM1 2453. |
854 |
SUMSIG=SUMSIG+DSIG(L) 2453.5 |
855 |
TSZM(J)=TSZM(J)+DQDY1(J,L)*DSIG(L) 2454. |
856 |
2366 TEMJ(J)=TEMJ(J)+TEMZ(J,L)*DSIG(L) 2454.5 |
857 |
TSZM(J)=TSZM(J)/SUMSIG 2455. |
858 |
TEMJ(J) = TEMJ(J)/SUMSIG 2455.5 |
859 |
TSZM(J)=TSZM(J)*MOMCOEF |
860 |
2365 CONTINUE 2456. |
861 |
if(HPRNT)then |
862 |
print *,' Eddypa 2365' |
863 |
print *,'TEMJ' |
864 |
print *,TEMJ |
865 |
print *,'TSZM' |
866 |
print *,TSZM |
867 |
endif |
868 |
C**** CALCULATE THE ADJUSTMENT COEFFICIENT FOR MOMENTUM BALANCING 2456.5 |
869 |
PFRAC=.0 2457. |
870 |
YSCALE=400000./RADIUS 2457.5 |
871 |
DO 2355 JH=1,2 2458.5 |
872 |
JB=2 2459. |
873 |
JE=JM 2459.5 |
874 |
IF(JH.EQ.1) JE=JEQ-1 2460. |
875 |
IF(JH.EQ.2) JB=JEQ+1 2460.5 |
876 |
YZERO=YCUTN 2461. |
877 |
IF(JH.EQ.1) YZERO=YCUTS 2461.5 |
878 |
TEM1=0. 2462. |
879 |
TEM2=0. 2462.5 |
880 |
DO 2353 J=JB,JE 2463. |
881 |
FKB=1. 2463.5 |
882 |
DELTY=VLAT(J)-YZERO 2464. |
883 |
c NLAT=ABS(DELTY/DLAT)+1. 2464.5 |
884 |
IF(JH.EQ.1.AND.DELTY.GT.0.) 2465. |
885 |
* FKB=EXP(-RADIUS*DELTY/DSCALE) |
886 |
IF(JH.EQ.2.AND.DELTY.LT.0.) 2466. |
887 |
* FKB=EXP(RADIUS*DELTY/DSCALE) |
888 |
SP=.5*(P(1,J)+P(1,J-1)) 2467. |
889 |
TEM1=TEM1+TEMJ(J)*SP*FKB 2467.5 |
890 |
TEM2=TEM2+TSZM(J)*SP 2468. |
891 |
2353 CONTINUE 2468.5 |
892 |
DO 2354 J=JB,JE 2469. |
893 |
SP=.5*(P(1,J)+P(1,J-1)) 2469.5 |
894 |
2354 TSZM(J)=TSZM(J)-PFRAC*TEM2/(SP*(JE-JB+1)) 2470. |
895 |
COEKB=-TEM2*(1.-PFRAC)/(TEM1+1.E-20) 2470.5 |
896 |
C IF(JH.EQ.1) CKS=1./COEKB 2471. |
897 |
C IF(JH.EQ.2) CKN=1./COEKB 2471.5 |
898 |
IF(DMOD(TAU,25.).EQ.0..AND.IABS(MRCH).EQ.1) 2472. |
899 |
* WRITE(6,909) JH,COEKB,TAU 2472.5 |
900 |
DO 2358 J=JB,JE 2473. |
901 |
FKB=1. 2473.5 |
902 |
DELTY=VLAT(J)-YZERO 2474. |
903 |
c NLAT=ABS(DELTY/DLAT)+1. 2474.5 |
904 |
IF(JH.EQ.1.AND.DELTY.GT.0.) 2475. |
905 |
* FKB=EXP(-RADIUS*DELTY/DSCALE) |
906 |
IF(JH.EQ.2.AND.DELTY.LT.0.) 2476. |
907 |
* FKB=EXP(RADIUS*DELTY/DSCALE) |
908 |
CONKB=COEKB*FKB 2477. |
909 |
2358 TEMJ(J)=DYU(J)*(CONKB*TEMJ(J)+1.000*TSZM(J)) 2477.5 |
910 |
2355 CONTINUE 2478. |
911 |
909 FORMAT(10X,'JH =',I5,8X,'KB =',E16.5,5X,'TAU =',F10.2) 2478.5 |
912 |
C TOTLJ=JM-3 2479. |
913 |
DO 2451 J=1,JM 2479.5 |
914 |
2451 FD(1,J)=PT(1,J)*DXYP(J) 2480. |
915 |
FD(1,1)=2.*FD(1,1) 2480.5 |
916 |
FD(1,JM)=2.*FD(1,JM) 2481. |
917 |
DO 2452 J=2,JM 2481.5 |
918 |
RFDU=2./(FD(1,J)+FD(1,J-1)) 2482. |
919 |
DO 2452 L=1,LZM1 2482.5 |
920 |
2452 UOLD(J,L)=UT(1,J,L)*RFDU 2483. |
921 |
c |
922 |
if(uvtr)then |
923 |
c |
924 |
C**** CALCULATE U'V' CONTRIBUTION IN U-MOMENTUM EQUATION 2483.5 |
925 |
C TDISG=1.-SIGE(LZ) 2484. |
926 |
C TEM2=.36-SIGE(LZ)*SIGE(LZ) 2484.5 |
927 |
C TEM3=.216-SIGE(LZ)**3. 2485. |
928 |
C XCONT=TDISG/(TEM2*.3-TEM3/3.) 2485.5 |
929 |
if(OLDUVVERT)then |
930 |
C Old |
931 |
VVSUM=0. 2486. |
932 |
DO 2384 L=1,LM 2486.5 |
933 |
2384 VVSUM=VVSUM+VVN(L)+VVS(L) 2487. |
934 |
FSUM=0. 2487.5 |
935 |
DO 2383 L=1,LM 2488. |
936 |
2383 FSUM=FSUM+(VVN(L)+VVS(L))*DSIG(L)/VVSUM 2488.5 |
937 |
else |
938 |
C New |
939 |
VVNEWAV=0. |
940 |
ZZD=250. |
941 |
CALL VFUNCT(VVD,ZZD) |
942 |
DO L=1,LM |
943 |
ZZZ=PHITAV(L)/GRAV+250. |
944 |
if(L.eq.LM)then |
945 |
ZZU=30000. |
946 |
else |
947 |
ZZU=0.5*(PHITAV(L)+PHITAV(L+1))/GRAV+250. |
948 |
endif |
949 |
CALL VFUNCT(VVNEW,ZZZ) |
950 |
CALL VFUNCT(VVU,ZZU) |
951 |
VVNEWA(L)=0.25*(VVD+2.*VVNEW+VVU) |
952 |
VVNEWAV=VVNEWAV+VVNEWA(L)*DSIG(L) |
953 |
VVD=VVU |
954 |
ZZD=ZZU |
955 |
ENDDO |
956 |
if(abs(1.-VVNEWAV).gt.1.e-1)then |
957 |
print *,' VVNEWAV=',VVNEWAV |
958 |
endif |
959 |
endif ! OLDUVVERT |
960 |
DO 2385 L=1,LZ 2489. |
961 |
if(OLDUVVERT)then |
962 |
TEM=(VVN(L)+VVS(L))/(VVSUM*FSUM) 2489.5 |
963 |
C IF(SIG(L).LT.0.6) TEM=XCONT*SIG(L)*(0.6-SIG(L)) 2490. |
964 |
else |
965 |
TEM=VVNEWA(L)/VVNEWAV |
966 |
endif ! OLDUVVERT |
967 |
C New Aexp 049.98 |
968 |
JB=2 |
969 |
JE=JMM1 |
970 |
DO 2360 J=JB,JE |
971 |
VU(J,L)=VU(J-1,L)*P(1,J-1)*(COSL(J-1)/COSL(J))**2.+ |
972 |
* TEMJ(J)*TEM*.5*(P(1,J)+P(1,J-1))/COSL(J)**2. |
973 |
if(HPRNT)then |
974 |
c if(J.eq.JPR-1.and.L.eq.4)then |
975 |
if(L.eq.4)then |
976 |
print *,' Eddypa 2360' |
977 |
print *,VU(J-1,L),TEMJ(J),P(1,J-1) |
978 |
c print *,TEMJ(J),TEM,P(1,J) |
979 |
endif |
980 |
endif |
981 |
2360 VU(J,L)=VU(J,L)/P(1,J) |
982 |
C New |
983 |
DO 2370 J=JB,JE 2495. |
984 |
DO 2370 I=1,IX 2495.5 |
985 |
if(HPRNT)then |
986 |
if(J.eq.JPR-1.and.L.eq.4)then |
987 |
print *,' Eddypa 2370' |
988 |
print *,DTE,DXP(J),VU(J,L),P(1,J) |
989 |
print *,FLUX,UT(I,J+1,L) |
990 |
endif |
991 |
endif |
992 |
FLUX=DTE*VU(J,L)*P (1,J)*DXP(J) 2496. |
993 |
UT(I,J+1,L)=UT(I,J+1,L)+FLUX 2496.5 |
994 |
2370 UT(I,J,L)=UT(I,J,L)-FLUX 2497. |
995 |
2385 CONTINUE 2497.5 |
996 |
c |
997 |
if(HPRNT)then |
998 |
print *,' Eddypa 2385' |
999 |
print *,'UT' |
1000 |
print *,(UT(1,JPR,L),L=1,LM) |
1001 |
endif ! PRNT |
1002 |
endif |
1003 |
c end uvtr |
1004 |
C 2498. |
1005 |
if(wvtr)then |
1006 |
c |
1007 |
C PRESCRIBE W'U' BASED ON THE SEMI-SPECTRAL MODEL 2498.5 |
1008 |
C 2499. |
1009 |
C**** ADJUSTMENT OF U TO MAINTAIN TOTAL KINETIC ENERGY 2499.5 |
1010 |
C 2500. |
1011 |
DO 2454 J=2,JM 2500.5 |
1012 |
RFDU=2./(FD(1,J)+FD(1,J-1)) 2501. |
1013 |
DO 2454 L=1,LZM1 2501.5 |
1014 |
2454 UNEW(J,L)=UT(1,J,L)*RFDU 2502. |
1015 |
DO 2456 L=1,LZM1 2502.5 |
1016 |
DO 2456 JH=1,2 2503. |
1017 |
JB=2 2503.5 |
1018 |
JE=JM 2504. |
1019 |
IF(JH.EQ.1) JE=JEQM1 2504.5 |
1020 |
IF(JH.EQ.2) JB=JEQP1 2505. |
1021 |
SUM1=0. 2505.5 |
1022 |
SUM2=0. 2506. |
1023 |
DO 2458 J=JB,JE 2506.5 |
1024 |
TEM=0.5*(FD(1,J)+FD(1,J-1)) 2507. |
1025 |
SUM1=SUM1+TEM*UOLD(J,L)*UOLD(J,L) 2507.5 |
1026 |
2458 SUM2=SUM2+TEM*UNEW(J,L)*UNEW(J,L) 2508. |
1027 |
DIFF=(SUM2-SUM1)/SUM2 2508.5 |
1028 |
DO 2460 J=JB,JE 2509. |
1029 |
TEM=0.5*(FD(1,J)+FD(1,J-1)) 2509.5 |
1030 |
TEM1=TEM*UNEW(J,L)*UNEW(J,L)*(1.-DIFF) 2510. |
1031 |
IF(TEM1.LE.0.) GO TO 2465 2510.5 |
1032 |
TEM2=SQRT(TEM1/TEM) 2511. |
1033 |
IF(UNEW(J,L).GE.0.) UNEW(J,L)=TEM2 2511.5 |
1034 |
IF(UNEW(J,L).LT.0.) UNEW(J,L)=-TEM2 2512. |
1035 |
GO TO 2460 2512.5 |
1036 |
2465 UNEW(J,L)=0. 2513. |
1037 |
2460 CONTINUE 2513.5 |
1038 |
2456 CONTINUE 2514. |
1039 |
DO 2462 J=2,JM 2514.5 |
1040 |
FDU=0.5*(FD(1,J)+FD(1,J-1)) 2515. |
1041 |
DO 2462 L=1,LZM1 2515.5 |
1042 |
2462 UT(1,J,L)=UNEW(J,L)*FDU 2516. |
1043 |
if(HPRNT)then |
1044 |
print *,' Eddypa 2462' |
1045 |
print *,'UT' |
1046 |
print *,(UT(1,JPR,L),L=1,LM) |
1047 |
endif ! PRNT |
1048 |
endif |
1049 |
c end wvtr |
1050 |
C 2516.5 |
1051 |
if(uvtr)then |
1052 |
C**** CALCULATE OTHER U'V' AND V'V' CONTRIBUTIONS IN MOMENTUM EQUATIONS 2517. |
1053 |
C 2517.5 |
1054 |
DO 2490 L=1,LZM1 2518. |
1055 |
DO 2490 J=2,JMM1 2518.5 |
1056 |
FLUX=DTE*P(1,J)*DXP(J)*.5*(VV(J,L)+VV(J+1,L)) 2519. |
1057 |
VT(1,J+1,L)=VT(1,J+1,L)+FLUX 2519.5 |
1058 |
2490 VT(1,J,L)=VT(1,J,L)-FLUX 2520. |
1059 |
DO 2492 L=1,LZM1 2520.5 |
1060 |
DO 2494 J=2,JM 2521. |
1061 |
TEM=DTE*(DXP(J)-DXP(J-1))*(P(1,J)+P(1,J-1))*.5 2521.5 |
1062 |
FLUXU=.5*TEM*(VU(J,L)+VU(J-1,L)) 2522. |
1063 |
FLUXV=TEM*VV(J,L) 2522.5 |
1064 |
UT(1,J,L)=UT(1,J,L)-FLUXU 2523. |
1065 |
2494 VT(1,J,L)=VT(1,J,L)+FLUXV 2523.5 |
1066 |
2492 CONTINUE 2524. |
1067 |
if(HPRNT)then |
1068 |
print *,' Eddypa 2492' |
1069 |
print *,'UT' |
1070 |
print *,(UT(1,JPR,L),L=1,LM) |
1071 |
endif ! PRNT |
1072 |
endif |
1073 |
c end uvtr |
1074 |
C 2524.5 |
1075 |
if(uttr)then |
1076 |
C PARAMETERIZATION OF MERIDIONAL EDDY HEAT TRANSPORTS 2525. |
1077 |
C 2525.5 |
1078 |
COEKC=1. 2526. |
1079 |
DO 2245 L=1,LZ 2526.5 |
1080 |
DO 2240 J=2,JM 2527. |
1081 |
COEKD=CKS 2527.5 |
1082 |
IF(J.GT.JHALF) COEKD=CKN 2528. |
1083 |
TEM=.5*(PHI(1,J,L)+PHI(1,J-1,L)-PHIS(1,J)-PHIS(1,J-1)) 2528.5 |
1084 |
* /GRAV/(DWV(J)+1.E-20) 2529. |
1085 |
EXPFAC=EXP(-TEM/COEKD) 2529.5 |
1086 |
COEJL(J,L)=COEJL(J,L)*EXPFAC 2530. |
1087 |
IF(L.EQ.1) COEJL(J,L)=COEJL(J,L)*.3935 2530.5 |
1088 |
IF(L.EQ.2) COEJL(J,L)=COEJL(J,L)*.8347 2531. |
1089 |
VTH(J,L)=-.5*COEJL(J,L)*DTDYJ(J) 2531.5 |
1090 |
#if ( defined CPL_CHEM ) |
1091 |
! |
1092 |
! --- Get the Kt |
1093 |
! |
1094 |
! xxx =(T(1,j,l)-T(1,j-1,l))/dyv(j) |
1095 |
xxx =dthdy(j,l) |
1096 |
fkt(j,l)=-vth(j,l)/(xxx+1.e-20) |
1097 |
c fkt(j,l)=dmax1(0.0,fkt(j,l)) |
1098 |
c |
1099 |
c 012896: increase midlatitude eddy |
1100 |
c here 0.0 can be changed to minimum |
1101 |
c diffussion constant too |
1102 |
c |
1103 |
c fkt(j,l)=max(2.e6,fkt(j,l)*1.5) |
1104 |
c fkt(j,l)=max(0.0,fkt(j,l)*1.5) |
1105 |
c |
1106 |
c 111596: climate model change vth to 1.5 time large already: |
1107 |
c |
1108 |
fkt(j,l)=max(0.0,fkt(j,l)) |
1109 |
! |
1110 |
#endif |
1111 |
|
1112 |
FLUXT=DTE*VTH(J,L)*( P(1,J)+ P(1,J-1))*DXU(J)*.5 2532. |
1113 |
TT(1,J,L)=TT(1,J,L)+FLUXT 2532.5 |
1114 |
TT(1,J-1,L)=TT(1,J-1,L)-FLUXT 2533. |
1115 |
2240 CONTINUE 2533.5 |
1116 |
2245 CONTINUE 2534. |
1117 |
endif |
1118 |
c end uttr |
1119 |
C |
1120 |
C 2534.5 |
1121 |
if(uqtr)then |
1122 |
C PARAMETERIZATION OF MERIDIONAL EDDY MOISTURE TRANSPORTS 2535. |
1123 |
C 2535.5 |
1124 |
COEKC=1. 2536. |
1125 |
DO 2345 L=1,LZ 2536.5 |
1126 |
DO 2340 J=2,JM 2537. |
1127 |
IF(MRCH.EQ.0 ) GO TO 2340 2537.5 |
1128 |
VQ(J,L)=VTH(J,L)*P1000*.5*(DQSDT(J,L)+DQSDT(J-1,L)) 2538. |
1129 |
IF(DSHDY(J,L)*DTHDY(J,L).LT.0.) VQ(J,L)=0. 2538.1 |
1130 |
FLUXSH=VQ(J,L)*DTE*.5*(P(1,J)+P(1,J-1))*DXV(J) 2538.5 |
1131 |
IF(FLUXSH.GT..5* QT(1,J-1,L)) FLUXSH=.5* QT(1,J-1,L) 2539. |
1132 |
IF(FLUXSH.LT.-.5* QT(1,J,L)) FLUXSH=-.5* QT(1,J,L) 2539.5 |
1133 |
QT(1,J,L)= QT(1,J,L)+FLUXSH 2540. |
1134 |
QT(1,J-1,L)= QT(1,J-1,L)-FLUXSH 2540.5 |
1135 |
2340 CONTINUE 2541. |
1136 |
2345 CONTINUE 2541.5 |
1137 |
endif |
1138 |
c end uqtr |
1139 |
C 2542. |
1140 |
if(wttr) then |
1141 |
C PARAMETERIZE W'T' 2542.5 |
1142 |
C 2543. |
1143 |
COEF=.25 2543.5 |
1144 |
DO 4150 J=2,JMM1 2544. |
1145 |
COEKD=1. 2544.5 |
1146 |
C IF(J.GT.JHALF) COEKD=CKN 2545. |
1147 |
FCONT=ABS(FO(J)) 2545.5 |
1148 |
c TEM2=.5*(HTJ(J)+HTJ(J+1)) 2546. |
1149 |
TEM2=HTJP(J) |
1150 |
TEM3=ABS(.5*(DTDYJ(J)+DTDYJ(J+1))) 2546.5 |
1151 |
TEM=TEM3**3.*TEM2*TEM2*GRAV*GRAV/ 2547. |
1152 |
* (BVF(J)*THSZM(J)*THSZM(J)*FCONT*FCONT+1.E-20) 2547.5 |
1153 |
IF ((DTDYJ(J)+DTDYJ(J+1))*FO(J).GT.0.) TEM=0. 2547.6 |
1154 |
RFAC=(GAMAJ(J)/GAMAD)*(DTDZJ(J)+GAMAJ(J))/(DTDZJ(J)+GAMAD) 2551. |
1155 |
IF(RFAC.LT.0.) RFAC=0. 2551.01 |
1156 |
RFAC=RFAC+(TEM3*GRAV/(THSZM(J)*FCONT*BVF(J)+1.E-20))**2 2551.1 |
1157 |
IF(RFAC.LT.0.01) RFAC=0.01 2551.5 |
1158 |
ALAM=.573*SQRT(RFAC)/(1.-.427*RFAC**.302) 2552. |
1159 |
c TEMH=DWV(J)+DWV(J+1) 2553. |
1160 |
TEMH=DWVPA(J) |
1161 |
DO 4150 L=1,LZM2 2552.5 |
1162 |
TEMP=(PHI(1,J,L)+PHI(1,J,L+1)-2.*PHIS(1,J))/ 2554. |
1163 |
* (GRAV*TEMH+1.E-20) 2554.5 |
1164 |
WTH(J,L)=-CONSTA*TEM*EXP(-TEMP/COEKD)*(P(1,J)*SIGE(L+1)+PTOP) 2555. |
1165 |
* /SCLH(J,L)/P(1,J)*TEMP*(1.-.25*TEMP)/COEKD 2555.5 |
1166 |
WTH(J,L)=WTH(J,L)*(1.+ALAM/RFAC)/(1.+ALAM) 2556. |
1167 |
C IF(WTH(J,L).GT.0.) WTH(J,L)=0. 2556.5 |
1168 |
4150 CONTINUE 2557. |
1169 |
DO 4160 J=2,JMM1 2557.5 |
1170 |
SP= P(1,J) 2558. |
1171 |
DO 4160 I=1,IS 2558.5 |
1172 |
FLUXDT=0. 2559. |
1173 |
DO 4155 L=1,LZM2 2559.5 |
1174 |
FLUXT=DTE*SP*DXYP(J)*WTH(J,L) 2560. |
1175 |
TT(I,J,L)=TT(I,J,L)+(FLUXT-FLUXDT)/DSIG(L) 2560.5 |
1176 |
FLUXDT=FLUXT 2561. |
1177 |
4155 CONTINUE 2561.5 |
1178 |
TT(I,J,LZM1)=TT(I,J,LZM1)-FLUXDT/DSIG(LZM1) 2562. |
1179 |
4160 CONTINUE 2562.5 |
1180 |
endif |
1181 |
c end wttr |
1182 |
C |
1183 |
C 2563. |
1184 |
if(wqtr) then |
1185 |
C PARAMETERIZE W'Q' 2563.5 |
1186 |
C 2564. |
1187 |
IF(MRCH.ne.0) then |
1188 |
DO 3170 J=2,JMM1 2564.5 |
1189 |
SP=P(1,J) 2565. |
1190 |
BETA=BETAP(J) |
1191 |
FCONT=ABS(FO(J)) 2565.6 |
1192 |
TEM3=ABS(.5*(DTDYJ(J)+DTDYJ(J+1))) 2565.7 |
1193 |
ccc GAMAP=.5*(GAMA(J)+GAMA(J+1)) 2566. |
1194 |
C FUNK=.55*(1.+GAMAP)*(1.+GAMAP)-.25 2566.5 |
1195 |
c FUND=SHJ(J)/(1.48*GAMAP+.48) 2567. |
1196 |
FUND=SHJ(J)/(2./RCZZ*(1.+GAMAP(j))-1.) |
1197 |
RFAC=(GAMAJ(J)/GAMAD)*(DTDZJ(J)+GAMAJ(J))/(DTDZJ(J)+GAMAD) 2567.5 |
1198 |
IF(RFAC.LT.0.) RFAC=0. 2567.51 |
1199 |
RFAC=RFAC+(TEM3*GRAV/(THSZM(J)*FCONT*BVF(J)+1.E-20))**2 2567.6 |
1200 |
c IF(RFAC.LT.0.01) RFAC=0.01 2568. |
1201 |
ALAM=.573*SQRT(RFAC)/(1.-.427*RFAC**.302) 2568.5 |
1202 |
TEMPB=BETA*SHJ(J)/(FO(J)*GAMAP(j)+1.E-20) 2569. |
1203 |
IF ((DTDYJ(J)+DTDYJ(J+1))*FO(J).GT.0.) TEMPB=0. 2569.1 |
1204 |
|
1205 |
#if ( defined CPL_CHEM ) |
1206 |
! |
1207 |
! --- Get parameters: |
1208 |
! |
1209 |
beta1(j)=tempb |
1210 |
beta2(j)=rfac |
1211 |
! |
1212 |
#endif |
1213 |
|
1214 |
DO 3170 I=1,IS 2569.5 |
1215 |
FLUXDT=0. 2570. |
1216 |
FLUXDS=0. 2570.5 |
1217 |
DO 3175 L=1,LZM2 2571. |
1218 |
ZHEI=.5*(PHI(1,J,L)+PHI(1,J,L+1)-2.*PHIS(1,J))/GRAV 2572. |
1219 |
ZOVD=ZHEI/FUND 2572.5 |
1220 |
TEMDN1=.25*(DSHDY(J,L)+DSHDY(J,L+1)+DSHDY(J+1,L)+DSHDY(J+1,L+1)) 2573. |
1221 |
IF(J.GT.JHALF.AND.TEMDN1.GT.-DQDYMN) TEMDN1=-DQDYMN 2573.5 |
1222 |
IF(J.LE.JHALF.AND.TEMDN1.LT.DQDYMN) TEMDN1=DQDYMN 2574. |
1223 |
TERM1=ZOVD+.25*ZOVD*ZOVD*TEMPB*DSHDZ(J,L+1)/TEMDN1 2574.5 |
1224 |
DSHDZE=CP*(DTDZJ(J)+GAMAJ(J))*(1.-GAMAJ(J)/GAMAD)/HLAT 2575. |
1225 |
IF(DSHDZE.LT.0.) DSHDZE=0. 2575.1 |
1226 |
TEMDN2=.25*(DSHDYS(J,L)+DSHDYS(J,L+1)+DSHDYS(J+1,L)+ 2575.5 |
1227 |
* DSHDYS(J+1,L+1)) 2576. |
1228 |
IF(J.GT.JHALF.AND.TEMDN2.GT.-DQDYMN) TEMDN2=-DQDYMN 2576.5 |
1229 |
IF(J.LE.JHALF.AND.TEMDN2.LT.DQDYMN) TEMDN2=DQDYMN 2577. |
1230 |
TERM2=ZOVD+.25*ZOVD*ZOVD*TEMPB*DSHDZE/RFAC/TEMDN2 2577.5 |
1231 |
TEM=(TEMPB/(1.+ALAM))*(TERM1+ALAM*TERM2/RFAC) 2578. |
1232 |
WQ(J,L)=-.25*(VQ(J,L)+VQ(J,L+1)+VQ(J+1,L)+VQ(J+1,L+1))*TEM* 2578.5 |
1233 |
* (P(1,J)*SIGE(L+1)+PTOP)/SCLH(J,L)/P(1,J) 2579. |
1234 |
|
1235 |
#if ( defined CPL_CHEM ) |
1236 |
! |
1237 |
! --- Get parameters: |
1238 |
! |
1239 |
beta3(j,l+1)=zovd |
1240 |
c 112696 added /p(1,j) |
1241 |
beta4(j,l+1)=(p(1,j)*sige(l+1)+ptop)/sclh(j,l)/p(1,j) |
1242 |
! |
1243 |
#endif |
1244 |
|
1245 |
FLUXS=DTE*SP*DXYP(J)*WQ(J,L) 2579.5 |
1246 |
IF(FLUXS.GT..5* QT(I,J,L+1)*DSIG(L+1)) FLUXS=.5* QT(I,J,L+1)* 2580. |
1247 |
* DSIG(L+1) 2580.5 |
1248 |
IF(FLUXS.LT.-.5* QT(I,J,L)*DSIG(L)) FLUXS=-.5* QT(I,J,L)*DSIG(L) 2581. |
1249 |
QT(I,J,L)= QT(I,J,L)+(FLUXS-FLUXDS)/DSIG(L) 2581.5 |
1250 |
FLUXDS=FLUXS 2582. |
1251 |
3175 CONTINUE 2582.5 |
1252 |
QT(I,J,LZM1)= QT(I,J,LZM1)-FLUXDS/DSIG(LZM1) |
1253 |
3170 CONTINUE 2583.5 |
1254 |
endif ! MRCH |
1255 |
endif |
1256 |
c end wqtr |
1257 |
C 2584. |
1258 |
c print *,' TPRIM from eddypa' |
1259 |
c print *,TPRIM2 |
1260 |
c do j=1,jm |
1261 |
c do l=1,lm |
1262 |
c print *,TPRIM2(J,L)*(EXPBYK(SIG(L)*P(1,J)+PTOP))**2. |
1263 |
c enddo |
1264 |
c enddo |
1265 |
|
1266 |
RETURN 2584.5 |
1267 |
END 2585. |
1268 |
subroutine FLTR(X,JM,LM,JIN) |
1269 |
dimension X(JM,LM),Y(JM,LM) |
1270 |
do l=1,LM |
1271 |
do j=1,JM |
1272 |
Y(j,l)=X(j,l) |
1273 |
enddo |
1274 |
enddo |
1275 |
do l=1,LM |
1276 |
do j=JIN+1,JM-1 |
1277 |
X(j,l)=0.25*(X(j-1,l)+2.0*X(j,l)+X(j+1,l)) |
1278 |
enddo |
1279 |
X(JIN,l)=(2.0*X(JIN,l)+X(JIN+1,l))/3.0 |
1280 |
X(JM,l)=(2.0*X(JM,l)+X(JM-1,l))/3.0 |
1281 |
enddo |
1282 |
return |
1283 |
end |