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

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

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


Revision 1.2 - (show annotations) (download)
Mon Apr 23 21:20:18 2007 UTC (18 years, 3 months ago) by jscott
Branch: MAIN
CVS Tags: HEAD
Changes since 1.1: +48 -4 lines
bring igsm atmos code up to date

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_CLM
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 CLM )
38 #if ( defined CPL_CHEM )
39 !
40 #include "chem_para"
41 #include "chem_com"
42 !
43 #endif
44
45 #include "BD2G04.COM"
46
47 #include "CLM.h"
48
49 COMMON/SPEC2/KM,KINC,COEK,C3LAND(IO0,JM0),C3OICE(IO0,JM0) 5808.1
50 * ,C3LICE(IO0,JM0),WMGE(IO0,JM0),TSSFC(1,JM0,4) 5808.2
51 COMMON U,V,T,P,Q 5809.
52 COMMON/WORK1/CONV(IM0,JM0,LM0),PK(IM0,JM0,LM0),PREC(IM0,JM0),
53 & TPREC(IM0,JM0), 5810.
54 * COSZ1(IO0,JM0) 5811.
55 COMMON/WORK2/UT(IM0,JM0,LM0),VT(IM0,JM0,LM0),DU1(IO0,JM0),
56 & DV1(IO0,JM0), 5812.
57 * RA(8),ID(8),UMS(8) 5813.
58 COMMON/WORK3/E0(IO0,JM0,4),E1(IO0,JM0,4),EVAPOR(IO0,JM0,4), 5814.
59 * TGRND(IO0,JM0,4) 5814.1
60 COMMON/RDATA/ROUGHL(IO0,JM0) 5815.
61 DIMENSION SINI(72),COSI(72) 5816.
62 LOGICAL POLE,PRNT,HPRNT
63 common/conprn/HPRNT
64 ! common/TSUR/TSURFC(JM0,0:13),TSURFT(JM0),TSURFD(JM0),DTEMSR(JM0)
65 #include "TSRF.COM"
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 QSAT(TM,PR,QLH)=3.797915*EXP(QLH*(7.93252E-6-2.166847E-3/TM))/PR 5836.
81 DLQSDT(TM,QLH)=QLH*2.166847E-3/(TM*TM)
82 DATA IFIRST/1/ 5838.
83 ROSNOW(X)=0.54*X/LOG(1.+0.54*X/275.)
84 ALSNOW(X)=2.8E-6*X**2
85 C**** 5839.
86 C**** FDATA 2 LAND COVERAGE (1) 5840.
87 C**** 3 RATIO OF LAND ICE COVERAGE TO LAND COVERAGE (1) 5841.
88 C**** 5842.
89 C**** ODATA 1 OCEAN TEMPERATURE (C) 5843.
90 C**** 2 RATIO OF OCEAN ICE COVERAGE TO WATER COVERAGE (1) 5844.
91 C**** 3 OCEAN ICE AMOUNT OF SECOND LAYER (KG/M**2) 5845.
92 C**** 5846.
93 C**** GDATA 1 OCEAN ICE SNOW AMOUNT (KG/M**2) 5847.
94 C**** 2 EARTH SNOW AMOUNT (KG/M**2) 5848.
95 C**** 3 OCEAN ICE TEMPERATURE OF FIRST LAYER (C) 5849.
96 C**** 4 EARTH TEMPERATURE OF FIRST LAYER (C) 5850.
97 C**** 5 EARTH WATER OF FIRST LAYER (KG/M**2) 5851.
98 C**** 6 EARTH ICE OF FIRST LAYER (KG/M**2) 5852.
99 C**** 7 OCEAN ICE TEMPERATURE OF SECOND LAYER (C) 5853.
100 C**** 8 EARTH TEMPERATURE OF SECOND LAYER (C) 5854.
101 C**** 9 EARTH WATER OF SECOND LAYER (KG/M**2) 5855.
102 C**** 10 EARTH ICE OF SECOND LAYER (KG/M**2) 5856.
103 C**** 12 LAND ICE SNOW AMOUNT (KG/M**2) 5857.
104 C**** 13 LAND ICE TEMPERATURE OF FIRST LAYER (C) 5858.
105 C**** 14 LAND ICE TEMPERATURE OF SECOND LAYER (C) 5859.
106 C**** 5860.
107 C**** BLDATA 1 COMPOSITE SURFACE WIND MAGNITUDE (M/S) 5861.
108 C**** 2 COMPOSITE SURFACE AIR TEMPERATURE (K) 5862.
109 C**** 3 COMPOSITE SURFACE AIR SPECIFIC HUMIDITY (1) 5863.
110 C**** 4 LAYER TO WHICH DRY CONVECTION MIXES (1) 5864.
111 C**** 5 FREE 5865.
112 C**** 6 COMPOSITE SURFACE U WIND 5866.
113 C**** 7 COMPOSITE SURFACE V WIND 5867.
114 C**** 8 COMPOSITE SURFACE MOMENTUM TRANSFER (TAU) 5868.
115 C**** 5869.
116 C**** VDATA 9 WATER FIELD CAPACITY OF FIRST LAYER (KG/M**2) 5870.
117 C**** 10 WATER FIELD CAPACITY OF SECOND LAYER (KG/M**2) 5871.
118 C**** 5872.
119 C**** ROUGHL LOG(ZGS/ROUGHNESS LENGTH) (LOGARITHM TO BASE 10) 5873.
120 C**** ROUGHL will be ROUGHNESS LENGTH
121 C**** 5874.
122 c print *,'surface TAU=',TAU
123 NSTEPS=NSURF*NSTEP/NDYN 5875.
124 IF(IFIRST.NE.1) GO TO 50 5876.
125 print *,' SURFACE FOR LAND AFTER CLM'
126 IFIRST=0 5877.
127 COEFSN=100./ROSNOW(10.)
128 COEFSN=1.
129 DTSURF=NDYN*DT/NSURF 5884.
130 NGRNDZ=NGRND
131 DTGRND=DTSURF/NGRNDZ 6143.
132 NCLMPERDAY=(24.*3600.)/DTSURF
133 print *,' DTSURF=',DTSURF
134 print *,' DTGRND=',DTGRND
135 print *,' NCLMPERDAY=',NCLMPERDAY
136 SHA=RGAS/KAPA 5886.
137 RVX=0. 5887.
138 ZS1CO=.5*DSIG(1)*RGAS/GRAV 5896.
139 P1000K=EXPBYK(1000.) 5897.
140 COEFS=GRAV/(100.*DSIG(1)) 5898.
141 COEF1=(1.-SIG(2))/DSIGO(1) 5899.
142 COEF2=(SIG(1)-1.)/DSIGO(1) 5900.
143 50 CONTINUE
144 C**** ZERO OUT ENERGY AND EVAPORATION FOR GROUND AND INITIALIZE TGRND 5906.
145 DO 70 J=1,JM 5907.
146 DO 70 I=1,IM 5908.
147 TGRND(I,J,3)=GDATA(I,J,13) 5910.
148 TGRND(I,J,4)=GDATA(I,J,4) 5911.
149 DO 70 K=3,4 5912.
150 EVAPOR(I,J,K)=0.
151 E1(I,J,K)=0. 5913.
152 E0(I,J,K)=0. 5913.
153 70 CONTINUE
154 IHOUR=1.5+TOFDAY 5914.
155 C**** 5915.
156 C**** OUTSIDE LOOP OVER TIME STEPS, EXECUTED NSURF TIMES EVERY HOUR 5916.
157 C**** 5917.
158 C**** ZERO OUT LAYER 1 WIND INCREMENTS 5922.
159 DO 60 J=1,JM 5923.
160 DUL1(J)=0.
161 60 DVL1(J)=0.
162 C**** 5927.
163 C**** OUTSIDE LOOP OVER J AND I, EXECUTED ONCE FOR EACH GRID POINT 5928.
164 C**** 5929.
165 DO 7000 J=1,JM 5930.
166 if(PRNT)then
167 if(ns.eq.1)then
168 write(78,*) ,' '
169 write(78,*) ,'TAU=',TAU
170 endif
171 write(78,*),'NS=',ns
172 endif
173 100 CONTINUE
174 BTRHDT=0. 5958.
175 BSHDT=0. 5961.
176 BEVHDT=0. 5964.
177 BT2=0. 5967.
178 BTAUL=0.
179 BTAUF=0.
180 IMAX=IM 5969.
181 DO 6000 I=1,IMAX 5970.
182 C**** 5971.
183 C**** DETERMINE SURFACE CONDITIONS 5972.
184 C**** 5973.
185 PLAND=FDATA(I,J,2) 5974.
186 if(PLAND.gt.0.0)then
187 SP=P(I,J) 5980.
188 PS=SP+PTOP 5981.
189 PSK=EXPBYK(PS) 5982.
190 P1=SIG(1)*SP+PTOP 5983.
191 P1K=EXPBYK(P1) 5984.
192 C surface fluxes from radiation
193 TRHT0=TRSURF(J,2)
194 SRHEAT=SRSURF(J,2)*COSZ1(I,J)
195 C surface fluxes from radiation
196
197 RMBYA=100.*SP*DSIG(1)/GRAV 6001.
198 C**** ZERO OUT QUANTITIES TO BE SUMMED OVER SURFACE TYPES 6002.
199 TAUS=0. 6008.
200 T2MS=0.
201 C**** 6032.
202 SNOW=snwdclm(i,j)
203 ACE1=h2oiclm(i,j,1)
204 WTR1=h2olclm(i,j,1)
205 IF(SNOW.GT.0.) THEN
206 ELHX=LHS
207 ELSE
208 PFROZN=ACE1/(WTR1+ACE1+1.E-20)
209 ELHX=LHE+LHM*PFROZN
210 ENDIF
211
212 SRHDT=SRHEAT*DTSURF
213 TKV=THV1*PSK 6137.
214 ZS1=ZS1CO*TKV*SP/PS 6138.
215 P1=SIG(1)*SP+PTOP 6139.
216 SHDT=0. 6144.
217 EVHDT=0. 6145.
218 TRHDT=0. 6146.
219 F1DT=0. 6147.
220 C**** LOOP OVER GROUND TIME STEPS 6148.
221 C**** CALCULATE FLUXES OF SENSIBLE HEAT, LATENT HEAT, THERMAL 6478.
222 C**** RADIATION, AND CONDUCTION HEAT (WATTS/M**2) 6479.
223 SHEAT=shfclm(i,j)
224 EVHEAT=lhfclm(i,j)
225 TG=(abs(lwuclm(i,j))/STBO)**(1./4.)
226 TG1=TG-TF
227 c print *,'From surf_clm TAU=',TAU,' J=',j
228 c print *,LHS,LHE,LHM,ELHX
229 c print *,shfclm(i,j),lhfclm(i,j)
230 c TRHEAT=TRHT0-STBO*(TG*TG)*(TG*TG)
231 TRHEAT=TRHT0+lwuclm(i,j)
232 SHDT=DTSURF*SHEAT 6505.
233 EVHDT=DTSURF*EVHEAT 6506.
234 TRHDT=DTSURF*TRHEAT
235 ! if(TAU.gt.3623.0.and.TAU.lt.3744.0) then
236 if(J.eq.-29)then
237 ! print *,'Form surf_clm J=',j,' TAU=',TAU
238 ! print *,'TG1=',TG1,' lwuclm(i,j)=',lwuclm(i,j)
239 ! print 'A5,3F7.2',' TGL ',TAU,TOFDAY,TG1
240 ! print *,'TRHT0=',TRHT0,' TRHEAT=',TRHEAT
241 ! print *,'Form surf_clm J=',j,' TAU=',TAU
242 ! print *,'shfclm(i,j)=',shfclm(i,j),
243 ! & 'lhfclm(i,j)=',lhfclm(i,j)
244 ! print *,'tref2mclm(i,j)=',tref2mclm(i,j),
245 ! & 'tgndclm(i,j)=',tgndclm(i,j)
246 ! print *,'TG1=',TG1,' lwuclm(i,j)=',lwuclm(i,j)
247 write (329) ,TAU,TOFDAY,TG1,TRHT0,lwuclm(i,j),
248 & shfclm(i,j),lhfclm(i,j)
249 endif
250 if(J.eq.-30)then
251 write (330) ,TAU,TOFDAY,TG1,TRHT0,lwuclm(i,j),
252 & shfclm(i,j),lhfclm(i,j)
253 endif
254 ! endif
255
256 c DQ1=EVHDT/(ELHX*RMBYA)
257 c EVAP=-EVHDT/ELHX
258 EVAP=vetclm(i,j)+sevclm(i,j)+cevclm(i,j)
259 !! EVAP=vetclm(i,j)
260 EVAPV=vetclm(i,j)
261 EVAPS=sevclm(i,j)
262 EVAPC=cevclm(i,j)
263 EVAP=EVAP*DTSURF
264 DQ1=-EVAP/RMBYA
265 c print *,EVHDT,SHDT,DQ1,EVAP
266
267 F0=SRHEAT+TRHEAT+SHEAT+EVHEAT 6487.
268
269 TAUL=tauxclm(i,j)
270 TAUF=tauyclm(i,j)
271 WR=SQRT(VSSL(J)**2+USSL(J)**2)/WSSL(J)
272 TAUL=WR*TAUL
273 TAUF=WR*TAUF
274 DUL1(J)=DUL1(J)+PLAND*DTGRND*TAUL*COEFS/SP
275 DVL1(J)=DVL1(J)+PLAND*DTGRND*TAUF*COEFS/SP
276 TAUYG(J)=TAUL
277 TAUXG(J)=TAUF
278 c print *,tauxclm(i,j),tauyclm(i,j)
279 c print *,PLAND,DTGRND,COEFS,SP
280 c print *,DUL1(J),DVL1(J)
281
282
283 c TH2M=tref2mclm(i,j)
284 c t2md4tem(j)=t2md4tem(j)+TH2M
285 T2M=tref2mclm(i,j)
286 t2md4tem(j)=t2md4tem(j)+T2M
287 nt2md4tem(j)=nt2md4tem(j)+1
288 F0DT=CORSR*SRHDT+TRHDT+SHDT+EVHDT 6510.
289 if (J.eq.-23) then
290 print *,' From surf_clm TAU=',TAU
291 print *,j,i,' T2M=',t2m
292 print *,tsoiclm(i,j,1),TG,tgndclm(i,j)
293 ! TG=(abs(lwuclm(i,j))/STBO)**(1./4.)
294 print *,SHEAT,EVHEAT
295 print *,TRHEAT,SRHEAT
296 endif
297 c print *,'From surface ',TAU,CORSR,SRHDT,TRHDT,SHDT,EVHDT
298 C**** ACCUMULATE SURFACE FLUXES AND PROGNOSTIC AND DIAGNOSTIC QUANTITIES6517.
299 do ITYPE=3,4
300 E0(I,J,ITYPE)=E0(I,J,ITYPE)+F0DT 6518.
301 E1(I,J,ITYPE)=E1(I,J,ITYPE)+F1DT 6519.
302 EVAPOR(I,J,ITYPE)=EVAPOR(I,J,ITYPE)+EVAP 6520.
303 TGRND(I,J,ITYPE)=TG1 6521.
304 enddo
305
306 DTH1=-SHDT*PLAND/(SHA*RMBYA*P1K)
307 DQQ1=-DQ1*PLAND
308
309
310 TAUS=TAUS+SQRT(TAUL**2+TAUF**2)*PLAND
311 T2MS=T2MS+T2M*PLAND
312 BSHDT=BSHDT+SHDT*PLAND
313 BEVHDT=BEVHDT+EVHDT*PLAND
314 BTRHDT=BTRHDT+TRHDT*PLAND
315 c BT2=BT2+(TH2M-TF)*PLAND
316 BT2=BT2+(T2M-TF)*PLAND
317 BTAUL=BTAUL+TAUL*PLAND
318 BTAUF=BTAUF+TAUF*PLAND
319
320 5000 CONTINUE
321 DT1L(J)=DTH1
322 DQ1L(J)=DQQ1
323 TAUSL(J)=TAUS
324 T2ML(J)=T2MS
325 TLANDD(J)=TLANDD(J)+(T2M-273.16)/NCLMPERDAY
326 C**** 6596.
327 C**** ACCUMULATE DIAGNOSTICS 6597.
328 C**** 6598.
329 endif
330 6000 CONTINUE
331 C**** QUANTITIES ACCUMULATED FOR SURFACE TYPE TABLES IN DIAG1 6663.
332 BLJ(J,9)=BTRHDT
333 BLJ(J,13)=BSHDT
334 BLJ(J,14)=BEVHDT
335 BLJ(J,32)=BTAUL
336 BLJ(J,33)=BTAUF
337 BLJ(J,38)=BTAUL
338 BLJ(J,26)=BT2
339 BJ(J,41)=BJ(J,41)+EVAPV*PLAND
340 BJ(J,42)=BJ(J,42)+EVAPS*PLAND
341 BJ(J,16)=BJ(J,16)+EVAPC*PLAND
342 if(J.eq.-23)then
343 print *,'TAU=',TAU,' EVHEAT=',EVHEAT
344 print *,'EVAP=',EVAP,' EVAP1=',-EVHDT/ELHX
345 endif
346 7000 CONTINUE 6677.
347 ! do J=1,JM
348 ! TLANDD(J)=TLANDD(J)+(T2ML(J)-273.16)/NCLMPERDAY
349 ! enddo
350 ! print *,' From surf_clm T2ML TAU=',TAU
351 ! print *,T2ML
352 ! print *,'TLANDD'
353 ! print *,TLANDD
354 ! print *,'NCLMPERDAY=',NCLMPERDAY
355 c write(935),TAU,ELHTG,SHTG,TAUXG,TAUYG
356 C**** 6678.
357 #endif
358 RETURN 6795.
359 990 FORMAT ('0PPBL',3I4,14F8.2) 6818.
360 991 FORMAT ('0SURFACE ',4I4,5F10.4,3F11.7) 6819.
361 992 FORMAT ('0',I2,10F10.4/23X,4F10.4,10X,2F10.4/ 6820.
362 * 33X,3F10.4,10X,2F10.4) 6821.
363 993 FORMAT ('0',I2,10F10.4/23X,7F10.4/33X,7F10.4) 6822.
364 994 FORMAT ('0',I2,11F10.4) 6823.
365 END 6824.

  ViewVC Help
Powered by ViewVC 1.1.22