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

Annotation of /MITgcm_contrib/jscott/igsm/src/surface.F

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


Revision 1.3 - (hide 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.2: +23 -14 lines
bring igsm atmos code up to date

1 jscott 1.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 SURFCE 5801.
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 jscott 1.3 #include "CLM.h"
48 jscott 1.1 #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     DIMENSION WMGMINO(JM0)
64     LOGICAL POLE,PRNT,HPRNT
65     common/conprn/HPRNT
66 jscott 1.3 ! common/TSUR/TSURFC(JM0,0:13),TSURFT(JM0),TSURFD(JM0),DTEMSR(JM0)
67     #include "TSRF.COM"
68 jscott 1.1 common/SURRAD/TRSURF(JM0,4),SRSURF(JM0,4)
69     c REAL*8 B,TGV,TKV,TSV0,TSV1,TSV 5818.
70     COMMON/CWMG/WMGEA(JM0),NWMGEA(JM0),CHAVER(JM0),DTAV(JM0),DQAV(JM0)
71     & ,Z0AV(JM0),WSAV(JM0),WS0AV(JM0),TAUAV(JM0)
72     C
73     #if ( defined OCEAN_3D || defined ML_2D)
74 jscott 1.2 #include "AGRID.h"
75 jscott 1.1 #endif
76     c
77     DATA RVAP/461.5/ 5819.
78     DATA SHV/0./,SHW/4185./,SHI/2060./,RHOW/1000./,RHOI/916.6/, 5820.
79     * ALAMI/2.1762/,STBO/.5672573E-7/,TF/273.16/,TFO/-1.56/ 5821.
80     DATA Z1I/.1/,Z2LI/2.9/,Z1E/.1/,Z2E/4./,RHOS/91.66/,ALAMS/.35/ 5822.
81     DIMENSION AROUGH(20),BROUGH(20),CROUGH(20),DROUGH(20),EROUGH(20) 5823.
82     DATA AROUGH/16.59,13.99,10.4,7.35,5.241,3.926,3.126,2.632,2.319, 5824.
83     *2.116,1.982,1.893,1.832,1.788,1.757,1.733,1.714,1.699,1.687,1.677/5825.
84     DATA BROUGH/3.245,1.733,0.8481,0.3899,0.1832,0.9026E-1,0.4622E-1, 5826.
85     * .241E-1,.1254E-1,.6414E-2,.3199E-2,.1549E-2,.7275E-3,.3319E-3, 5827.
86     * .1474E-3,.6392E-4,.2713E-4,.1130E-4,.4630E-5,.1868E-5/ 5828.
87     DATA CROUGH/5.111,3.088,1.682,.9239,.5626,.3994,.3282,.3017,.299 5829.
88     *,.3114,.3324,.3587,.3881,.4186,.4492,.4792,.5082,.5361,.5627, 5830.
89     * .5882/ 5831.
90     DATA DROUGH/1.24,1.02,0.806,0.682,0.661,0.771,0.797,0.895,0.994, 5832.
91     * 1.09,1.18,1.27,1.35,1.43,1.50,1.58,1.65,1.71,1.78,1.84/ 5833.
92     DATA EROUGH/0.128,0.130,0.141,0.174,0.238,0.330,0.438,0.550,0.660,5834.
93     * 0.766,0.866,0.962,1.05,1.14,1.22,1.30,1.37,1.45,1.52,1.58/ 5835.
94     QSAT(TM,PR,QLH)=3.797915*EXP(QLH*(7.93252E-6-2.166847E-3/TM))/PR 5836.
95     DLQSDT(TM,QLH)=QLH*2.166847E-3/(TM*TM)
96     c TLOG(Z0)=ALOG(.36*SQRTT/(FMAG*Z0))+2.302585*ROUGH-.08 5837.
97     DATA IFIRST/1/ 5838.
98     ROSNOW(X)=0.54*X/LOG(1.+0.54*X/275.)
99     ALSNOW(X)=2.8E-6*X**2
100     C**** 5839.
101     C**** FDATA 2 LAND COVERAGE (1) 5840.
102     C**** 3 RATIO OF LAND ICE COVERAGE TO LAND COVERAGE (1) 5841.
103     C**** 5842.
104     C**** ODATA 1 OCEAN TEMPERATURE (C) 5843.
105     C**** 2 RATIO OF OCEAN ICE COVERAGE TO WATER COVERAGE (1) 5844.
106     C**** 3 OCEAN ICE AMOUNT OF SECOND LAYER (KG/M**2) 5845.
107     C**** 5846.
108     C**** GDATA 1 OCEAN ICE SNOW AMOUNT (KG/M**2) 5847.
109     C**** 2 EARTH SNOW AMOUNT (KG/M**2) 5848.
110     C**** 3 OCEAN ICE TEMPERATURE OF FIRST LAYER (C) 5849.
111     C**** 4 EARTH TEMPERATURE OF FIRST LAYER (C) 5850.
112     C**** 5 EARTH WATER OF FIRST LAYER (KG/M**2) 5851.
113     C**** 6 EARTH ICE OF FIRST LAYER (KG/M**2) 5852.
114     C**** 7 OCEAN ICE TEMPERATURE OF SECOND LAYER (C) 5853.
115     C**** 8 EARTH TEMPERATURE OF SECOND LAYER (C) 5854.
116     C**** 9 EARTH WATER OF SECOND LAYER (KG/M**2) 5855.
117     C**** 10 EARTH ICE OF SECOND LAYER (KG/M**2) 5856.
118     C**** 12 LAND ICE SNOW AMOUNT (KG/M**2) 5857.
119     C**** 13 LAND ICE TEMPERATURE OF FIRST LAYER (C) 5858.
120     C**** 14 LAND ICE TEMPERATURE OF SECOND LAYER (C) 5859.
121     C**** 5860.
122     C**** BLDATA 1 COMPOSITE SURFACE WIND MAGNITUDE (M/S) 5861.
123     C**** 2 COMPOSITE SURFACE AIR TEMPERATURE (K) 5862.
124     C**** 3 COMPOSITE SURFACE AIR SPECIFIC HUMIDITY (1) 5863.
125     C**** 4 LAYER TO WHICH DRY CONVECTION MIXES (1) 5864.
126     C**** 5 SURFACE MOMENTUM TRANSFER (TAU) OCEAN 5865.
127     C**** 6 COMPOSITE SURFACE U WIND 5866.
128     C**** 7 COMPOSITE SURFACE V WIND 5867.
129     C**** 8 COMPOSITE SURFACE MOMENTUM TRANSFER (TAU) 5868.
130     C**** 5869.
131     C**** VDATA 9 WATER FIELD CAPACITY OF FIRST LAYER (KG/M**2) 5870.
132     C**** 10 WATER FIELD CAPACITY OF SECOND LAYER (KG/M**2) 5871.
133     C**** 5872.
134     C**** ROUGHL LOG(ZGS/ROUGHNESS LENGTH) (LOGARITHM TO BASE 10) 5873.
135     C**** ROUGHL will be ROUGHNESS LENGTH
136     C**** 5874.
137     c print *,'surface TAU=',TAU
138     NSTEPS=NSURF*NSTEP/NDYN 5875.
139     IF(IFIRST.NE.1) GO TO 50 5876.
140     print *,' SURFACE CORSR=',CORSR
141     print *,' ZGS=30 m for LAND '
142     IFIRST=0 5877.
143     WMGMINL = 5.0
144     print *,'WMGMIN 4 LAND=',WMGMINL
145     print *,'over land WMG=max(WMG0,WMGMIN)'
146     WMGM0=8.0
147     WMGM45=25.
148     print *,' WMGM0=', WMGM0,' WMGM45=',WMGM45
149     WMGMAV=0.5*( WMGM0+WMGM45)
150     DWGM=0.5*( WMGM0-WMGM45)
151     do j = 1,jm0
152     rhrad = 3.14159*(-90.+4.*(j-1))/180.
153     WMGMINO(j) = WMGMAV+DWGM*cos(4.*rhrad)
154     enddo
155     print *,'WMGMIN 4 OCEAN is a function of latitude'
156     print 258,(WMGMINO(J),J=1,JM)
157     print *,' WMGE'
158     print 258,(WMGE(1,J),J=1,JM)
159     258 format(12f5.1)
160     ! print *,'ODATA(1,7,2)=',ODATA(1,7,2)
161     COEFSN=100./ROSNOW(10.)
162     COEFSN=1.
163     print *,' COEFSN=',COEFSN
164     do 2567 J=1,JM
165     NWMGEA(J)=0
166     WMGEA(J)=0.
167     CHAVER(J)=0.
168     DTAV(J)=0.
169     DQAV(J)=0.
170     Z0AV(J)=0.
171     WSAV(J)=0.
172     WS0AV(J)=0.
173     TAUAV(J)=0.
174     2567 CONTINUE
175     READ (519) ((ROUGHL(I,J),I=1,IO),J=1,JM) 5878.
176     c DO 10 J=2,JMM1 5878.01
177     C *************
178     DO 10 J=1,JM
179     C *************
180     ILAND=0.
181     SUM1=0. 5878.02
182     CONT1=0. 5878.03
183     CONT2=0.
184     DO 11 I=1,IO 5878.04
185     PLAND=C3LAND(I,J) 5878.05
186     CONT1=CONT1+PLAND 5878.06
187     ROUGHL(I,J)=10**(log10(30.)-ROUGHL(I,J))
188     C**** ROUGHL IS NOW ROUGHNESS LENGTH
189     11 SUM1=SUM1+PLAND*ROUGHL(I,J) 5878.07
190     IF(CONT1.LE.0.) GO TO 10 5878.08
191     SUM1=SUM1/CONT1 5878.09
192     DO 12 I=1,IO 5878.1
193     12 ROUGHL(I,J)=SUM1 5878.11
194     10 CONTINUE 5878.12
195     C SRCORX=1. 5878.13
196     CIAX=0.3
197     print *,' surfacen CIAX=',CIAX
198     print *,' QS=Q1, TS=T1'
199     print *,' WS=sqrt(0.75*W1+WGEM) '
200     print *,' ROUGHL'
201     print *,(ROUGHL(1,J),J=1,jm)
202     REWIND 519 5879.
203     LBLMM1=LBLM-1 5880.
204     IQ1=IM/4+1 5881.
205     IQ2=IM/2+1 5882.
206     IQ3=3*IM/4+1 5883.
207     DTSURF=NDYN*DT/NSURF 5884.
208     print *,' DTSURF=',DTSURF
209     DTSRCE=DT*NDYN 5885.
210     SHA=RGAS/KAPA 5886.
211     RVX=0. 5887.
212     ACE1I=Z1I*RHOI 5888.
213     HC1I=ACE1I*SHI 5889.
214     HC2LI=Z2LI*RHOI*SHI 5890.
215     HC1DE=Z1E*1129950. 5891.
216     HC2DE=Z2E*1129950.+3.5*.125*RHOW*3100. 5892.
217     Z1IBYL=Z1I/ALAMI 5893.
218     Z2LI3L=Z2LI/(3.*ALAMI) 5894.
219     BYRSL=1./(RHOS*ALAMS) 5895.
220     ZS1CO=.5*DSIG(1)*RGAS/GRAV 5896.
221     P1000K=EXPBYK(1000.) 5897.
222     COEFS=GRAV/(100.*DSIG(1)) 5898.
223     COEF1=(1.-SIG(2))/DSIGO(1) 5899.
224     COEF2=(SIG(1)-1.)/DSIGO(1) 5900.
225     DO 20 I=1,IM 5901.
226     SINI(I)=SIN((I-1)*TWOPI/FIM) 5902.
227     20 COSI(I)=COS((I-1)*TWOPI/FIM) 5903.
228     50 S0=S0X*1367./RSDIST 5904.
229     BYS0=RSDIST/1367. 5905.
230     C**** ZERO OUT ENERGY AND EVAPORATION FOR GROUND AND INITIALIZE TGRND 5906.
231     DO 70 J=1,JM 5907.
232     DO 70 I=1,IM 5908.
233     TGRND(I,J,2)=GDATA(I,J,3) 5909.
234     TGRND(I,J,3)=GDATA(I,J,13) 5910.
235     TGRND(I,J,4)=GDATA(I,J,4) 5911.
236     DO 70 K=1,12 5912.
237     70 E0(I,J,K)=0. 5913.
238     IHOUR=1.5+TOFDAY 5914.
239     C**** 5915.
240     C**** OUTSIDE LOOP OVER TIME STEPS, EXECUTED NSURF TIMES EVERY HOUR 5916.
241     C**** 5917.
242     DO 9000 NS=1,NSURF 5918.
243     MODDSF=MOD(NSTEPS+NS-1,NDASF) 5919.
244     IF(MODDSF.EQ.0) IDACC(3)=IDACC(3)+1 5920.
245     MODD6=MOD(IDAY+NS,NSURF) 5921.
246     C**** ZERO OUT LAYER 1 WIND INCREMENTS 5922.
247     DO 60 J=1,JM 5923.
248     DO 60 I=1,IM 5924.
249     DU1(I,J)=0. 5925.
250     60 DV1(I,J)=0. 5926.
251     C**** 5927.
252     C**** OUTSIDE LOOP OVER J AND I, EXECUTED ONCE FOR EACH GRID POINT 5928.
253     C**** 5929.
254     JPR=-7
255     DO 7000 J=1,JM 5930.
256     PRNT=j.eq.8
257     PRNT=.FALSE.
258     if(PRNT)then
259     if(ns.eq.1)then
260     write(78,*) ,' '
261     write(78,*) ,'TAU=',TAU
262     endif
263     write(78,*),'NS=',ns
264     endif
265     HEMI=1. 5931.
266     IF(J.LE.JM/2) HEMI=-1. 5932.
267     FCOR=2.*OMEGA*SINP(J) 5933.
268     FMAG=FCOR*HEMI 5934.
269     ROOT2F=SQRT(2.*FMAG) 5935.
270     IF(J.EQ.1) GO TO 80 5936.
271     IF(J.EQ.JM) GO TO 90 5937.
272     WMG0=.5*(WMGE(1,J)+WMGE(1,J+1))+.001 5937.5
273     POLE=.FALSE. 5938.
274     IMAX=IM 5939.
275     GO TO 100 5940.
276     C**** CONDITIONS AT THE SOUTH POLE 5941.
277     80 POLE=.TRUE. 5942.
278     IMAX=1 5943.
279     JVPO=2 5944.
280     RAPO=2.*RAPVN(1) 5945.
281     U1=.25*(U(1,2,1)+V(IQ1,2,1)-U(IQ2,2,1)-V(IQ3,2,1)) 5946.
282     V1=.25*(V(1,2,1)-U(IQ1,2,1)-V(IQ2,2,1)+U(IQ3,2,1)) 5947.
283     WMG0=WMGE(1,2)+.001 5947.5
284     GO TO 100 5948.
285     C**** CONDITIONS AT THE NORTH POLE 5949.
286     90 POLE=.TRUE. 5950.
287     IMAX=1 5951.
288     JVPO=JM 5952.
289     RAPO=2.*RAPVS(JM) 5953.
290     U1=.25*(U(1,JM,1)-V(IQ1,JM,1)-U(IQ2,JM,1)+V(IQ3,JM,1)) 5954.
291     V1=.25*(V(1,JM,1)+U(IQ1,JM,1)-V(IQ2,JM,1)-U(IQ3,JM,1)) 5955.
292     WMG0=WMGE(1,JM)+.001 5955.5
293     C**** ZERO OUT SURFACE DIAGNOSTICS WHICH WILL BE SUMMED OVER LONGITUDE 5956.
294     100 ATRHDT=0. 5957.
295     BTRHDT=0. 5958.
296     CTRHDT=0. 5959.
297     ASHDT=0. 5960.
298     BSHDT=0. 5961.
299     CSHDT=0. 5962.
300     AEVHDT=0. 5963.
301     BEVHDT=0. 5964.
302     CEVHDT=0. 5965.
303     ATS=0. 5966.
304     BTS=0. 5967.
305     CTS=0. 5968.
306     AT2=0. 5966.
307     BT2=0. 5967.
308     CT2=0. 5968.
309     ATAUL=0.
310     ATAUF=0.
311     BTAUL=0.
312     BTAUF=0.
313     CTAUL=0.
314     CTAUF=0.
315     AWS=0.
316     BWS=0.
317     CWS=0.
318     AWMG=0.
319     BWMG=0.
320     CWMG=0.
321     ACH=0.
322     BCH=0.
323     CCH=0.
324     IM1=IM 5969.
325     #if ( defined CLM )
326     if(NS.eq.1)then
327 jscott 1.3 i=1
328     tsl4clm(i,j)=0.0
329     qs4clm(i,j)=0.0
330     ps4clm(i,j)=0.0
331     ws4clm(i,j)=0.0
332     us4clm(i,j)=0.0
333     vs4clm(i,j)=0.0
334 jscott 1.1 endif
335     #endif
336     DO 6000 I=1,IMAX 5970.
337     C**** 5971.
338     C**** DETERMINE SURFACE CONDITIONS 5972.
339     C**** 5973.
340     PLAND=FDATA(I,J,2) 5974.
341     PWATER=1.-PLAND 5975.
342     PLICE=FDATA(I,J,3)*PLAND 5976.
343     PEARTH=PLAND-PLICE 5977.
344     POICE=ODATA(I,J,2)*PWATER 5978.
345     POCEAN=PWATER-POICE 5979.
346     if(POCEAN.LE.1.E-5)then
347     POCEAN=0.
348     POICE=PWATER
349     endif
350     TTOFR=PEARTH+PLICE+POICE+POCEAN
351     if(abs(TTOFR-1).gt.1.e-3)then
352     print *,' From surface TTOFR=',TTOFR
353     print *,' J=',J,' PLAND=',PLAND,' POCEAN=',POCEAN
354     print *,'POICE=',POICE,' ODATA(I,J,2)=',ODATA(I,J,2)
355     stop
356     end if
357     SP=P(I,J) 5980.
358     PS=SP+PTOP 5981.
359     PSK=EXPBYK(PS) 5982.
360     P1=SIG(1)*SP+PTOP 5983.
361     P1K=EXPBYK(P1) 5984.
362     WSOLD=BLDATA(I,J,1) 5985.
363     USOLD=BLDATA(I,J,6) 5986.
364     VSOLD=BLDATA(I,J,7) 5987.
365     TOLD=BLDATA(I,J,8) 5988.
366     SQRTT=SQRT(TOLD) 5989.
367     GKBYFW=.1296*GRAV/(FCOR*FMAG*WSOLD+1.E-20) 5990.
368     COSWS=GKBYFW*USOLD 5991.
369     SINWS=GKBYFW*VSOLD 5992.
370     IF(POLE) GO TO 1200 5993.
371     U1=.25*(U(IM1,J,1)+U(I,J,1)+U(IM1,J+1,1)+U(I,J+1,1)) 5994.
372     V1=.25*(V(IM1,J,1)+V(I,J,1)+V(IM1,J+1,1)+V(I,J+1,1)) 5995.
373     if(J.eq.JPR.or.J.eq.-12)then
374     print *,' J=',J
375     print *,'TAU=',TAU,' NS=',NS,' ITYPE=',ITYPE
376     print *,'U(I,J,1)=',U(I,J,1),' V(I,J,1)=',V(I,J,1)
377     print *,'U(I,J+1,1)=',U(I,J+1,1),' V(I,J+1,1)=',V(I,J+1,1)
378     print *,'U(IM1,J,1)=',U(IM1,J,1),' V(IM1,J,1)=',V(IM1,J,1)
379     print *,'U(IM1,J+1,1)=',U(IM1,J+1,1),
380     & ' V(IM1,J+1,1)=',V(IM1,J+1,1)
381     endif
382     1200 TH1=T(I,J,1) 5996.
383     Q1=Q(I,J,1) 5997.
384     DTH1=0.0
385     DQQ1=0.0
386     THV1=TH1*(1.+Q1*RVX) 5998.
387     c SRHEAT=SRHR(I,J,1)*COSZ1(I,J)*SRCOR 5999.
388     c SRHDT=SRHEAT*DTSURF 6000.
389     RMBYA=100.*SP*DSIG(1)/GRAV 6001.
390     C**** ZERO OUT QUANTITIES TO BE SUMMED OVER SURFACE TYPES 6002.
391     USS=0. 6003.
392     VSS=0. 6004.
393     WSS=0. 6005.
394     TSS=0. 6006.
395     QSS=0. 6007.
396     TAUS=0. 6008.
397     SINAPS=0. 6009.
398     COSAPS=0. 6010.
399     JR=J
400     DXYPJ=DXYP(J) 6012.
401     TG1S=0. 6013.
402     QGS=0. 6014.
403     BETAS=0. 6015.
404     TRHDTS=0. 6016.
405     SHDTS=0. 6017.
406     EVHDTS=0. 6018.
407     UGS=0. 6019.
408     VGS=0. 6020.
409     WGS=0. 6021.
410     USRS=0. 6022.
411     VSRS=0. 6023.
412     RIS1S=0. 6024.
413     RIGSS=0. 6025.
414     CDMS=0. 6026.
415     CDHS=0. 6027.
416     DGSS=0. 6028.
417     EDS1S=0. 6029.
418     PPBLS=0. 6030.
419     EVAPS=0. 6031.
420     C**** 6032.
421     IF(POCEAN.LE.0.) GO TO 2200 6033.
422     C**** 6034.
423     C**** OCEAN 6035.
424     C**** 6036.
425     ITYPE=1 6037.
426     PTYPE=POCEAN 6038.
427     c formula charnoka
428     TOCEAN=BLDATA(I,J,5)
429     ROUGH=MAX(0.018*TOCEAN/GRAV,1.5E-4)
430     c ROUGH=MAX(0.035*TOCEAN/GRAV,2.5E-4) used in 008.03
431     ZGS=10. 6041.
432     ! WMGMIN=8.
433     WMGMIN=WMGMINO(J)
434     NGRNDZ=1 6043.
435     TG1=ODATA(I,J,1) 6044.
436     BETA=1. 6045.
437     ELHX=LHE 6046.
438     TRHT0=TRSURF(J,1)
439     SRHEAT=SRSURF(J,1)*COSZ1(I,J)*SRCOR
440     GO TO 3000 6047.
441     C**** 6048.
442     2200 IF(POICE.LE.0.) GO TO 2400 6049.
443     C**** 6050.
444     C**** OCEAN ICE 6051.
445     C**** 6052.
446     ITYPE=2 6053.
447     PTYPE=POICE 6054.
448     NGRNDZ=NGRND 6055.
449     SNOW=GDATA(I,J,1) 6056.
450     TG1=TGRND(I,J,2) 6057.
451     TG2=GDATA(I,J,7) 6058.
452     ACE2=ODATA(I,J,3) 6059.
453     Z2=ACE2/RHOI 6060.
454     Z2BY4L=Z2/(4.*ALAMI) 6061.
455     if (SNOW.gt.10.)then
456     RHOS0=ROSNOW(SNOW)
457     else
458     RHOS0=275.
459     endif
460     RHOS=COEFSN*RHOS0
461     ALAMS=ALSNOW(RHOS0)
462     BYRSL=1./(RHOS*ALAMS)
463     c Z1BY6L=(Z1IBYL+SNOW*BYRSL)*.1666667 6062.
464     c CDTERM=1.5*TG2-.5*TFO 6063.
465     CDTERM=TG2
466     c CDENOM=1./(2.*Z1BY6L+Z2BY4L) 6064.
467     Z1BY2L=(Z1IBYL+SNOW*BYRSL)*0.5
468     CDENOM=1./(Z1BY2L+2.*Z2BY4L)
469     ROUGH=10**(log10(10.)-4.37)
470     ZGS=10. 6067.
471     ! WMGMIN=8.
472     WMGMIN=WMGMINO(J)
473     HC1=HC1I+SNOW*SHI 6069.
474     BETA=1. 6070.
475     ELHX=LHS 6071.
476     TRHT0=TRSURF(J,3)
477     SRHEAT=SRSURF(J,3)*COSZ1(I,J)*SRCOR
478     GO TO 3000 6072.
479     C**** 6073.
480     2400 IF(PLAND.LE.0.) GO TO 5000 6074.
481     NGRNDZ=NGRND 6075.
482     ROUGH=ROUGHL(I,J) 6076.
483     ZGS=30. 6078.
484     WMGMIN=WMGMINL
485     TRHT0=TRSURF(J,2)
486     SRHEAT=SRSURF(J,2)*COSZ1(I,J)*SRCOR
487     IF(PLICE.LE.0.) GO TO 2600 6080.
488     C**** 6081.
489     C**** LAND ICE 6082.
490     C**** 6083.
491     ITYPE=3 6084.
492     PTYPE=PLICE 6085.
493     SNOW=GDATA(I,J,12) 6086.
494     TG1=TGRND(I,J,3) 6087.
495     TG2=GDATA(I,J,14) 6088.
496     if (SNOW.gt.10.)then
497     RHOS0=ROSNOW(SNOW)
498     else
499     RHOS0=275.
500     endif
501     RHOS=COEFSN*RHOS0
502     ALAMS=ALSNOW(RHOS0)
503     BYRSL=1./(RHOS*ALAMS)
504     c Z1BY6L=(Z1IBYL+SNOW*BYRSL)*.1666667 6089.
505     CDTERM=TG2 6090.
506     c CDENOM=1./(2.*Z1BY6L+Z2LI3L) 6091.
507     Z1BY2L=(Z1IBYL+SNOW*BYRSL)*0.5
508     CDENOM=1./(Z1BY2L+3.*Z2LI3L/2.)
509     HC1=HC1I+SNOW*SHI 6092.
510     BETA=1. 6093.
511     ELHX=LHS 6094.
512     GO TO 3000 6095.
513     C**** 6096.
514     2600 IF(PEARTH.LE.0.) GO TO 5000 6097.
515     C**** 6098.
516     C**** EARTH 6099.
517     C**** 6100.
518     ITYPE=4 6101.
519     PTYPE=PEARTH 6102.
520     SNOW=GDATA(I,J,2) 6103.
521     TG1=TGRND(I,J,4) 6104.
522     WTR1=GDATA(I,J,5) 6105.
523     ACE1=GDATA(I,J,6) 6106.
524     TG2=GDATA(I,J,8) 6107.
525     WTR2=GDATA(I,J,9) 6108.
526     ACE2=GDATA(I,J,10) 6109.
527     WFC1=VDATA(I,J,9) 6110.
528     WFC2=VDATA(I,J,10) 6111.
529     WTR1DRY=0.025*WFC1
530     HC1=HC1DE+WTR1*SHW+(ACE1+SNOW)*SHI 6112.
531     ALAM1D=2.+.5*(1.+2.*WTR1/WFC1) 6113.
532     ALAM2D=4. 6114.
533     RMULCH=1. 6115.
534     IF((SINP(J).GT..5).AND.(JDAY-91)*(273-JDAY).LT.0) RMULCH=.25 6116.
535     IF((SINP(J).LT.-.5).AND.(JDAY-91)*(273-JDAY).GE.0) RMULCH=.25 6117.
536     ALAM1V=RMULCH*(.4185+1.2555*WTR1/WFC1+ALAMI*ACE1/(Z1E*RHOI)) 6118.
537     ALAM3V=.8370 6119.
538     IF(TG2.LT.0.) ALAM3V=.4185+ALAMI*.15 6120.
539     ALAM2V=.125*(.4185+1.2555*WTR2/WFC2+ALAMI*ACE2/(5.*Z1E*RHOI)) 6121.
540     * +.875*ALAM3V 6122.
541     ALAM1E=VDATA(I,J,1)*ALAM1D+(1.-VDATA(I,J,1))*ALAM1V 6123.
542     ALAM2E=VDATA(I,J,1)*ALAM2D+(1.-VDATA(I,J,1))*ALAM2V 6124.
543     if (SNOW.gt.10.)then
544     RHOS0=ROSNOW(SNOW)
545     else
546     RHOS0=275.
547     endif
548     RHOS=COEFSN*RHOS0
549     ALAMS=ALSNOW(RHOS0)
550     BYRSL=1./(RHOS*ALAMS)
551     c Z1BY6L=(Z1E/ALAM1E+SNOW*BYRSL)*.1666667 6125.
552     Z1BY2L=(Z1E/ALAM1E+SNOW*BYRSL)*0.5
553     CDTERM=TG2 6126.
554     c CDENOM=1./(2.*Z1BY6L+Z2E/(3.*ALAM2E)) 6127.
555     CDENOM=1./(Z1BY2L+Z2E/(2.*ALAM2E))
556     BETA=1. 6128.
557     ELHX=LHS 6129.
558     IF(SNOW.GT.0.) GO TO 3000 6130.
559     BETA=(WTR1+ACE1)/WFC1 6131.
560     BETA=max(((WTR1+ACE1-WTR1DRY)/WFC1),0.0)
561     PFROZN=ACE1/(WTR1+ACE1+1.E-20) 6132.
562     ELHX=LHE+LHM*PFROZN 6133.
563     HC2E=HC2DE+WTR2*SHW+ACE2*SHI
564     C**** 6134.
565     C**** BOUNDARY LAYER INTERACTION 6135.
566     C**** 6136.
567     3000 continue
568     SRHDT=SRHEAT*DTSURF
569     TKV=THV1*PSK 6137.
570     ZS1=ZS1CO*TKV*SP/PS 6138.
571     P1=SIG(1)*SP+PTOP 6139.
572     DTGRND=DTSURF/NGRNDZ 6143.
573     SHDT=0. 6144.
574     EVHDT=0. 6145.
575     TRHDT=0. 6146.
576     F1DT=0. 6147.
577     C**** LOOP OVER GROUND TIME STEPS 6148.
578     DO 3600 NG=1,NGRNDZ 6149.
579     TG=TG1+TF 6150.
580     QG=QSAT(TG,PS,ELHX) 6151.
581     TGV=TG*(1.+QG*RVX) 6152.
582     W1=SQRT(U1*U1+V1*V1)
583     WS0=W1
584     c WS=SQRT(W1*W1+0.8*WMG)
585     ! WMG=WMG0+WMGMIN
586     ! 07/17/2006
587     if(ITYPE.le.2)then
588     WMG=WMG0+WMGMIN
589     else
590     WMG=max(WMG0,WMGMIN)
591     endif
592     ! 07/17/2006
593     WS=SQRT((0.75*W1)**2+WMG)
594     if(J.eq.JPR)then
595     print *,' '
596     print *,'TAU=',TAU,' NS=',NS,' ITYPE=',ITYPE
597     print *,'TG=',TG,' QG=',QG
598     print *,'RVX=',RVX,' TG1=',TG1
599     endif
600 jscott 1.3 ! if(TAU.gt.3623.0.and.TAU.lt.3744.0) then
601     ! if(ITYPE.eq.4)then
602     ! if(J.eq.29)then
603     ! write (329) ,TAU,TOFDAY,TG1
604     ! endif
605     ! endif
606     ! endif
607 jscott 1.1
608     #if ( defined CPL_OCEANCO2 || defined OCEAN_3D)
609     if(ITYPE.eq.1)then
610     NWMGEA(J)=NWMGEA(J)+1
611     WSAV(J)=WSAV(J)+WS
612     end if
613     #endif
614    
615     WG=WS
616     THS=TH1
617     QS=Q1
618     TSV=THS*PSK
619     Z0=ROUGH
620     ROUGH=log10(ZGS/ROUGH)
621     CDN=.0231/(ROUGH*ROUGH)
622     c if(ITYPE.eq.1)then
623     c CDN=.00075+.000067*WSOLD
624     c ROUGH=7.126-1.068*LOG(WSOLD+1.E-12)
625     c endif
626     LR=ROUGH*2.-.5
627     IF(LR.GT.20) LR=20
628     IF(LR.LT.1) LR=1
629     RIGS=ZGS*GRAV*(TSV-TGV)/(TGV*WS*WS)
630     SINAP=0.
631     COSAP=1.
632     IF(RIGS.LE.0) THEN
633     C surface layer has unstable stratification
634     CIA=TWOPI*0.0625/(1.+WS*CIAX)
635     DM=SQRT((1.-AROUGH(LR)*RIGS)*(1.-BROUGH(LR)*RIGS)/
636     * (1.-CROUGH(LR)*RIGS))
637     DH=1.35*SQRT((1.-DROUGH(LR)*RIGS)/(1.-EROUGH(LR)*RIGS))
638     ELSE
639     C surface layer has stable stratification
640     CIA=TWOPI*(0.09375-0.03125/(1.+4*RIGS**2))/(1.+WS*CIAX)
641     DM=1./(1.+(11.238+89.9*RIGS)*RIGS)
642     DH=1.35/(1.+1.93*RIGS)
643     END IF
644     CDH=CDN*DM*DH
645     if(J.eq.JPR)then
646     print *,'TAU=',TAU,' NS=',NS,' ITYPE=',ITYPE
647     print *,'WS=',WS,' ZGS=',ZGS
648     print *,'DM=',DM,' DH=',DH
649     print *,'RIGS=',RIGS,' TGV=',TGV
650     endif
651     USR=COS(CIA)
652     VSR=SIN(CIA)*HEMI
653     UG=U1
654     VG=V1
655     US=(USR*UG-VSR*VG)
656     VS=(VSR*UG+USR*VG)
657     RCDHWS=CDH*WS*100.*PS/(RGAS*TSV)
658     if(J.eq.JPR)then
659     c print *,' '
660     print *,'TAU=',TAU,' NS=',NS,' ITYPE=',ITYPE
661     print *,'CDH=',CDH,' RGAS=',RGAS
662     print *,'PS=',PS,' TSV=',TSV
663     print *,'WS=',WS,' RCDHWS=',RCDHWS
664     endif
665     TS=TSV/(1.+QS*RVX) 6467.
666     QSATS=QSAT(TS,PS,ELHX) 6468.
667     c dLQS/dTs
668     DLQSDTS=DLQSDT(TS,ELHX)
669     c dLQS/dTs
670     IF(QS.LE.QSATS) GO TO 3500 6469.
671     DQSSDT=QSATS*ELHX/(RVAP*TS*TS) 6470.
672     X=(QS-QSATS)/(DQSSDT+(SHA/ELHX)) 6471.
673     TS=TS+X 6472.
674     QS=QS+X*(SHA/ELHX) 6473.
675     C**** CALCULATE RHOS*CDM*WS AND RHOS*CDH*WS 6474.
676     3500 CDM=CDN*DM 6475.
677     RCDMWS=CDM*WS*100.*PS/(RGAS*TS) 6476.
678     C**** CALCULATE FLUXES OF SENSIBLE HEAT, LATENT HEAT, THERMAL 6478.
679     C**** RADIATION, AND CONDUCTION HEAT (WATTS/M**2) 6479.
680     SHEAT=SHA*RCDHWS*(TS-TG) 6480.
681     BETAUP=BETA 6481.
682     IF(QS.GT.QG) BETAUP=1. 6482.
683     EVHEAT=(LHE+TG1*SHV)*BETAUP*RCDHWS*(QS-QG) 6483.
684     c TRHEAT=TRHR(I,J,1)-STBO*(TG*TG)*(TG*TG) 6484.
685     TRHEAT=TRHT0-STBO*(TG*TG)*(TG*TG)
686     #if ( defined CLM )
687     if(NS.eq.1)then
688     if(ITYPE.EQ.4.or.ITYPE.EQ.3)then
689 jscott 1.3 tsl4clm(i,j)=tsl4clm(i,j)+TS*PTYPE/PLAND
690     qs4clm(i,j)=qs4clm(i,j)+QS*PTYPE/PLAND
691     ps4clm(i,j)=ps4clm(i,j)+PS*PTYPE/PLAND
692     ws4clm(i,j)=ws4clm(i,j)+WS*PTYPE/PLAND
693     us4clm(i,j)=us4clm(i,j)+US*PTYPE/PLAND
694     vs4clm(i,j)=vs4clm(i,j)+VS*PTYPE/PLAND
695 jscott 1.1 endif
696     endif
697     #endif
698     if(J.eq.JPR)then
699     c print *,' '
700     print *,'TAU=',TAU,' NS=',NS,' ITYPE=',ITYPE
701     print *,'TRHT0=',TRHT0,' STBO=',STBO
702     print *,'TG=',TG,' TS=',TS
703     print *,'TRHEAT=',TRHEAT
704     print *,'SHA=',SHA,' RCDHWS=',RCDHWS
705     print *,'SHEAT=',SHEAT
706     endif
707     TG1OLD=TG1
708     SHEATOLD=SHEAT
709     #if ( defined OCEAN_3D )
710     IF(ITYPE.EQ.1 .or. ITYPE.EQ.2) GO TO 3620
711     #else
712     IF(ITYPE.EQ.1) GO TO 3620 6485.
713     #endif
714     C**** CALCULATE FLUXES USING IMPLICIT TIME STEP FOR NON-OCEAN POINTS 6486.
715     F0=SRHEAT+TRHEAT+SHEAT+EVHEAT 6487.
716     c F1=(TG1-CDTERM-F0*Z1BY6L)*CDENOM 6488.
717     F1=(TG1-CDTERM)*CDENOM
718     DSHDTG=-RCDHWS*SHA
719     DQGDTG=QG*ELHX/(RVAP*TG*TG) 6490.
720     DEVDTG=-RCDHWS*LHE*BETAUP*DQGDTG
721     DTRDTG=-4.*STBO*TG*TG*TG 6492.
722     DF0DTG=DSHDTG+DEVDTG+DTRDTG 6493.
723     c DFDTG=DF0DTG-(1.-DF0DTG*Z1BY6L)*CDENOM 6493.5
724     DFDTG=DF0DTG-CDENOM
725     c DF1DTG=(1.-DF0DTG*Z1BY6L)*CDENOM
726     DF1DTG=CDENOM
727     DTG=(F0-F1)*DTGRND/(HC1-DTGRND*DFDTG) 6494.
728     SHDT=SHDT+DTGRND*(SHEAT+DTG*DSHDTG) 6495.
729     EVHDT=EVHDT+DTGRND*(EVHEAT+DTG*DEVDTG) 6496.
730     TRHDT=TRHDT+DTGRND*(TRHEAT+DTG*DTRDTG) 6497.
731     TG1=TG1+DTG
732     c F1DT=F1DT+DTGRND*(TG1-CDTERM-(F0+DTG*DF0DTG)*Z1BY6L)*CDENOM 6498.
733     F1DT=F1DT+DTGRND*(TG1-CDTERM)*CDENOM
734     DU1(I,J)=DU1(I,J)+PTYPE*DTGRND*RCDMWS*US*COEFS/SP 6499.
735     DV1(I,J)=DV1(I,J)+PTYPE*DTGRND*RCDMWS*VS*COEFS/SP 6500.
736     c TG1=TG1+DTG 6501.
737     3600 CONTINUE 6502.
738     GO TO 3700 6503.
739     C**** CALCULATE FLUXES USING EXPLICIT TIME STEP FOR OCEAN POINTS 6504.
740     3620 SHDT=DTSURF*SHEAT 6505.
741     EVHDT=DTSURF*EVHEAT 6506.
742     TRHDT=DTSURF*TRHEAT 6507.
743     DU1(I,J)=DU1(I,J)+PTYPE*DTSURF*RCDMWS*US*COEFS/SP 6508.
744     DV1(I,J)=DV1(I,J)+PTYPE*DTSURF*RCDMWS*VS*COEFS/SP 6509.
745     3700 CONTINUE
746     EPS=1.D-8
747     c print *,'FROM SURFACE NS=',NS
748     c print *,'J=',J,' ITYPE=',ITYPE
749     c print *,RCDMWS,WS
750     WWS=max(W1,1.D-4)
751     c RO=SP*100/(RGAS*TG)
752     c print *,'RO=',RO
753     c USTAR=SQRT(RCDMWS*WS/RO)
754     c TSTAR=SHEATOLD/(0.35*1007.*RO*USTAR)
755     c ALPHAH=DH
756     c TT2M=TG+TSTAR/ALPHAH*LOG(2.0/Z0)
757     c TT2M=TG+TSTAR/ALPHAH*LOG(ZGS/Z0)
758     c print *,'RIGS=',RIGS,' Z0=',Z0
759     c print *,'CDN=',CDN
760     c print *,'H=',SHDT/DTSURF,' TGM=',RCDMWS*WS
761     c print *,' SHEATOLD=',SHEATOLD
762     c print *,' USTAR=',USTAR,' TSTAR=',TSTAR
763     c print *,' ALPHAH=',ALPHAH,' TT2M=',TT2M
764     c print *,' TT2M=',TT2M
765     ZTEM=ZGS
766     ZTEM=2.0
767     c print *,'ZTEM=',ZTEM
768     CALL TZM(T2M,TH2M,ZTEM,Z0,ZGS,SP,TG,TS,RIGS,WS,
769     & -SHEATOLD,RCDMWS*WS,LR,EPS)
770     c print *,'FROM SURFACE'
771     c print *,'TS=',TS,' TG=',TG
772     c print *,' T2M=',T2M,' TH2M=',TH2M
773     F0DT=CORSR*SRHDT+TRHDT+SHDT+EVHDT 6510.
774     if(J.eq.JPR)then
775     print *,'TAU=',TAU,' NS=',NS,' ITYPE=',ITYPE
776     print *,'DTSURF=',DTSURF,' CORSR=',CORSR
777     print *,'SRHDT=',SRHDT,' TRHDT=',TRHDT
778     print *,'SHDT=',SHDT,' EVHDT=',EVHDT
779     print *,'F0DT=',F0DT
780     print *,'US=',US,' VS=',VS
781     print *,'COEFS=',COEFS,' SP=',SP
782     endif
783     c print *,'From surface ',TAU,CORSR,SRHDT,TRHDT,SHDT,EVHDT
784     C**** CALCULATE EVAPORATION 6511.
785     CCC DQ1=EVHDT/((LHE+TG1*SHV)*RMBYA) 6512.
786     DQ1=EVHDT/(ELHX*RMBYA)
787     IF(DQ1*PTYPE.LE.Q1) GO TO 3720 6513.
788     DQ1=Q1/PTYPE 6514.
789     CCC EVHDT=DQ1*(LHE+TG1*SHV)*RMBYA 6515.
790     EVHDT=DQ1*ELHX*RMBYA
791     3720 EVAP=-DQ1*RMBYA 6516.
792     C**** ACCUMULATE SURFACE FLUXES AND PROGNOSTIC AND DIAGNOSTIC QUANTITIES6517.
793     E0(I,J,ITYPE)=E0(I,J,ITYPE)+F0DT 6518.
794     E1(I,J,ITYPE)=E1(I,J,ITYPE)+F1DT 6519.
795     EVAPOR(I,J,ITYPE)=EVAPOR(I,J,ITYPE)+EVAP 6520.
796     if(PRNT)then
797     c write(78,*) ,' '
798     c write(78,*) ,'TAU=',TAU
799     write(78,*) ,'J=',j,' ITYPE=',ITYPE,' PTYPE=',PTYPE
800     write(78,*) ,'TS=',TS,' TG=',TG,' QS=',QS
801     write(78,*) ,'TG1=',TG1,' TG1OLD=',TG1OLD
802     write(78,*) ,'TG2=',TG2
803     write(78,*) ,'SHEAT=',SHEAT,' EVHEAT=',EVHEAT
804     write(78,*) ,'TRHEAT=',TRHEAT,' SRHEAT=',SRHEAT
805     write(78,*) ,'EVAP mm/day=',24.*3600.*EVAP/DTSURF
806     write(78,*) ,'EVAP=',EVAP,
807     & ' F0DT=',F0DT/DTSURF,' F1DT=',F1DT/DTSURF
808     endif
809     #if ( defined OCEAN_3D || defined ML_2D )
810     C For ocean model
811     c if(NS.eq.2)then
812     #if ( defined ML_2D )
813     if(ITYPE.eq.1)then
814     #endif
815     C DNetHeat by DTG
816     DSHDTG=-RCDHWS*SHA
817     DQGDTG=QG*ELHX/(RVAP*TG*TG)
818     DEVDTG=-RCDHWS*LHE*BETAUP*DQGDTG
819     DTRDTG=-4.*STBO*TG*TG*TG
820     DF0DTG=DSHDTG+DEVDTG+DTRDTG
821     if(EVHEAT.lt.0.0)then
822     DEVDTGEQ=EVHEAT*DLQSDTS
823     else
824     DEVDTGEQ=0.0
825     endif
826     C DNetHeat by DTG
827     #if ( defined OCEAN_3D )
828     if(ITYPE.eq.1)then
829     #endif
830     dhfodtg(j)=dhfodtg(j)+DF0DTG
831     devodtg(j)=devodtg(j)-DEVDTG/LHE
832     dhfodtgeq(j)=dhfodtgeq(j)+DEVDTGEQ
833     devodtgeq(j)=devodtgeq(j)-DEVDTGEQ/LHE
834     evao(j)=evao(j)+EVAP
835     hfluxo(j)=hfluxo(j)+F0DT
836     naveo(j)=naveo(j)+1
837     endif
838     if(ITYPE.eq.2)then
839     evai(j)=evai(j)+EVAP
840     hfluxi(j)=hfluxi(j)+F0DT
841     dhfidtg(j)=dhfidtg(j)+DF0DTG
842     devidtg(j)=devidtg(j)-DEVDTG/LHE
843     dhfidtgeq(j)=dhfidtgeq(j)+DEVDTGEQ
844     devidtgeq(j)=devidtgeq(j)-DEVDTGEQ/LHE
845     c tairi(j)=tairi(j)+TS
846     navei(j)=navei(j)+1
847     endif
848     c endif ! NS
849     tauu(j)=tauu(j)+RCDMWS*US*PTYPE
850     tauv(j)=tauv(j)+RCDMWS*VS*PTYPE
851     C For ocean model
852     #endif
853     TGRND(I,J,ITYPE)=TG1 6521.
854     TSSFC(I,J,ITYPE)=TS 6521.5
855    
856     c TH1=TH1-SHDT*PTYPE/(SHA*RMBYA*P1K) 6522.
857     c Q1=Q1-DQ1*PTYPE 6523.
858    
859     DTH1=DTH1-SHDT*PTYPE/(SHA*RMBYA*P1K)
860     DQQ1=DQQ1-DQ1*PTYPE
861    
862     USS=USS+US*PTYPE 6524.
863     VSS=VSS+VS*PTYPE 6525.
864     WSS=WSS+WS*PTYPE 6526.
865     TSS=TSS+TS*PTYPE 6527.
866     QSS=QSS+QS*PTYPE 6528.
867     TAUS=TAUS+CDM*WS*W1*PTYPE 6529.
868     SINAPS=SINAPS+SINAP*PTYPE 6530.
869     COSAPS=COSAPS+COSAP*PTYPE 6531.
870     TG1S=TG1S+TG1*PTYPE 6532.
871     QGS=QGS+QG*PTYPE 6533.
872     BETAS=BETAS+BETA*PTYPE 6534.
873     TRHDTS=TRHDTS+TRHDT*PTYPE 6535.
874     SHDTS=SHDTS+SHDT*PTYPE 6536.
875     EVHDTS=EVHDTS+EVHDT*PTYPE 6537.
876     UGS=UGS+UG*PTYPE 6538.
877     VGS=VGS+VG*PTYPE 6539.
878     WGS=WGS+WG*PTYPE 6540.
879     USRS=USRS+USR*PTYPE 6541.
880     VSRS=VSRS+VSR*PTYPE 6542.
881     c RIS1S=RIS1S+RIS1*PTYPE 6543.
882     RIGSS=RIGSS+RIGS*PTYPE 6544.
883     CDMS=CDMS+CDM*PTYPE 6545.
884     CDHS=CDHS+CDH*PTYPE 6546.
885     c DGSS=DGSS+DGS*PTYPE 6547.
886     c EDS1S=EDS1S+EDS1*PTYPE 6548.
887     c PPBLS=PPBLS+PPBL*PTYPE 6549.
888     EVAPS=EVAPS+EVAP*PTYPE 6550.
889     GO TO (4000,4100,4400,4600),ITYPE 6551.
890     C**** 6552.
891     C**** OCEAN 6553.
892     C**** 6554.
893     4000 ASHDT=ASHDT+SHDT*POCEAN 6555.
894     AEVHDT=AEVHDT+EVHDT*POCEAN 6556.
895     ATRHDT=ATRHDT+TRHDT*POCEAN 6557.
896     ATS=ATS+(TS-TF)*POCEAN 6558.
897     AT2=AT2+(TH2M-TF)*POCEAN
898     BLDATA(I,J,5)=CDM*WS*W1
899     ATAUL=ATAUL+RCDMWS*US*POCEAN
900     ATAUF=ATAUF+RCDMWS*VS*POCEAN
901     AWS=AWS+WS*POCEAN
902     AWMG=AWMG+SQRT(WMG)*POCEAN
903     ACH=ACH+RCDHWS*POCEAN
904     GO TO 2200 6559.
905     C**** 6560.
906     C**** OCEAN ICE 6561.
907     C**** 6562.
908     4100 CSHDT=CSHDT+SHDT*POICE 6563.
909     CEVHDT=CEVHDT+EVHDT*POICE 6564.
910     CTRHDT=CTRHDT+TRHDT*POICE 6565.
911     CTS=CTS+(TS-TF)*POICE 6566.
912     CT2=CT2+(TH2M-TF)*POICE 6566.
913     CTAUL=CTAUL+RCDMWS*US*POICE
914     CTAUF=CTAUF+RCDMWS*VS*POICE
915     CWS=CWS+WS*POICE
916     CWMG=CWMG+SQRT(WMG)*POICE
917     CCH=CCH+RCDHWS*POICE
918     GO TO 2400 6567.
919     C**** 6568.
920     C**** LAND ICE 6569.
921     C**** 6570.
922     4400 BSHDT=BSHDT+SHDT*PLICE 6571.
923     BEVHDT=BEVHDT+EVHDT*PLICE 6572.
924     BTRHDT=BTRHDT+TRHDT*PLICE 6573.
925     BTS=BTS+(TS-TF)*PLICE 6574.
926     BT2=BT2+(TH2M-TF)*PLICE
927     BTAUL=BTAUL+RCDMWS*US*PLICE
928     BTAUF=BTAUF+RCDMWS*VS*PLICE
929     BWS=BWS+WS*PLICE
930     BWMG=BWMG+SQRT(WMG)*PLICE
931     BCH=BCH+RCDHWS*PLICE
932     GO TO 2600 6575.
933     C**** 6576.
934     C**** EARTH 6577.
935     C**** 6578.
936     4600 BSHDT=BSHDT+SHDT*PEARTH 6579.
937     BEVHDT=BEVHDT+EVHDT*PEARTH 6580.
938     BTRHDT=BTRHDT+TRHDT*PEARTH 6581.
939     BTS=BTS+(TS-TF)*PEARTH 6582.
940     BT2=BT2+(TH2M-TF)*PEARTH
941     BTAUL=BTAUL+RCDMWS*US*PEARTH
942     BTAUF=BTAUF+RCDMWS*VS*PEARTH
943     BWS=BWS+WS*PEARTH
944     BWMG=BWMG+SQRT(WMG)*PEARTH
945     BCH=BCH+RCDHWS*PEARTH
946     C**** NON-OCEAN POINTS WHICH ARE NOT MELTING OR FREEZING WATER USE 6583.
947     C**** IMPLICIT TIME STEPS 6584.
948     C**** 6585.
949     C**** UPDATE SURFACE AND FIRST LAYER QUANTITIES 6586.
950     C**** 6587.
951     5000 CONTINUE
952     T(I,J,1)=TH1 6588.
953     & +DTH1
954     Q(I,J,1)=Q1 6589.
955     & +DQQ1
956     BLDATA(I,J,1)=WSS 6590.
957     BLDATA(I,J,2)=TSS 6591.
958     BLDATA(I,J,3)=QSS 6592.
959     BLDATA(I,J,6)=USS 6593.
960     BLDATA(I,J,7)=VSS 6594.
961     BLDATA(I,J,8)=TAUS 6595.
962     c print *,j,T(I,J,1),Q(I,J,1)
963     c print *,(TGRND(I,J,k),k=1,4)
964     c print *,(EVAPOR(I,J,k),k=1,4)
965     c print *,(E0(I,J,k),k=1,4)
966     c print *,(E1(I,J,k),k=1,4)
967     c print *,j,DU1(1,j),DV1(1,j)
968     C**** 6596.
969     C**** ACCUMULATE DIAGNOSTICS 6597.
970     C**** 6598.
971     C**** QUANTITIES ACCUMULATED FOR REGIONS IN DIAG1 6599.
972     IF(JR.EQ.JM) GO TO 5700 6600.
973     DJ(JR,9)=DJ(JR,9)+TRHDTS*DXYPJ 6601.
974     DJ(JR,13)=DJ(JR,13)+SHDTS*DXYPJ 6602.
975     DJ(JR,14)=DJ(JR,14)+EVHDTS*DXYPJ 6603.
976     DJ(JR,19)=DJ(JR,19)+EVAPS*DXYPJ 6604.
977     IF(MODDSF.NE.0) GO TO 5700 6605.
978     DJ(JR,23)=DJ(JR,23)+(TSS-TF)*DXYPJ 6606.
979     5700 CONTINUE
980     6000 IM1=I 6662.
981     C**** QUANTITIES ACCUMULATED FOR SURFACE TYPE TABLES IN DIAG1 6663.
982     AJ(J,9)=AJ(J,9)+ATRHDT 6664.
983     BJ(J,9)=BJ(J,9)+BTRHDT 6665.
984     CJ(J,9)=CJ(J,9)+CTRHDT 6666.
985     AJ(J,13)=AJ(J,13)+ASHDT 6667.
986     BJ(J,13)=BJ(J,13)+BSHDT 6668.
987     CJ(J,13)=CJ(J,13)+CSHDT 6669.
988     AJ(J,14)=AJ(J,14)+AEVHDT 6670.
989     BJ(J,14)=BJ(J,14)+BEVHDT 6671.
990     CJ(J,14)=CJ(J,14)+CEVHDT 6672.
991     AJ(J,32)=AJ(J,32)+ATAUL
992     BJ(J,32)=BJ(J,32)+BTAUL
993     CJ(J,32)=CJ(J,32)+CTAUL
994     AJ(J,33)=AJ(J,33)+ATAUF
995     BJ(J,33)=BJ(J,33)+BTAUF
996     CJ(J,33)=CJ(J,33)+CTAUF
997     AJ(J,37)=AJ(J,37)+AWS
998     BJ(J,37)=BJ(J,37)+BWS
999     CJ(J,37)=CJ(J,37)+CWS
1000     AJ(J,28)=AJ(J,28)+AWMG
1001     BJ(J,28)=BJ(J,28)+BWMG
1002     CJ(J,28)=CJ(J,28)+CWMG
1003     AJ(J,38)=AJ(J,38)+ATAUL/NSURF
1004     BJ(J,38)=BJ(J,38)+BTAUL/NSURF
1005     CJ(J,38)=CJ(J,38)+CTAUL/NSURF
1006     IF(MODDSF.NE.0) GO TO 7000 6673.
1007     AJ(J,23)=AJ(J,23)+ATS 6674.
1008     BJ(J,23)=BJ(J,23)+BTS 6675.
1009     CJ(J,23)=CJ(J,23)+CTS 6676.
1010     AJ(J,26)=AJ(J,26)+AT2 6674.
1011     BJ(J,26)=BJ(J,26)+BT2 6675.
1012     CJ(J,26)=CJ(J,26)+CT2 6676.
1013     c print *,j,'ATS=',ATS,' AT2=',AT2
1014     c print *,'BLDATA'
1015     c print *,(BLDATA(1,j,k),k=1,3)
1016     c print *,(BLDATA(1,j,k),k=6,8)
1017    
1018     7000 CONTINUE 6677.
1019     C**** 6678.
1020     C**** ADD IN SURFACE FRICTION TO FIRST LAYER WIND 6679.
1021     C**** 6680.
1022     DO 7600 I=1,IM 6681.
1023     U(I,2,1)=U(I,2,1)-2.*(DU1(1,1)*COSI(I)-DV1(1,1)*SINI(I))*RAPVN(1) 6682.
1024     V(I,2,1)=V(I,2,1)-2.*(DV1(1,1)*COSI(I)+DU1(1,1)*SINI(I))*RAPVN(1) 6683.
1025     U(I,JM,1)=U(I,JM,1) 6684.
1026     * -2.*(DU1(1,JM)*COSI(I)+DV1(1,JM)*SINI(I))*RAPVS(JM) 6685.
1027     7600 V(I,JM,1)=V(I,JM,1) 6686.
1028     * -2.*(DV1(1,JM)*COSI(I)-DU1(1,JM)*SINI(I))*RAPVS(JM) 6687.
1029     DO 7700 J=2,JMM1 6688.
1030     I=IM 6689.
1031     DO 7700 IP1=1,IM 6690.
1032     if(J.eq.JPR.or.J.eq.-12)then
1033     print *,' J=',J,' before'
1034     print *,'U(I,J,1)=',U(I,J,1),' V(I,J,1)=',V(I,J,1)
1035     print *,'U(I,J+1,1)=',U(I,J+1,1),' V(I,J+1,1)=',V(I,J+1,1)
1036     print *,'DU1(I,J)=',DU1(I,J),' DU1(IP1,J)=',DU1(IP1,J)
1037     endif
1038     U(I,J,1)=U(I,J,1)-(DU1(I,J)+DU1(IP1,J))*RAPVS(J) 6691.
1039     V(I,J,1)=V(I,J,1)-(DV1(I,J)+DV1(IP1,J))*RAPVS(J) 6692.
1040     U(I,J+1,1)=U(I,J+1,1)-(DU1(I,J)+DU1(IP1,J))*RAPVN(J) 6693.
1041     V(I,J+1,1)=V(I,J+1,1)-(DV1(I,J)+DV1(IP1,J))*RAPVN(J) 6694.
1042     if(J.eq.JPR.or.J.eq.-12)then
1043     print *,' J=',J,' after'
1044     print *,'U(I,J,1)=',U(I,J,1),' V(I,J,1)=',V(I,J,1)
1045     print *,'U(I,J+1,1)=',U(I,J+1,1),' V(I,J+1,1)=',V(I,J+1,1)
1046     print *,'DU1(I,J)=',DU1(I,J),' DU1(IP1,J)=',DU1(IP1,J)
1047     endif
1048     7700 I=IP1 6695.
1049     c print *,'U V'
1050     c do j=1,jm
1051     c print *,j,U(I,J,1),v(I,J,1)
1052     c enddo
1053     C**** 6696.
1054     C**** DRY CONVECTION ORIGINATING FROM THE FIRST LAYER 6697.
1055     C**** 6698.
1056     C**** LOAD U,V INTO UT,VT. UT,VT WILL BE FIXED DURING DRY CONVECTION 6699.
1057     C**** WHILE U,V WILL BE UPDATED. 6700.
1058     DO 8050 L=1,LM 6701.
1059     DO 8050 J=2,JM 6702.
1060     DO 8050 I=1,IM 6703.
1061     UT(I,J,L)=U(I,J,L) 6704.
1062     8050 VT(I,J,L)=V(I,J,L) 6705.
1063     C**** OUTSIDE LOOPS OVER J AND I 6706.
1064     DO 8500 J=1,JM 6707.
1065     POLE=.FALSE. 6708.
1066     IF(J.EQ.1.OR.J.EQ.JM) POLE=.TRUE. 6709.
1067     IMAX=IM 6710.
1068     IF(POLE) IMAX=IM 6711.
1069     DO 8120 K=1,2 6712.
1070     RA(K)=RAPVS(J) 6713.
1071     8120 RA(K+2)=RAPVN(J) 6714.
1072     IM1=IM 6715.
1073     DO 8500 I=1,IMAX 6716.
1074     BLDATA(I,J,4)=1. 6717.
1075     IF(T(I,J,1)*(1.+Q(I,J,1)*RVX).LE. 6718.
1076     * T(I,J,2)*(1.+Q(I,J,2)*RVX)) GO TO 8500 6719.
1077     C**** MIX HEAT AND MOISTURE THROUGHOUT THE BOUNDARY LAYER 6720.
1078     PKMS=PK(I,J,1)*DSIG(1)+PK(I,J,2)*DSIG(2) 6721.
1079     THPKMS=T(I,J,1)*(PK(I,J,1)*DSIG(1))+T(I,J,2)*(PK(I,J,2)*DSIG(2)) 6722.
1080     QMS=Q(I,J,1)*DSIG(1)+Q(I,J,2)*DSIG(2) 6723.
1081     TVMS=T(I,J,1)*(1.+Q(I,J,1)*RVX)*(PK(I,J,1)*DSIG(1)) 6724.
1082     * +T(I,J,2)*(1.+Q(I,J,2)*RVX)*(PK(I,J,2)*DSIG(2)) 6725.
1083     THETA=TVMS/PKMS 6726.
1084    
1085     #if ( defined CPL_CHEM )
1086     !
1087     ! --- 03/23/95
1088     !
1089     cfc11ms = cfc11(i,j,1)*dsig(1) + cfc11(i,j,2)*dsig(2)
1090    
1091     cfc12ms = cfc12(i,j,1)*dsig(1) + cfc12(i,j,2)*dsig(2)
1092    
1093     xn2oms = xn2o (i,j,1)*dsig(1) + xn2o (i,j,2)*dsig(2)
1094    
1095     o3ms = o3 (i,j,1)*dsig(1) + o3 (i,j,2)*dsig(2)
1096    
1097     coms = co (i,j,1)*dsig(1) + co (i,j,2)*dsig(2)
1098    
1099     zco2ms = zco2 (i,j,1)*dsig(1) + zco2 (i,j,2)*dsig(2)
1100    
1101     xnoms = xno (i,j,1)*dsig(1) + xno (i,j,2)*dsig(2)
1102    
1103     xno2ms = xno2 (i,j,1)*dsig(1) + xno2 (i,j,2)*dsig(2)
1104    
1105     xn2o5ms = xn2o5(i,j,1)*dsig(1) + xn2o5(i,j,2)*dsig(2)
1106    
1107     hno3ms = hno3 (i,j,1)*dsig(1) + hno3 (i,j,2)*dsig(2)
1108    
1109     ch4ms = ch4 (i,j,1)*dsig(1) + ch4 (i,j,2)*dsig(2)
1110    
1111     ch2oms = ch2o (i,j,1)*dsig(1) + ch2o (i,j,2)*dsig(2)
1112    
1113     so2ms = so2 (i,j,1)*dsig(1) + so2 (i,j,2)*dsig(2)
1114    
1115     h2so4ms = h2so4(i,j,1)*dsig(1) + h2so4(i,j,2)*dsig(2)
1116    
1117     ! === if hfc, pfc, and sf6 are included:
1118     #ifdef INC_3GASES
1119     ! === 032698
1120     hfc134ams = hfc134a(i,j,1)*dsig(1)
1121     & + hfc134a(i,j,2)*dsig(2)
1122    
1123     pfcms = pfc(i,j,1)*dsig(1)
1124     & + pfc(i,j,2)*dsig(2)
1125    
1126     sf6ms = sf6(i,j,1)*dsig(1)
1127     & + sf6(i,j,2)*dsig(2)
1128     ! ===
1129     #endif
1130    
1131     bcms = bcarbon (i,j,1)*dsig(1) + bcarbon (i,j,2)*dsig(2)
1132     ocms = ocarbon (i,j,1)*dsig(1) + ocarbon (i,j,2)*dsig(2)
1133    
1134     c 062295
1135     c h2o2ms = h2o2 (i,j,1)*dsig(1) + h2o2 (i,j,2)*dsig(2)
1136    
1137     !
1138     #endif
1139    
1140     DO 8140 L=3,LM 6727.
1141     IF(THETA.LT.T(I,J,L)*(1.+Q(I,J,L)*RVX)) GO TO 8160 6728.
1142     PKMS=PKMS+(PK(I,J,L)*DSIG(L)) 6729.
1143     THPKMS=THPKMS+T(I,J,L)*(PK(I,J,L)*DSIG(L)) 6730.
1144     QMS=QMS+Q(I,J,L)*DSIG(L) 6731.
1145     TVMS=TVMS+T(I,J,L)*(1.+Q(I,J,L)*RVX)*(PK(I,J,L)*DSIG(L)) 6732.
1146    
1147     #if ( defined CPL_CHEM )
1148     !
1149     ! --- 03/23/95
1150     !
1151     cfc11ms = cfc11ms + cfc11(i,j,l)*dsig(l)
1152    
1153     cfc12ms = cfc12ms + cfc12(i,j,l)*dsig(l)
1154    
1155     xn2oms = xn2oms + xn2o (i,j,l)*dsig(l)
1156    
1157     o3ms = o3ms + o3 (i,j,l)*dsig(l)
1158    
1159     coms = coms + co (i,j,l)*dsig(l)
1160    
1161     zco2ms = zco2ms + zco2 (i,j,l)*dsig(l)
1162    
1163     xnoms = xnoms + xno (i,j,l)*dsig(l)
1164    
1165     xno2ms = xno2ms + xno2 (i,j,l)*dsig(l)
1166    
1167     xn2o5ms = xn2o5ms + xn2o5(i,j,l)*dsig(l)
1168    
1169     hno3ms = hno3ms + hno3 (i,j,l)*dsig(l)
1170    
1171     ch4ms = ch4ms + ch4 (i,j,l)*dsig(l)
1172    
1173     ch2oms = ch2oms + ch2o (i,j,l)*dsig(l)
1174    
1175     so2ms = so2ms + so2 (i,j,l)*dsig(l)
1176    
1177     h2so4ms = h2so4ms + h2so4(i,j,l)*dsig(l)
1178    
1179     ! === if hfc, pfc, and sf6 are included:
1180     #ifdef INC_3GASES
1181     ! === 032698
1182     hfc134ams = hfc134ams
1183     & + hfc134a(i,j,l)*dsig(l)
1184    
1185     pfcms = pfcms
1186     & + pfc(i,j,l)*dsig(l)
1187    
1188     sf6ms = sf6ms
1189     & + sf6(i,j,l)*dsig(l)
1190     ! ===
1191     #endif
1192    
1193     bcms = bcms + bcarbon (i,j,l)*dsig(l)
1194     ocms = ocms + ocarbon (i,j,l)*dsig(l)
1195    
1196     c 062295
1197     c h2o2ms = h2o2ms + h2o2 (i,j,l)*dsig(l)
1198     !
1199     #endif
1200    
1201     8140 THETA=TVMS/PKMS 6733.
1202     L=LM+1 6734.
1203     8160 LMAX=L-1 6735.
1204     RDSIGS=1./(SIGE(1)-SIGE(LMAX+1)) 6736.
1205     THM=THPKMS/PKMS 6737.
1206     QMS=QMS*RDSIGS 6738.
1207    
1208     #if ( defined CPL_CHEM )
1209     !
1210     ! --- 03/23/95
1211     !
1212     cfc11ms = cfc11ms*rdsigs
1213    
1214     cfc12ms = cfc12ms*rdsigs
1215    
1216     xn2oms = xn2oms *rdsigs
1217    
1218     o3ms = o3ms *rdsigs
1219    
1220     coms = coms *rdsigs
1221    
1222     zco2ms = zco2ms *rdsigs
1223    
1224     xnoms = xnoms *rdsigs
1225    
1226     xno2ms = xno2ms *rdsigs
1227    
1228     xn2o5ms = xn2o5ms*rdsigs
1229    
1230     hno3ms = hno3ms *rdsigs
1231    
1232     ch4ms = ch4ms *rdsigs
1233    
1234     ch2oms = ch2oms *rdsigs
1235    
1236     so2ms = so2ms *rdsigs
1237    
1238     h2so4ms = h2so4ms*rdsigs
1239    
1240     ! === if hfc, pfc, and sf6 are included:
1241     #ifdef INC_3GASES
1242     ! === 032698
1243     hfc134ams = hfc134ams*rdsigs
1244    
1245     pfcms = pfcms*rdsigs
1246    
1247     sf6ms = sf6ms*rdsigs
1248     ! ===
1249     #endif
1250    
1251     bcms = bcms*rdsigs
1252     ocms = ocms*rdsigs
1253    
1254     c 062295
1255     c h2o2ms = h2o2ms*rdsigs
1256     c
1257     !
1258     #endif
1259    
1260     BLDATA(I,J,4)=LMAX 6739.
1261     DO 8180 L=1,LMAX 6740.
1262     AJL(J,L,12)=AJL(J,L,12)+(THM-T(I,J,L))*PK(I,J,L)*P(I,J) 6741.
1263     T(I,J,L)=THM 6742.
1264    
1265     #if ( defined CPL_CHEM )
1266     !
1267     ! --- 03/23/95
1268     !
1269     cfc11(i,j,l) = cfc11ms
1270    
1271     cfc12(i,j,l) = cfc12ms
1272    
1273     xn2o(i,j,l) = xn2oms
1274    
1275     o3(i,j,l) = o3ms
1276    
1277     co(i,j,l) = coms
1278    
1279     zco2(i,j,l) = zco2ms
1280    
1281     xno(i,j,l) = xnoms
1282    
1283     xno2(i,j,l) = xno2ms
1284    
1285     xn2o5(i,j,l) = xn2o5ms
1286    
1287     hno3(i,j,l) = hno3ms
1288    
1289     ch4(i,j,l) = ch4ms
1290    
1291     ch2o(i,j,l) = ch2oms
1292    
1293     so2(i,j,l) = so2ms
1294    
1295     h2so4(i,j,l) = h2so4ms
1296    
1297     ! === if hfc, pfc, and sf6 are included:
1298     #ifdef INC_3GASES
1299     ! === 032698
1300     hfc134a(i,j,l) = hfc134ams
1301    
1302     pfc(i,j,l) = pfcms
1303    
1304     sf6(i,j,l) = sf6ms
1305     ! ===
1306     #endif
1307    
1308     bcarbon(i,j,l) = bcms
1309     ocarbon(i,j,l) = ocms
1310    
1311     c 062295
1312     c h2o2(i,j,l) = h2o2ms
1313     c
1314     !
1315     #endif
1316    
1317     8180 Q(I,J,L)=QMS 6743.
1318     IF(POLE) GO TO 8300 6744.
1319     C**** MIX MOMENTUM THROUGHOUT THE BOUNDARY LAYER AT NON-POLAR GRID BOXES6745.
1320     ID(1)=I+(J-1)*IM 6748.
1321     ID(2)=ID(1)+IM*JM*LM 6749.
1322     ID(3)=I+J*IM 6752.
1323     ID(4)=ID(3)+IM*JM*LM 6753.
1324     if(J.eq.JPR)then
1325     print *,'ID for J=',j
1326     print *,(ID(k),k=1,4)
1327     print *,'RA for J=',j
1328     print *,(RA(k),k=1,4)
1329     endif
1330     DO 8240 K=1,4 6754.
1331     UMS(K)=0. 6755.
1332     DO 8220 L=1,LMAX 6756.
1333     8220 UMS(K)=UMS(K)+UT(ID(K),1,L)*DSIG(L) 6757.
1334     8240 UMS(K)=UMS(K)*RDSIGS 6758.
1335     DO 8260 L=1,LMAX 6759.
1336     AJL(J,L,38)=AJL(J,L,38)+(UMS(1)-UT(I,J,L))*.5* 6760.
1337     * P(I,J)*RA(1) 6761.
1338     AJL(J+1,L,38)=AJL(J+1,L,38)+(UMS(3)- 6762.
1339     * UT(I,J+1,L))*P(I,J)*RA(3)*.5 6763.
1340     DO 8260 K=1,4 6764.
1341     if(J.eq.JPR)then
1342     print *,'L=',L,' K=',K
1343     print *,'ID(K)=',ID(K),' RA(K)=',RA(K)
1344     print *,'UMS(K)=',UMS(K),' UT(ID(K),1,L)=',UT(ID(K),1,L)
1345     endif
1346     8260 U(ID(K),1,L)=U(ID(K),1,L)+(UMS(K)-UT(ID(K),1,L))*RA(K) 6765.
1347     GO TO 8400 6766.
1348     C**** MIX MOMENTUM THROUGHOUT THE BOUNDARY LAYER AT POLAR GRID BOXES 6767.
1349     8300 JVPO=2 6768.
1350     IF(J.EQ.JM) JVPO=JM 6769.
1351     RAPO=2.*RAPVN(1) 6770.
1352     DO 8360 IPO=1,IM 6771.
1353     UMSPO=0. 6772.
1354     VMSPO=0. 6773.
1355     DO 8320 L=1,LMAX 6774.
1356     UMSPO=UMSPO+UT(IPO,JVPO,L)*DSIG(L) 6775.
1357     8320 VMSPO=VMSPO+VT(IPO,JVPO,L)*DSIG(L) 6776.
1358     UMSPO=UMSPO*RDSIGS 6777.
1359     VMSPO=VMSPO*RDSIGS 6778.
1360     DO 8340 L=1,LMAX 6779.
1361     U(IPO,JVPO,L)=U(IPO,JVPO,L)+(UMSPO-UT(IPO,JVPO,L))*RAPO 6780.
1362     V(IPO,JVPO,L)=V(IPO,JVPO,L)+(VMSPO-VT(IPO,JVPO,L))*RAPO 6781.
1363     8340 AJL(JVPO,L,38)=AJL(JVPO,L,38) 6782.
1364     * +(UMSPO-UT(IPO,JVPO,L))*P(1,J)*RAPO 6783.
1365     8360 CONTINUE 6784.
1366     C**** ACCUMULATE BOUNDARY LAYER DIAGNOSTICS 6785.
1367     8400 IF(MODD6.NE.0) GO TO 8500 6786.
1368     DO 8420 KR=1,4 6787.
1369     IF(I.EQ.IJD6(1,KR).AND.J.EQ.IJD6(2,KR)) GO TO 8440 6788.
1370     8420 CONTINUE 6789.
1371     GO TO 8500 6790.
1372     8440 ADAILY(IHOUR,47,KR)=ADAILY(IHOUR,47,KR)+1. 6791.
1373     ADAILY(IHOUR,48,KR)=ADAILY(IHOUR,48,KR)+LMAX 6792.
1374     8500 IM1=I 6793.
1375     do j=1,jm
1376     I=1
1377     if(J.eq.JPR.or.J.eq.-12)then
1378     print *,' J=',J,' after dry convection'
1379     print *,'U(I,J,1)=',U(I,J,1),' V(I,J,1)=',V(I,J,1)
1380     print *,'U(I,J+1,1)=',U(I,J+1,1),' V(I,J+1,1)=',V(I,J+1,1)
1381     endif
1382     enddo
1383     9000 CONTINUE 6794.
1384     do 9001 J=1,JM
1385     TSURFD(J)=TSURFD(J)+(BLDATA(1,J,2)-273.16)/24.
1386     9001 continue
1387     c write (935) ,ps4clm,
1388     c & tsl4clm,
1389     c & qs4clm,ws4clm
1390     c & ,us4clm,vs4clm
1391     RETURN 6795.
1392     990 FORMAT ('0PPBL',3I4,14F8.2) 6818.
1393     991 FORMAT ('0SURFACE ',4I4,5F10.4,3F11.7) 6819.
1394     992 FORMAT ('0',I2,10F10.4/23X,4F10.4,10X,2F10.4/ 6820.
1395     * 33X,3F10.4,10X,2F10.4) 6821.
1396     993 FORMAT ('0',I2,10F10.4/23X,7F10.4/33X,7F10.4) 6822.
1397     994 FORMAT ('0',I2,11F10.4) 6823.
1398     END 6824.

  ViewVC Help
Powered by ViewVC 1.1.22