/[MITgcm]/MITgcm_contrib/jscott/igsm/src/eddypa.F
ViewVC logotype

Contents of /MITgcm_contrib/jscott/igsm/src/eddypa.F

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


Revision 1.1 - (show annotations) (download)
Fri Aug 11 19:35:30 2006 UTC (18 years, 11 months ago) by jscott
Branch: MAIN
CVS Tags: HEAD
atm2d package

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

  ViewVC Help
Powered by ViewVC 1.1.22