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

Contents of /MITgcm_contrib/jscott/igsm/src/surf_land.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:32 2006 UTC (18 years, 11 months ago) by jscott
Branch: MAIN
atm2d package

1
2 #include "ctrparam.h"
3
4 ! ==========================================================
5 !
6 ! SURFACE.F: THIS SUBROUTINE CALCULATES THE SURFACE FLUXES
7 ! WHICH INCLUDE SENSIBLE HEAT, EVAPORATION,
8 ! THERMAL RADIATION, AND MOMENTUM DRAG. IT ALSO
9 ! CALCULATES INSTANTANEOUSLY SURFACE TEMPERATURE,
10 ! SURFACE SPECIFIC HUMIDITY, AND SURFACE WIND
11 ! COMPONENTS.
12 !
13 ! ----------------------------------------------------------
14 !
15 ! Author of Chemistry Modules: Chien Wang
16 !
17 ! ----------------------------------------------------------
18 !
19 ! Revision History:
20 !
21 ! When Who What
22 ! ---- ---------- -------
23 ! 073100 Chien Wang repack based on CliChem3 and add cpp
24 ! 092301 Chien Wang add bc and oc
25 !
26 ! ==========================================================
27
28 SUBROUTINE SURF_LAND
29
30 C**** 5802.
31 C**** THIS SUBROUTINE CALCULATES THE SURFACE FLUXES WHICH INCLUDE 5803.
32 C**** SENSIBLE HEAT, EVAPORATION, THERMAL RADIATION, AND MOMENTUM 5804.
33 C**** DRAG. IT ALSO CALCULATES INSTANTANEOUSLY SURFACE TEMPERATURE, 5805.
34 C**** SURFACE SPECIFIC HUMIDITY, AND SURFACE WIND COMPONENTS. 5806.
35 C**** 5807.
36
37 #if ( defined CPL_CHEM )
38 !
39 #include "chem_para"
40 #include "chem_com"
41 !
42 #endif
43
44 #include "BD2G04.COM"
45
46 #if ( defined CLM )
47 #include "CLM.COM"
48 #endif
49
50 COMMON/SPEC2/KM,KINC,COEK,C3LAND(IO0,JM0),C3OICE(IO0,JM0) 5808.1
51 * ,C3LICE(IO0,JM0),WMGE(IO0,JM0),TSSFC(1,JM0,4) 5808.2
52 COMMON U,V,T,P,Q 5809.
53 COMMON/WORK1/CONV(IM0,JM0,LM0),PK(IM0,JM0,LM0),PREC(IM0,JM0),
54 & TPREC(IM0,JM0), 5810.
55 * COSZ1(IO0,JM0) 5811.
56 COMMON/WORK2/UT(IM0,JM0,LM0),VT(IM0,JM0,LM0),DU1(IO0,JM0),
57 & DV1(IO0,JM0), 5812.
58 * RA(8),ID(8),UMS(8) 5813.
59 COMMON/WORK3/E0(IO0,JM0,4),E1(IO0,JM0,4),EVAPOR(IO0,JM0,4), 5814.
60 * TGRND(IO0,JM0,4) 5814.1
61 COMMON/RDATA/ROUGHL(IO0,JM0) 5815.
62 DIMENSION SINI(72),COSI(72) 5816.
63 LOGICAL POLE,PRNT,HPRNT
64 common/conprn/HPRNT
65 common/TSUR/TSURFC(JM0,0:13),TSURFT(JM0),TSURFD(JM0),DTEMSR(JM0)
66 common/SURRAD/TRSURF(JM0,4),SRSURF(JM0,4)
67 c REAL*8 B,TGV,TKV,TSV0,TSV1,TSV 5818.
68 COMMON/CWMG/WMGEA(JM0),NWMGEA(JM0),CHAVER(JM0),DTAV(JM0),DQAV(JM0)
69 & ,Z0AV(JM0),WSAV(JM0),WS0AV(JM0),TAUAV(JM0)
70 C
71 COMMON/SURFLAND/ DUL1(JM0),DVL1(JM0),DT1L(JM0),DQ1L(JM0),
72 & WSSL(JM0),T2ML(JM0),
73 & TSSL(JM0),QSSL(JM0),USSL(JM0),VSSL(JM0),TAUSL(JM0),BLJ(JM0,50)
74 & ,ELHTG(JM0),SHTG(JM0),TAUXG(JM0),TAUYG(JM0)
75 c
76 DATA RVAP/461.5/ 5819.
77 DATA SHV/0./,SHW/4185./,SHI/2060./,RHOW/1000./,RHOI/916.6/, 5820.
78 * ALAMI/2.1762/,STBO/.5672573E-7/,TF/273.16/,TFO/-1.56/ 5821.
79 DATA Z1I/.1/,Z2LI/2.9/,Z1E/.1/,Z2E/4./,RHOS/91.66/,ALAMS/.35/ 5822.
80 DIMENSION AROUGH(20),BROUGH(20),CROUGH(20),DROUGH(20),EROUGH(20) 5823.
81 DATA AROUGH/16.59,13.99,10.4,7.35,5.241,3.926,3.126,2.632,2.319, 5824.
82 *2.116,1.982,1.893,1.832,1.788,1.757,1.733,1.714,1.699,1.687,1.677/5825.
83 DATA BROUGH/3.245,1.733,0.8481,0.3899,0.1832,0.9026E-1,0.4622E-1, 5826.
84 * .241E-1,.1254E-1,.6414E-2,.3199E-2,.1549E-2,.7275E-3,.3319E-3, 5827.
85 * .1474E-3,.6392E-4,.2713E-4,.1130E-4,.4630E-5,.1868E-5/ 5828.
86 DATA CROUGH/5.111,3.088,1.682,.9239,.5626,.3994,.3282,.3017,.299 5829.
87 *,.3114,.3324,.3587,.3881,.4186,.4492,.4792,.5082,.5361,.5627, 5830.
88 * .5882/ 5831.
89 DATA DROUGH/1.24,1.02,0.806,0.682,0.661,0.771,0.797,0.895,0.994, 5832.
90 * 1.09,1.18,1.27,1.35,1.43,1.50,1.58,1.65,1.71,1.78,1.84/ 5833.
91 DATA EROUGH/0.128,0.130,0.141,0.174,0.238,0.330,0.438,0.550,0.660,5834.
92 * 0.766,0.866,0.962,1.05,1.14,1.22,1.30,1.37,1.45,1.52,1.58/ 5835.
93 QSAT(TM,PR,QLH)=3.797915*EXP(QLH*(7.93252E-6-2.166847E-3/TM))/PR 5836.
94 DLQSDT(TM,QLH)=QLH*2.166847E-3/(TM*TM)
95 c TLOG(Z0)=ALOG(.36*SQRTT/(FMAG*Z0))+2.302585*ROUGH-.08 5837.
96 DATA IFIRST/1/ 5838.
97 ROSNOW(X)=0.54*X/LOG(1.+0.54*X/275.)
98 ALSNOW(X)=2.8E-6*X**2
99 C**** 5839.
100 C**** FDATA 2 LAND COVERAGE (1) 5840.
101 C**** 3 RATIO OF LAND ICE COVERAGE TO LAND COVERAGE (1) 5841.
102 C**** 5842.
103 C**** ODATA 1 OCEAN TEMPERATURE (C) 5843.
104 C**** 2 RATIO OF OCEAN ICE COVERAGE TO WATER COVERAGE (1) 5844.
105 C**** 3 OCEAN ICE AMOUNT OF SECOND LAYER (KG/M**2) 5845.
106 C**** 5846.
107 C**** GDATA 1 OCEAN ICE SNOW AMOUNT (KG/M**2) 5847.
108 C**** 2 EARTH SNOW AMOUNT (KG/M**2) 5848.
109 C**** 3 OCEAN ICE TEMPERATURE OF FIRST LAYER (C) 5849.
110 C**** 4 EARTH TEMPERATURE OF FIRST LAYER (C) 5850.
111 C**** 5 EARTH WATER OF FIRST LAYER (KG/M**2) 5851.
112 C**** 6 EARTH ICE OF FIRST LAYER (KG/M**2) 5852.
113 C**** 7 OCEAN ICE TEMPERATURE OF SECOND LAYER (C) 5853.
114 C**** 8 EARTH TEMPERATURE OF SECOND LAYER (C) 5854.
115 C**** 9 EARTH WATER OF SECOND LAYER (KG/M**2) 5855.
116 C**** 10 EARTH ICE OF SECOND LAYER (KG/M**2) 5856.
117 C**** 12 LAND ICE SNOW AMOUNT (KG/M**2) 5857.
118 C**** 13 LAND ICE TEMPERATURE OF FIRST LAYER (C) 5858.
119 C**** 14 LAND ICE TEMPERATURE OF SECOND LAYER (C) 5859.
120 C**** 5860.
121 C**** BLDATA 1 COMPOSITE SURFACE WIND MAGNITUDE (M/S) 5861.
122 C**** 2 COMPOSITE SURFACE AIR TEMPERATURE (K) 5862.
123 C**** 3 COMPOSITE SURFACE AIR SPECIFIC HUMIDITY (1) 5863.
124 C**** 4 LAYER TO WHICH DRY CONVECTION MIXES (1) 5864.
125 C**** 5 FREE 5865.
126 C**** 6 COMPOSITE SURFACE U WIND 5866.
127 C**** 7 COMPOSITE SURFACE V WIND 5867.
128 C**** 8 COMPOSITE SURFACE MOMENTUM TRANSFER (TAU) 5868.
129 C**** 5869.
130 C**** VDATA 9 WATER FIELD CAPACITY OF FIRST LAYER (KG/M**2) 5870.
131 C**** 10 WATER FIELD CAPACITY OF SECOND LAYER (KG/M**2) 5871.
132 C**** 5872.
133 C**** ROUGHL LOG(ZGS/ROUGHNESS LENGTH) (LOGARITHM TO BASE 10) 5873.
134 C**** ROUGHL will be ROUGHNESS LENGTH
135 C**** 5874.
136 c print *,'surface TAU=',TAU
137 NSTEPS=NSURF*NSTEP/NDYN 5875.
138 IF(IFIRST.NE.1) GO TO 50 5876.
139 print *,' SURFACE FOR LAND ONLY'
140 print *,' ZGS=30 m for LAND '
141 WMGMIN=8.
142 WMGMIN=4.
143 print *,'WMGMIN 4 LAND=',WMGMIN
144 IFIRST=0 5877.
145 print *,' WMGE'
146 print 258,(WMGE(1,J),J=1,JM)
147 258 format(12f5.1)
148 ! print *,'ODATA(1,7,2)=',ODATA(1,7,2)
149 COEFSN=100./ROSNOW(10.)
150 COEFSN=1.
151 print *,' COEFSN=',COEFSN
152 do 2567 J=1,JM
153 NWMGEA(J)=0
154 WMGEA(J)=0.
155 CHAVER(J)=0.
156 DTAV(J)=0.
157 DQAV(J)=0.
158 Z0AV(J)=0.
159 WSAV(J)=0.
160 WS0AV(J)=0.
161 TAUAV(J)=0.
162 2567 CONTINUE
163 READ (519) ((ROUGHL(I,J),I=1,IO),J=1,JM) 5878.
164 DO 10 J=1,JM
165 ILAND=0.
166 SUM1=0. 5878.02
167 CONT1=0. 5878.03
168 CONT2=0.
169 DO 11 I=1,IO 5878.04
170 PLAND=C3LAND(I,J) 5878.05
171 CONT1=CONT1+PLAND 5878.06
172 ROUGHL(I,J)=10**(log10(30.)-ROUGHL(I,J))
173 C**** ROUGHL IS NOW ROUGHNESS LENGTH
174 11 SUM1=SUM1+PLAND*ROUGHL(I,J) 5878.07
175 IF(CONT1.LE.0.) GO TO 10 5878.08
176 SUM1=SUM1/CONT1 5878.09
177 DO 12 I=1,IO 5878.1
178 12 ROUGHL(I,J)=SUM1 5878.11
179 10 CONTINUE 5878.12
180 C SRCORX=1. 5878.13
181 CIAX=0.3
182 print *,' surfacen CIAX=',CIAX
183 print *,' QS=Q1, TS=T1'
184 print *,' WS=sqrt(0.75*W1+WGEM) '
185 print *,' ROUGHL'
186 print *,(ROUGHL(1,J),J=1,jm)
187 REWIND 519 5879.
188 LBLMM1=LBLM-1 5880.
189 IQ1=IM/4+1 5881.
190 IQ2=IM/2+1 5882.
191 IQ3=3*IM/4+1 5883.
192 DTSURF=NDYN*DT/NSURF 5884.
193 print *,' DTSURF=',DTSURF
194 DTSRCE=DT*NDYN 5885.
195 SHA=RGAS/KAPA 5886.
196 RVX=0. 5887.
197 ACE1I=Z1I*RHOI 5888.
198 HC1I=ACE1I*SHI 5889.
199 HC2LI=Z2LI*RHOI*SHI 5890.
200 HC1DE=Z1E*1129950. 5891.
201 HC2DE=Z2E*1129950.+3.5*.125*RHOW*3100. 5892.
202 Z1IBYL=Z1I/ALAMI 5893.
203 Z2LI3L=Z2LI/(3.*ALAMI) 5894.
204 BYRSL=1./(RHOS*ALAMS) 5895.
205 ZS1CO=.5*DSIG(1)*RGAS/GRAV 5896.
206 P1000K=EXPBYK(1000.) 5897.
207 COEFS=GRAV/(100.*DSIG(1)) 5898.
208 COEF1=(1.-SIG(2))/DSIGO(1) 5899.
209 COEF2=(SIG(1)-1.)/DSIGO(1) 5900.
210 DO 20 I=1,IM 5901.
211 SINI(I)=SIN((I-1)*TWOPI/FIM) 5902.
212 20 COSI(I)=COS((I-1)*TWOPI/FIM) 5903.
213 50 S0=S0X*1367./RSDIST 5904.
214 BYS0=RSDIST/1367. 5905.
215 C**** ZERO OUT ENERGY AND EVAPORATION FOR GROUND AND INITIALIZE TGRND 5906.
216 DO 70 J=1,JM 5907.
217 DO 70 I=1,IM 5908.
218 c TGRND(I,J,2)=GDATA(I,J,3) 5909.
219 TGRND(I,J,3)=GDATA(I,J,13) 5910.
220 TGRND(I,J,4)=GDATA(I,J,4) 5911.
221 DO 70 K=3,4 5912.
222 EVAPOR(I,J,K)=0.
223 E1(I,J,K)=0. 5913.
224 E0(I,J,K)=0. 5913.
225 70 CONTINUE
226 IHOUR=1.5+TOFDAY 5914.
227 C**** 5915.
228 C**** OUTSIDE LOOP OVER TIME STEPS, EXECUTED NSURF TIMES EVERY HOUR 5916.
229 C**** 5917.
230 DO 9000 NS=1,NSURF 5918.
231 MODDSF=MOD(NSTEPS+NS-1,NDASF) 5919.
232 C IF(MODDSF.EQ.0) IDACC(3)=IDACC(3)+1 5920.
233 MODD6=MOD(IDAY+NS,NSURF) 5921.
234 C**** ZERO OUT LAYER 1 WIND INCREMENTS 5922.
235 DO 60 J=1,JM 5923.
236 DUL1(J)=0.
237 60 DVL1(J)=0.
238 C**** 5927.
239 C**** OUTSIDE LOOP OVER J AND I, EXECUTED ONCE FOR EACH GRID POINT 5928.
240 C**** 5929.
241 JPR=-7
242 DO 7000 J=1,JM 5930.
243 ELHTG(J)=0.0
244 SHTG(J)=0.0
245 TAUXG(J)=0.0
246 TAUYG(J)=0.0
247 PRNT=j.eq.8
248 PRNT=.FALSE.
249 if(PRNT)then
250 if(ns.eq.1)then
251 write(78,*) ,' '
252 write(78,*) ,'TAU=',TAU
253 endif
254 write(78,*),'NS=',ns
255 endif
256 HEMI=1. 5931.
257 IF(J.LE.JM/2) HEMI=-1. 5932.
258 FCOR=2.*OMEGA*SINP(J) 5933.
259 FMAG=FCOR*HEMI 5934.
260 ROOT2F=SQRT(2.*FMAG) 5935.
261 IF(J.EQ.1) GO TO 80 5936.
262 IF(J.EQ.JM) GO TO 90 5937.
263 WMG0=.5*(WMGE(1,J)+WMGE(1,J+1))+.001 5937.5
264 POLE=.FALSE. 5938.
265 IMAX=IM 5939.
266 GO TO 100 5940.
267 C**** CONDITIONS AT THE SOUTH POLE 5941.
268 80 POLE=.TRUE. 5942.
269 IMAX=1 5943.
270 JVPO=2 5944.
271 RAPO=2.*RAPVN(1) 5945.
272 U1=.25*(U(1,2,1)+V(IQ1,2,1)-U(IQ2,2,1)-V(IQ3,2,1)) 5946.
273 V1=.25*(V(1,2,1)-U(IQ1,2,1)-V(IQ2,2,1)+U(IQ3,2,1)) 5947.
274 WMG0=WMGE(1,2)+.001 5947.5
275 GO TO 100 5948.
276 C**** CONDITIONS AT THE NORTH POLE 5949.
277 90 POLE=.TRUE. 5950.
278 IMAX=1 5951.
279 JVPO=JM 5952.
280 RAPO=2.*RAPVS(JM) 5953.
281 U1=.25*(U(1,JM,1)-V(IQ1,JM,1)-U(IQ2,JM,1)+V(IQ3,JM,1)) 5954.
282 V1=.25*(V(1,JM,1)+U(IQ1,JM,1)-V(IQ2,JM,1)-U(IQ3,JM,1)) 5955.
283 WMG0=WMGE(1,JM)+.001 5955.5
284 C**** ZERO OUT SURFACE DIAGNOSTICS WHICH WILL BE SUMMED OVER LONGITUDE 5956.
285 100 ATRHDT=0. 5957.
286 BTRHDT=0. 5958.
287 CTRHDT=0. 5959.
288 ASHDT=0. 5960.
289 BSHDT=0. 5961.
290 CSHDT=0. 5962.
291 AEVHDT=0. 5963.
292 BEVHDT=0. 5964.
293 CEVHDT=0. 5965.
294 ATS=0. 5966.
295 BTS=0. 5967.
296 CTS=0. 5968.
297 AT2=0. 5966.
298 BT2=0. 5967.
299 CT2=0. 5968.
300 ATAUL=0.
301 ATAUF=0.
302 BTAUL=0.
303 BTAUF=0.
304 CTAUL=0.
305 CTAUF=0.
306 AWS=0.
307 BWS=0.
308 CWS=0.
309 AWMG=0.
310 BWMG=0.
311 CWMG=0.
312 ACH=0.
313 BCH=0.
314 CCH=0.
315 IM1=IM 5969.
316 #if ( defined CLM )
317 if(NS.eq.1)then
318 tsl4clm(j)=0.0
319 qs4clm(j)=0.0
320 ps4clm(j)=0.0
321 ws4clm(j)=0.0
322 us4clm(j)=0.0
323 vs4clm(j)=0.0
324 endif
325 #endif
326 DO 6000 I=1,IMAX 5970.
327 C**** 5971.
328 C**** DETERMINE SURFACE CONDITIONS 5972.
329 C**** 5973.
330 PLAND=FDATA(I,J,2) 5974.
331 PWATER=1.-PLAND 5975.
332 PLICE=FDATA(I,J,3)*PLAND 5976.
333 PEARTH=PLAND-PLICE 5977.
334 POICE=ODATA(I,J,2)*PWATER 5978.
335 POCEAN=PWATER-POICE 5979.
336 if(POCEAN.LE.1.E-5)then
337 POCEAN=0.
338 POICE=PWATER
339 endif
340 TTOFR=PEARTH+PLICE+POICE+POCEAN
341 if(abs(TTOFR-1).gt.1.e-3)then
342 print *,' From surface TTOFR=',TTOFR
343 print *,' J=',J,' PLAND=',PLAND,' POCEAN=',POCEAN
344 print *,'POICE=',POICE,' ODATA(I,J,2)=',ODATA(I,J,2)
345 stop
346 end if
347 SP=P(I,J) 5980.
348 PS=SP+PTOP 5981.
349 PSK=EXPBYK(PS) 5982.
350 P1=SIG(1)*SP+PTOP 5983.
351 P1K=EXPBYK(P1) 5984.
352 WSOLD=BLDATA(I,J,1) 5985.
353 USOLD=BLDATA(I,J,6) 5986.
354 VSOLD=BLDATA(I,J,7) 5987.
355 TOLD=BLDATA(I,J,8) 5988.
356 SQRTT=SQRT(TOLD) 5989.
357 GKBYFW=.1296*GRAV/(FCOR*FMAG*WSOLD+1.E-20) 5990.
358 COSWS=GKBYFW*USOLD 5991.
359 SINWS=GKBYFW*VSOLD 5992.
360 IF(POLE) GO TO 1200 5993.
361 U1=.25*(U(IM1,J,1)+U(I,J,1)+U(IM1,J+1,1)+U(I,J+1,1)) 5994.
362 V1=.25*(V(IM1,J,1)+V(I,J,1)+V(IM1,J+1,1)+V(I,J+1,1)) 5995.
363 if(J.eq.JPR.or.J.eq.-12)then
364 print *,' J=',J
365 print *,'TAU=',TAU,' NS=',NS,' ITYPE=',ITYPE
366 print *,'U(I,J,1)=',U(I,J,1),' V(I,J,1)=',V(I,J,1)
367 print *,'U(I,J+1,1)=',U(I,J+1,1),' V(I,J+1,1)=',V(I,J+1,1)
368 print *,'U(IM1,J,1)=',U(IM1,J,1),' V(IM1,J,1)=',V(IM1,J,1)
369 print *,'U(IM1,J+1,1)=',U(IM1,J+1,1),
370 & ' V(IM1,J+1,1)=',V(IM1,J+1,1)
371 endif
372 1200 TH1=T(I,J,1) 5996.
373 Q1=Q(I,J,1) 5997.
374 DTH1=0.0
375 DQQ1=0.0
376 THV1=TH1*(1.+Q1*RVX) 5998.
377 c SRHEAT=SRHR(I,J,1)*COSZ1(I,J)*SRCOR 5999.
378 c SRHDT=SRHEAT*DTSURF 6000.
379 RMBYA=100.*SP*DSIG(1)/GRAV 6001.
380 C**** ZERO OUT QUANTITIES TO BE SUMMED OVER SURFACE TYPES 6002.
381 USS=0. 6003.
382 VSS=0. 6004.
383 WSS=0. 6005.
384 TSS=0. 6006.
385 T2MS=0.
386 QSS=0. 6007.
387 TAUS=0. 6008.
388 SINAPS=0. 6009.
389 COSAPS=0. 6010.
390 JR=J
391 DXYPJ=DXYP(J) 6012.
392 TG1S=0. 6013.
393 QGS=0. 6014.
394 BETAS=0. 6015.
395 TRHDTS=0. 6016.
396 SHDTS=0. 6017.
397 EVHDTS=0. 6018.
398 UGS=0. 6019.
399 VGS=0. 6020.
400 WGS=0. 6021.
401 USRS=0. 6022.
402 VSRS=0. 6023.
403 RIS1S=0. 6024.
404 RIGSS=0. 6025.
405 CDMS=0. 6026.
406 CDHS=0. 6027.
407 DGSS=0. 6028.
408 EDS1S=0. 6029.
409 PPBLS=0. 6030.
410 EVAPS=0. 6031.
411 C**** 6032.
412 2400 IF(PLAND.LE.0.) GO TO 5000 6074.
413 NGRNDZ=NGRND 6075.
414 ROUGH=ROUGHL(I,J) 6076.
415 ZGS=30. 6078.
416 ! WMGMIN=5.
417 TRHT0=TRSURF(J,2)
418 SRHEAT=SRSURF(J,2)*COSZ1(I,J)*SRCOR
419 IF(PLICE.LE.0.) GO TO 2600 6080.
420 C**** 6081.
421 C**** LAND ICE 6082.
422 C**** 6083.
423 ITYPE=3 6084.
424 PTYPE=PLICE 6085.
425 SNOW=GDATA(I,J,12) 6086.
426 TG1=TGRND(I,J,3) 6087.
427 TG2=GDATA(I,J,14) 6088.
428 if (SNOW.gt.10.)then
429 RHOS0=ROSNOW(SNOW)
430 else
431 RHOS0=275.
432 endif
433 RHOS=COEFSN*RHOS0
434 ALAMS=ALSNOW(RHOS0)
435 BYRSL=1./(RHOS*ALAMS)
436 c Z1BY6L=(Z1IBYL+SNOW*BYRSL)*.1666667 6089.
437 CDTERM=TG2 6090.
438 c CDENOM=1./(2.*Z1BY6L+Z2LI3L) 6091.
439 Z1BY2L=(Z1IBYL+SNOW*BYRSL)*0.5
440 CDENOM=1./(Z1BY2L+3.*Z2LI3L/2.)
441 HC1=HC1I+SNOW*SHI 6092.
442 BETA=1. 6093.
443 ELHX=LHS 6094.
444 GO TO 3000 6095.
445 C**** 6096.
446 2600 IF(PEARTH.LE.0.) GO TO 5000 6097.
447 C**** 6098.
448 C**** EARTH 6099.
449 C**** 6100.
450 ITYPE=4 6101.
451 PTYPE=PEARTH 6102.
452 SNOW=GDATA(I,J,2) 6103.
453 TG1=TGRND(I,J,4) 6104.
454 WTR1=GDATA(I,J,5) 6105.
455 ACE1=GDATA(I,J,6) 6106.
456 TG2=GDATA(I,J,8) 6107.
457 WTR2=GDATA(I,J,9) 6108.
458 ACE2=GDATA(I,J,10) 6109.
459 WFC1=VDATA(I,J,9) 6110.
460 WFC2=VDATA(I,J,10) 6111.
461 WTR1DRY=0.025*WFC1
462 HC1=HC1DE+WTR1*SHW+(ACE1+SNOW)*SHI 6112.
463 ALAM1D=2.+.5*(1.+2.*WTR1/WFC1) 6113.
464 ALAM2D=4. 6114.
465 RMULCH=1. 6115.
466 IF((SINP(J).GT..5).AND.(JDAY-91)*(273-JDAY).LT.0) RMULCH=.25 6116.
467 IF((SINP(J).LT.-.5).AND.(JDAY-91)*(273-JDAY).GE.0) RMULCH=.25 6117.
468 ALAM1V=RMULCH*(.4185+1.2555*WTR1/WFC1+ALAMI*ACE1/(Z1E*RHOI)) 6118.
469 ALAM3V=.8370 6119.
470 IF(TG2.LT.0.) ALAM3V=.4185+ALAMI*.15 6120.
471 ALAM2V=.125*(.4185+1.2555*WTR2/WFC2+ALAMI*ACE2/(5.*Z1E*RHOI)) 6121.
472 * +.875*ALAM3V 6122.
473 ALAM1E=VDATA(I,J,1)*ALAM1D+(1.-VDATA(I,J,1))*ALAM1V 6123.
474 ALAM2E=VDATA(I,J,1)*ALAM2D+(1.-VDATA(I,J,1))*ALAM2V 6124.
475 if (SNOW.gt.10.)then
476 RHOS0=ROSNOW(SNOW)
477 else
478 RHOS0=275.
479 endif
480 RHOS=COEFSN*RHOS0
481 ALAMS=ALSNOW(RHOS0)
482 BYRSL=1./(RHOS*ALAMS)
483 c Z1BY6L=(Z1E/ALAM1E+SNOW*BYRSL)*.1666667 6125.
484 Z1BY2L=(Z1E/ALAM1E+SNOW*BYRSL)*0.5
485 CDTERM=TG2 6126.
486 c CDENOM=1./(2.*Z1BY6L+Z2E/(3.*ALAM2E)) 6127.
487 CDENOM=1./(Z1BY2L+Z2E/(2.*ALAM2E))
488 BETA=1. 6128.
489 ELHX=LHS 6129.
490 IF(SNOW.GT.0.) GO TO 3000 6130.
491 BETA=(WTR1+ACE1)/WFC1 6131.
492 BETA=max(((WTR1+ACE1-WTR1DRY)/WFC1),0.0)
493 PFROZN=ACE1/(WTR1+ACE1+1.E-20) 6132.
494 ELHX=LHE+LHM*PFROZN 6133.
495 HC2E=HC2DE+WTR2*SHW+ACE2*SHI
496 C**** 6134.
497 C**** BOUNDARY LAYER INTERACTION 6135.
498 C**** 6136.
499 3000 continue
500 SRHDT=SRHEAT*DTSURF
501 TKV=THV1*PSK 6137.
502 ZS1=ZS1CO*TKV*SP/PS 6138.
503 P1=SIG(1)*SP+PTOP 6139.
504 DTGRND=DTSURF/NGRNDZ 6143.
505 SHDT=0. 6144.
506 EVHDT=0. 6145.
507 TRHDT=0. 6146.
508 F1DT=0. 6147.
509 C**** LOOP OVER GROUND TIME STEPS 6148.
510 DO 3600 NG=1,NGRNDZ 6149.
511 TG=TG1+TF 6150.
512 QG=QSAT(TG,PS,ELHX) 6151.
513 TGV=TG*(1.+QG*RVX) 6152.
514 ! UG=U1
515 ! VG=V1
516 UG=0.75*U1 ! 07/14/2006
517 VG=0.75*V1 ! 07/14/2006
518 W1=SQRT(U1*U1+V1*V1)
519 W1=SQRT(UG*UG+VG*VG) ! 07/14/2006
520 WS0=W1
521 WMG=WMG0+WMGMIN
522 WS=SQRT(W1**2+WMG)
523 if(J.eq.JPR)then
524 print *,' '
525 print *,'TAU=',TAU,' NS=',NS,' ITYPE=',ITYPE
526 print *,'TG=',TG,' QG=',QG
527 print *,'RVX=',RVX,' TG1=',TG1
528 endif
529
530 WG=WS
531 THS=TH1
532 QS=Q1
533 TSV=THS*PSK
534 Z0=ROUGH
535 ROUGH=log10(ZGS/ROUGH)
536 CDN=.0231/(ROUGH*ROUGH)
537 LR=ROUGH*2.-.5
538 IF(LR.GT.20) LR=20
539 IF(LR.LT.1) LR=1
540 RIGS=ZGS*GRAV*(TSV-TGV)/(TGV*WS*WS)
541 SINAP=0.
542 COSAP=1.
543 IF(RIGS.LE.0) THEN
544 C surface layer has unstable stratification
545 CIA=TWOPI*0.0625/(1.+WS*CIAX)
546 DM=SQRT((1.-AROUGH(LR)*RIGS)*(1.-BROUGH(LR)*RIGS)/
547 * (1.-CROUGH(LR)*RIGS))
548 DH=1.35*SQRT((1.-DROUGH(LR)*RIGS)/(1.-EROUGH(LR)*RIGS))
549 ELSE
550 C surface layer has stable stratification
551 CIA=TWOPI*(0.09375-0.03125/(1.+4*RIGS**2))/(1.+WS*CIAX)
552 DM=1./(1.+(11.238+89.9*RIGS)*RIGS)
553 DH=1.35/(1.+1.93*RIGS)
554 END IF
555 CDH=CDN*DM*DH
556 if(J.eq.JPR)then
557 print *,'TAU=',TAU,' NS=',NS,' ITYPE=',ITYPE
558 print *,'WS=',WS,' ZGS=',ZGS
559 print *,'DM=',DM,' DH=',DH
560 print *,'RIGS=',RIGS,' TGV=',TGV
561 endif
562 USR=COS(CIA)
563 VSR=SIN(CIA)*HEMI
564 US=(USR*UG-VSR*VG)
565 VS=(VSR*UG+USR*VG)
566 RCDHWS=CDH*WS*100.*PS/(RGAS*TSV)
567 if(J.eq.JPR)then
568 c print *,' '
569 print *,'TAU=',TAU,' NS=',NS,' ITYPE=',ITYPE
570 print *,'CDH=',CDH,' RGAS=',RGAS
571 print *,'PS=',PS,' TSV=',TSV
572 print *,'WS=',WS,' RCDHWS=',RCDHWS
573 endif
574 TS=TSV/(1.+QS*RVX) 6467.
575 QSATS=QSAT(TS,PS,ELHX) 6468.
576 IF(QS.LE.QSATS) GO TO 3500 6469.
577 DQSSDT=QSATS*ELHX/(RVAP*TS*TS) 6470.
578 X=(QS-QSATS)/(DQSSDT+(SHA/ELHX)) 6471.
579 TS=TS+X 6472.
580 QS=QS+X*(SHA/ELHX) 6473.
581 3500 CONTINUE
582
583 #if ( defined CLM )
584 if(ITYPE.EQ.4.or.ITYPE.EQ.3)then
585 tsl4clm(j)=tsl4clm(j)+TS*PTYPE/PLAND
586 qs4clm(j)=qs4clm(j)+QS*PTYPE/PLAND
587 ps4clm(j)=ps4clm(j)+PS*PTYPE/PLAND
588 ws4clm(j)=ws4clm(j)+WS*PTYPE/PLAND
589 us4clm(j)=us4clm(j)+US*PTYPE/PLAND
590 vs4clm(j)=vs4clm(j)+VS*PTYPE/PLAND
591 endif
592 #endif
593
594 C**** CALCULATE RHOS*CDM*WS AND RHOS*CDH*WS 6474.
595 CDM=CDN*DM 6475.
596 RCDMWS=CDM*WS*100.*PS/(RGAS*TS) 6476.
597 C**** CALCULATE FLUXES OF SENSIBLE HEAT, LATENT HEAT, THERMAL 6478.
598 C**** RADIATION, AND CONDUCTION HEAT (WATTS/M**2) 6479.
599 SHEAT=SHA*RCDHWS*(TS-TG) 6480.
600 BETAUP=BETA 6481.
601 IF(QS.GT.QG) BETAUP=1. 6482.
602 EVHEAT=(LHE+TG1*SHV)*BETAUP*RCDHWS*(QS-QG) 6483.
603 c TRHEAT=TRHR(I,J,1)-STBO*(TG*TG)*(TG*TG) 6484.
604 TRHEAT=TRHT0-STBO*(TG*TG)*(TG*TG)
605 if(J.eq.JPR)then
606 c print *,' '
607 print *,'TAU=',TAU,' NS=',NS,' ITYPE=',ITYPE
608 print *,'TRHT0=',TRHT0,' STBO=',STBO
609 print *,'TG=',TG,' TS=',TS
610 print *,'TRHEAT=',TRHEAT
611 print *,'SHA=',SHA,' RCDHWS=',RCDHWS
612 print *,'SHEAT=',SHEAT
613 endif
614 TG1OLD=TG1
615 SHEATOLD=SHEAT
616 C**** CALCULATE FLUXES USING IMPLICIT TIME STEP FOR NON-OCEAN POINTS 6486.
617 F0=SRHEAT+TRHEAT+SHEAT+EVHEAT 6487.
618 c F1=(TG1-CDTERM-F0*Z1BY6L)*CDENOM 6488.
619 F1=(TG1-CDTERM)*CDENOM
620 DSHDTG=-RCDHWS*SHA
621 DQGDTG=QG*ELHX/(RVAP*TG*TG) 6490.
622 DEVDTG=-RCDHWS*LHE*BETAUP*DQGDTG
623 DTRDTG=-4.*STBO*TG*TG*TG 6492.
624 DF0DTG=DSHDTG+DEVDTG+DTRDTG 6493.
625 c DFDTG=DF0DTG-(1.-DF0DTG*Z1BY6L)*CDENOM 6493.5
626 DFDTG=DF0DTG-CDENOM
627 c DF1DTG=(1.-DF0DTG*Z1BY6L)*CDENOM
628 DF1DTG=CDENOM
629 DTG=(F0-F1)*DTGRND/(HC1-DTGRND*DFDTG) 6494.
630 SHDT=SHDT+DTGRND*(SHEAT+DTG*DSHDTG) 6495.
631 EVHDT=EVHDT+DTGRND*(EVHEAT+DTG*DEVDTG) 6496.
632 TRHDT=TRHDT+DTGRND*(TRHEAT+DTG*DTRDTG) 6497.
633 TG1=TG1+DTG
634 c F1DT=F1DT+DTGRND*(TG1-CDTERM-(F0+DTG*DF0DTG)*Z1BY6L)*CDENOM 6498.
635 F1DT=F1DT+DTGRND*(TG1-CDTERM)*CDENOM
636 DUL1(J)=DUL1(J)+PTYPE*DTGRND*RCDMWS*US*COEFS/SP
637 DVL1(J)=DVL1(J)+PTYPE*DTGRND*RCDMWS*VS*COEFS/SP
638 c TG1=TG1+DTG 6501.
639 3600 CONTINUE 6502.
640 EPS=1.D-8
641 c print *,'FROM SURFACE NS=',NS
642 c print *,'J=',J,' ITYPE=',ITYPE
643 c print *,RCDMWS,WS
644 WWS=max(W1,1.D-4)
645 c RO=SP*100/(RGAS*TG)
646 c print *,'RO=',RO
647 c USTAR=SQRT(RCDMWS*WS/RO)
648 c TSTAR=SHEATOLD/(0.35*1007.*RO*USTAR)
649 c ALPHAH=DH
650 c TT2M=TG+TSTAR/ALPHAH*LOG(2.0/Z0)
651 c TT2M=TG+TSTAR/ALPHAH*LOG(ZGS/Z0)
652 c print *,'RIGS=',RIGS,' Z0=',Z0
653 c print *,'CDN=',CDN
654 c print *,'H=',SHDT/DTSURF,' TGM=',RCDMWS*WS
655 c print *,' SHEATOLD=',SHEATOLD
656 c print *,' USTAR=',USTAR,' TSTAR=',TSTAR
657 c print *,' ALPHAH=',ALPHAH,' TT2M=',TT2M
658 c print *,' TT2M=',TT2M
659 ZTEM=ZGS
660 ZTEM=2.0
661 c print *,'ZTEM=',ZTEM
662 CALL TZM(T2M,TH2M,ZTEM,Z0,ZGS,SP,TG,TS,RIGS,WS,
663 & -SHEATOLD,RCDMWS*WS,LR,EPS)
664 c print *,'FROM SURFACE'
665 c print *,'TS=',TS,' TG=',TG
666 c print *,' T2M=',T2M,' TH2M=',TH2M
667 F0DT=CORSR*SRHDT+TRHDT+SHDT+EVHDT 6510.
668 if(J.eq.JPR)then
669 print *,'TAU=',TAU,' NS=',NS,' ITYPE=',ITYPE
670 print *,'DTSURF=',DTSURF,' CORSR=',CORSR
671 print *,'SRHDT=',SRHDT,' TRHDT=',TRHDT
672 print *,'SHDT=',SHDT,' EVHDT=',EVHDT
673 print *,'F0DT=',F0DT
674 print *,'US=',US,' VS=',VS
675 print *,'COEFS=',COEFS,' SP=',SP
676 endif
677 c print *,'From surface ',TAU,CORSR,SRHDT,TRHDT,SHDT,EVHDT
678 C**** CALCULATE EVAPORATION 6511.
679 CCC DQ1=EVHDT/((LHE+TG1*SHV)*RMBYA) 6512.
680 DQ1=EVHDT/(ELHX*RMBYA)
681 IF(DQ1*PTYPE.LE.Q1) GO TO 3720 6513.
682 DQ1=Q1/PTYPE 6514.
683 CCC EVHDT=DQ1*(LHE+TG1*SHV)*RMBYA 6515.
684 EVHDT=DQ1*ELHX*RMBYA
685 3720 EVAP=-DQ1*RMBYA 6516.
686 C**** ACCUMULATE SURFACE FLUXES AND PROGNOSTIC AND DIAGNOSTIC QUANTITIES6517.
687 E0(I,J,ITYPE)=E0(I,J,ITYPE)+F0DT 6518.
688 E1(I,J,ITYPE)=E1(I,J,ITYPE)+F1DT 6519.
689 EVAPOR(I,J,ITYPE)=EVAPOR(I,J,ITYPE)+EVAP 6520.
690 if(PRNT)then
691 c write(78,*) ,' '
692 c write(78,*) ,'TAU=',TAU
693 write(78,*) ,'J=',j,' ITYPE=',ITYPE,' PTYPE=',PTYPE
694 write(78,*) ,'TS=',TS,' TG=',TG,' QS=',QS
695 write(78,*) ,'TG1=',TG1,' TG1OLD=',TG1OLD
696 write(78,*) ,'TG2=',TG2
697 write(78,*) ,'SHEAT=',SHEAT,' EVHEAT=',EVHEAT
698 write(78,*) ,'TRHEAT=',TRHEAT,' SRHEAT=',SRHEAT
699 write(78,*) ,'EVAP mm/day=',24.*3600.*EVAP/DTSURF
700 write(78,*) ,'EVAP=',EVAP,
701 & ' F0DT=',F0DT/DTSURF,' F1DT=',F1DT/DTSURF
702 endif
703 TGRND(I,J,ITYPE)=TG1 6521.
704 TSSFC(I,J,ITYPE)=TS 6521.5
705 c TH1=TH1-SHDT*PTYPE/(SHA*RMBYA*P1K) 6522.
706 c Q1=Q1-DQ1*PTYPE 6523.
707
708 DTH1=DTH1-SHDT*PTYPE/(SHA*RMBYA*P1K)
709 DQQ1=DQQ1-DQ1*PTYPE
710
711
712 USS=USS+US*PTYPE 6524.
713 VSS=VSS+VS*PTYPE 6525.
714 WSS=WSS+WS*PTYPE 6526.
715 TSS=TSS+TS*PTYPE 6527.
716 T2MS=T2MS+T2M*PTYPE
717 QSS=QSS+QS*PTYPE 6528.
718 ELHTG(J)=ELHTG(J)+EVHEAT*PTYPE
719 SHTG(J)=SHTG(J)+SHEAT*PTYPE
720 TAUXG(J)=TAUXG(J)+RCDMWS*US*PTYPE
721 TAUYG(J)=TAUYG(J)+RCDMWS*VS*PTYPE
722 ! TAUS=TAUS+CDM*WS*W1*PTYPE 6529.
723 TAUS=TAUS+CDM*WS*WS*PTYPE ! 7/14/2006
724 SINAPS=SINAPS+SINAP*PTYPE 6530.
725 COSAPS=COSAPS+COSAP*PTYPE 6531.
726 TG1S=TG1S+TG1*PTYPE 6532.
727 QGS=QGS+QG*PTYPE 6533.
728 BETAS=BETAS+BETA*PTYPE 6534.
729 TRHDTS=TRHDTS+TRHDT*PTYPE 6535.
730 SHDTS=SHDTS+SHDT*PTYPE 6536.
731 EVHDTS=EVHDTS+EVHDT*PTYPE 6537.
732 UGS=UGS+UG*PTYPE 6538.
733 VGS=VGS+VG*PTYPE 6539.
734 WGS=WGS+WG*PTYPE 6540.
735 USRS=USRS+USR*PTYPE 6541.
736 VSRS=VSRS+VSR*PTYPE 6542.
737 c RIS1S=RIS1S+RIS1*PTYPE 6543.
738 RIGSS=RIGSS+RIGS*PTYPE 6544.
739 CDMS=CDMS+CDM*PTYPE 6545.
740 CDHS=CDHS+CDH*PTYPE 6546.
741 c DGSS=DGSS+DGS*PTYPE 6547.
742 c EDS1S=EDS1S+EDS1*PTYPE 6548.
743 c PPBLS=PPBLS+PPBL*PTYPE 6549.
744 EVAPS=EVAPS+EVAP*PTYPE 6550.
745 GO TO (5000,5000,4400,4600),ITYPE 6551.
746 C**** 6552.
747 C**** LAND ICE 6569.
748 C**** 6570.
749 4400 BSHDT=BSHDT+SHDT*PLICE 6571.
750 BEVHDT=BEVHDT+EVHDT*PLICE 6572.
751 BTRHDT=BTRHDT+TRHDT*PLICE 6573.
752 BTS=BTS+(TS-TF)*PLICE 6574.
753 c BT2=BT2+(TH2M-TF)*PLICE
754 BT2=BT2+(T2M-TF)*PLICE
755 BTAUL=BTAUL+RCDMWS*US*PLICE
756 BTAUF=BTAUF+RCDMWS*VS*PLICE
757 BWS=BWS+WS*PLICE
758 BWMG=BWMG+SQRT(WMG)*PLICE
759 BCH=BCH+RCDHWS*PLICE
760 GO TO 2600 6575.
761 C**** 6576.
762 C**** EARTH 6577.
763 C**** 6578.
764 4600 BSHDT=BSHDT+SHDT*PEARTH 6579.
765 BEVHDT=BEVHDT+EVHDT*PEARTH 6580.
766 BTRHDT=BTRHDT+TRHDT*PEARTH 6581.
767 BTS=BTS+(TS-TF)*PEARTH 6582.
768 c BT2=BT2+(TH2M-TF)*PEARTH
769 BT2=BT2+(T2M-TF)*PEARTH
770 BTAUL=BTAUL+RCDMWS*US*PEARTH
771 BTAUF=BTAUF+RCDMWS*VS*PEARTH
772 BWS=BWS+WS*PEARTH
773 BWMG=BWMG+SQRT(WMG)*PEARTH
774 BCH=BCH+RCDHWS*PEARTH
775
776 C**** NON-OCEAN POINTS WHICH ARE NOT MELTING OR FREEZING WATER USE 6583.
777 C**** IMPLICIT TIME STEPS 6584.
778 C**** 6585.
779 C**** UPDATE SURFACE AND FIRST LAYER QUANTITIES 6586.
780 C**** 6587.
781 5000 CONTINUE
782 DT1L(J)=DTH1
783 DQ1L(J)=DQQ1
784 WSSL(J)=WSS
785 TSSL(J)=TSS
786 T2ML(J)=T2MS
787 QSSL(J)=QSS
788 USSL(J)=USS
789 VSSL(J)=VSS
790 TAUSL(J)=TAUS
791 ELHTG(J)=ELHTG(J)/PLAND
792 SHTG(J)=SHTG(J)/PLAND
793 TAUXG(J)=TAUXG(J)/PLAND
794 TAUYG(J)=TAUYG(J)/PLAND
795 C**** 6596.
796 C**** ACCUMULATE DIAGNOSTICS 6597.
797 C**** 6598.
798 6000 IM1=I 6662.
799 C**** QUANTITIES ACCUMULATED FOR SURFACE TYPE TABLES IN DIAG1 6663.
800 BLJ(J,9)=BTRHDT
801 BLJ(J,13)=BSHDT
802 BLJ(J,14)=BEVHDT
803 BLJ(J,32)=BTAUL
804 BLJ(J,33)=BTAUF
805 BLJ(J,37)=BWS
806 BLJ(J,28)=BWMG
807 BLJ(J,38)=BTAUL
808 IF(MODDSF.NE.0) GO TO 7000 6673.
809 BLJ(J,23)=BTS
810 BLJ(J,26)=BT2
811 7000 CONTINUE 6677.
812 9000 CONTINUE
813 c write(935),TAU,ELHTG,SHTG,TAUXG,TAUYG
814 C**** 6678.
815 RETURN 6795.
816 990 FORMAT ('0PPBL',3I4,14F8.2) 6818.
817 991 FORMAT ('0SURFACE ',4I4,5F10.4,3F11.7) 6819.
818 992 FORMAT ('0',I2,10F10.4/23X,4F10.4,10X,2F10.4/ 6820.
819 * 33X,3F10.4,10X,2F10.4) 6821.
820 993 FORMAT ('0',I2,10F10.4/23X,7F10.4/33X,7F10.4) 6822.
821 994 FORMAT ('0',I2,11F10.4) 6823.
822 END 6824.

  ViewVC Help
Powered by ViewVC 1.1.22