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

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

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


Revision 1.5 - (hide annotations) (download)
Tue May 8 15:55:36 2007 UTC (18 years, 2 months ago) by jscott
Branch: MAIN
Changes since 1.4: +16 -21 lines
read from atm-ocean common block in atmosphere.F

1 jscott 1.5 c source 2007 sokolov users 76203 Apr 25 15:29 atmosphere.F
2 jscott 1.1
3     #include "ctrparam.h"
4    
5     ! ==========================================================
6     !
7     ! Atmosphere.F: Former main program of the MIT Global Climate and
8     ! Biogeochemistry Model.
9     !
10     ! ----------------------------------------------------------
11     !
12     ! Repacking Note: This version is combined with main.f of
13     ! several originally separated model versions such as
14     ! CliChem 3.0, CliChemNem, as well as MODEL24x11.
15     !
16     ! The chemistry module, ocean CO2 module,
17     ! and TEM module are all controlled by cpp now.
18     !
19     ! Chien Wang
20     ! MIT Joint Program for the Science and Policy
21     ! of Global Change
22     !
23     ! ----------------------------------------------------------
24     !
25     ! Revision History:
26     !
27     ! When Who What
28     ! ---- ---------- -------
29     ! 073100 Chien Wang Created from CliChem3.0 and MODEL24x11
30     ! 100500 Andrei Converted to subroutine. When called first
31     ! time does only initialization.
32     ! Creates monthly data for TEM
33     !
34     ! ==========================================================
35    
36    
37     SUBROUTINE ATMOSPHERE(DTATM,mndriver)
38     C**** MD2G04 BD2G04 MD2G04 01/02/93 0.1
39     C**** OPT(3) 0.2
40     C**** 0.3
41     C**** Dynamics and physics programs for 2-D model. 0.4
42     C**** Like D2G04 but run on work station. 0.5
43    
44     #include "BD2G04.COM" 1.
45     #include "ODIFF.COM"
46     #include "run.COM" 1.
47    
48     #include "chem_para"
49     #include "chem_com"
50    
51     #if ( defined CPL_CHEM )
52     !
53     #include "chem_tmp"
54     integer hrcnt, cnt3hr(8) ! for ozone impact
55     real sfc3hro3(nlat,8)
56     creal*4 sfc3hro3(nlat,8)
57     !
58     #endif
59    
60     COMMON/INTA/COE1(01,01,01),COE2(01,01,01) 1.1
61     c DIMENSION ACO(36,01),BSI(36,01),CCO(36,01),DSI(36,01) 1.2
62     COMMON/SPEC1/ 1.3
63     * XA(IM0,JM0),XB(IM0,JM0),YA(IM0,JM0),YB(IM0,JM0),ZA(IM0,JM0)
64     &,ZB(IM0,JM0) 1.4
65     COMMON/SPEC2/KM,KINC,COEK,C3LAND(IO0,JM0),C3OICE(IO0,JM0) 1.5
66     * ,C3LICE(IO0,JM0),WMGE(IO0,JM0),TSSFC(1,JM0,4) 1.6
67     COMMON/EPARA/VTH(JM0,LM0),WTH(JM0,LM0),VU(JM0,LM0),VV(JM0,LM0)
68     &,DQSDT(JM0,LM0) 1.7
69     * ,DWV(JM0),PHIT(JM0,LM0),TPRIM2(JM0,LM0),WU(JM0,LM0),CKS,CKN 1.8
70     * ,WQ(JM0,LM0),VQ(JM0,LM0),MRCHT 1.9
71     COMMON U,V,T,P,Q 2.
72     COMMON/WORK1/WORKX(IM0,JM0,LM0),UT(IM0,JM0,LM0),VT(IM0,JM0,LM0), 3.
73     * TT(IM0,JM0,LM0),PT(IM0,JM0),QT(IM0,JM0,LM0) 4.
74     COMMON/WORK2/UX(IM0,JM0,LM0),VX(IM0,JM0,LM0),TX(IM0,JM0,LM0)
75     &,PX(IM0,JM0) 5.
76     COMMON/OLDZO/ZMLOLD(IO0,JM0)
77     C COMMON/KEYS/KEYNR(42,50) 8.
78     CHARACTER C*4,CYEAR*4,CMMND*80 8.1
79     DIMENSION C(39),JC(100),RC(161) 8.2
80     EQUIVALENCE (JC(1),IM),(C(1),XLABEL(1),LABEL1),(RC(1),TAU) 8.3
81     CHARACTER*8 LABSSW,LABEL1,OFFSSW/'XXXXXXXX'/ 8.4
82     LOGICAL EVENT,wr25,HPRNT,TRANSR,LHORDIF 9.
83     & ,CONTRR,OBSFOR,FIRST,NOCLM
84     common/conprn/HPRNT,JPR,LPR
85     common/wrcom/wr25,TRANSR,CONTRR,OBSFOR
86     COMMON/CWMG/WMGEA(JM0),NWMGEA(JM0),RIGA(JM0),DTAV(JM0),DQAV(JM0),
87     *Z0AV(JM0),WSAV(JM0),WS0AV(JM0),TAUAV(JM0)
88     COMMON/OCN/TG3M(1,JM0,12),RTGO(1,JM0,lmo),STG3(1,JM0),DTG3(1,JM0)
89     common/SURRAD/TRSURF(JM0,4),SRSURF(JM0,4)
90     dimension RTGOAV(JM0,lmo)
91     common/tprmtg/tprmg(JM0),ntprmg(JM0)
92     common/aexpc/AEXP,ISTRT1,ISTRTCHEM,LYEAREM
93     common/mixlr/Z1OAV(JM0),NZ1OAV(JM0)
94     common/flxio/FLIO(JM0),NFLIO(JM0)
95     common/surps/srps(JM0+3),nsrps
96     c character *19 buf
97     c character *8 buf1
98     character *120 file1,file2,plotfl,nwrfl,qffile,clfile
99     CHARACTER*4 AMONTH(12),JMONTHPR
100     common/files/file1,file2,flotfl,nwrfl,qffile,clfile
101     common/PRNT1/NCOMP
102     common/Dscale/DWAV0(JM0)
103     COMMON/CO2TRND/ALFFOR,CO2TR
104     COMMON/FRMIC/ FRMDICE(JM0)
105     common/ BACKGRGHG/GHGBGR(5)
106     DATA AMONTH/'JAN','FEB','MAR','APR','MAY','JUNE','JULY'
107     & ,'AUG','SEP','OCT','NOV','DEC'/
108     dimension fluxnep(jm0)
109     logical odifcarbon
110    
111     #if ( defined CLM )
112 jscott 1.4 # include "CLM.h"
113     # if ( defined CPL_TEM )
114 jscott 1.1 C For TEM
115 jscott 1.4 # include "TEM.h"
116     # endif
117 jscott 1.1 #endif
118 jscott 1.4
119 jscott 1.1 #if ( defined CPL_OCEANCO2 )
120     #include "OCM.COM"
121     common /Garyflux/pC_atm(jm0),wind_amp,fluxco2(jm0)
122     # if ( defined ML_2D )
123     common/Garyclim/tggary(jm0),wsgary(jm0),areaml(jm0),arsrf(jm0)
124 jscott 1.4 common/Garydiff/depthml(jm0),edzcar(jm0),dzg(lmo),dzog(lmo-1),
125 jscott 1.1 &Rco2(jm0,lmo),edohd(lmo),zg(lmo),focean(jm0)
126     common /Garychem/Hg(jm0)
127     common/qfl/QFLUX(JM0,0:13),ZOAV(JM0),QFLUXT(JM0)
128     common /Garyvdif/iyearocm,vdfocm,acvdfc
129 jscott 1.4 common /Garyvlog/odifcarbon,ocarcont
130     common /Garykvct/cfkvct,edzcart(jm0)
131 jscott 1.1 # endif
132     #endif
133    
134 jscott 1.4 INTEGER dtatm, mndriver !routine arguments jrs
135 jscott 1.1
136     #if ( defined OCEAN_3D || defined ML_2D )
137 jscott 1.2 #include "AGRID.h"
138 jscott 1.4 Cjrs elimated COM file/moved elsewhere#include "HRD4OCN.COM"
139 jscott 1.1 dimension oimeltt(jm0),dhdtav(jm0),devdtav(jm0)
140     #endif
141    
142     C **** CLEAR SKY
143     common/clrsk/CLEAR(JM0),NCLR(JM0),AJCLR(JM0,12),BJCLR(JM0,12),
144     * CJCLR(JM0,12)
145 jscott 1.4 ! common/TSUR/TSURFC(JM0,0:13),TSURFT(JM0),TSURFD(JM0),DTSURF(JM0)
146     #include "TSRF.COM"
147     real TSURFW(JM0),TLANDW(JM0)
148 jscott 1.1 integer CLEAR
149    
150     common /ATCO2/atm_co2(jm0)
151    
152     #if ( defined CPL_NEM )
153     C For Emission
154     c === 031097
155     real ECH4COR(JM0),ECH4PY(JM0,12), ECH4OUT(JM0),EPJT(JM0)
156     &,ECH4CTR(JM0)
157     real EN2OCOR(JM0),EN2OPY(JM0,12), EN2OOUT(JM0),EPJTN2O(JM0)
158     &,EN2OCTR(JM0)
159     c ECH4CHIEN and EN2OCHIEN are used in chemistry model
160     common/EMFORCHIEN/ECH4CHIEN(JM0),EN2OCHIEN(JM0)
161     C For Emission
162     #endif
163     dimension NDAYMN(12)
164     data NDAYMN /31,28,31,30,31,30,31,31,30,31,30,31/
165     INTEGER JDOFM(13)
166     DATA JDOFM/0,31,59,90,120,151,181,212,243,273,304,334,365/
167     DATA FIRST/.TRUE./
168     SAVE
169    
170     C AJCLR
171     C 1 SW INC AT P0 RD (AJ(1))
172     C 2 SW ABS BELOW P0 RD (AJ(2))
173     C 3 SW ABS BELOW P1 RD (AJ(3))
174     C 4 SW ABS AT Z0 RD (AJ(6))
175     C 5 SW INC AT Z0 RD (AJ(5))
176     C 6 LW INC AT Z0 RD (AJ(67))
177     C 7 NET LW AT Z0 SF (AJ(9))
178     C 8 NET LW AT P0 RD (AJ(7))
179     C 9 NET LW AT P1 RD (AJ(8))
180     C 10 NET RAD AT P0 DG (AJ(10))
181     C 11 NET RAD AT P1 DG (AJ(11))
182     C 12 NET RAD AT Z0 DG (AJ(12))
183     C **** CLEAR SKY
184     INTFX(XTAU)=INT(XTAU*XINT+.5) 10.
185     EVENT(XTAU)=MOD(ITAU,INTFX(XTAU)).LT.IDTHR 11.
186    
187     !
188     ! --- assign input and output files
189     ! Note: Due to historical reasons, no all files are
190     ! assigned here - in case you want to search
191     ! something use
192     ! grep -i "needed characters" *.F
193     !
194     ! You have my sympathy.
195     !
196     ! Chien 080400
197     !
198    
199     !#include "assign.inc"
200     c
201     if(FIRST) then
202     c
203     CALL CLOCKS (MNOW) 12.
204     C CALL ERRSET (206,1,0,0,0,301) 13.
205     C CALL ERRSET (208,256,-1) 14.
206     MBEGIN=MNOW 14.1
207     c CALL HARMO(36,1,24,ACO,BSI,CCO,DSI,1) 14.5
208     IPFLAG=0 15.
209     C CALL ENQJOB 16.
210     CALL INPUT 17.
211    
212     #if ( defined CPL_CHEM )
213     !
214     ! --- Chemistry model
215     ! --- Set year and month index:
216     !
217     myyear = 1 !year from starting point
218     myyear = JYEAR-1976 !year from starting point
219     ! myyear = JYEAR-1891 !year from starting point
220     print *,'Emissioms for ',nchemyr,' year'
221     myyearlast = min(LYEAREM-1976,nchemyr) !last year of emission
222     ! myyearlast = min(LYEAREM-1891,nchemyr) !last year of emission
223     ! myyear = min(myyear,nchemyr)
224     myyear = min(myyear,myyearlast)
225     mymonth = 1 !month
226    
227     ihaha = 1
228     ievenodd = 0 ! even hour 0
229     ! odd hour 1
230    
231     call chembudget (p)
232     !
233     ! --- Set cfcnsf = 0.0
234     !
235     do k=1,nlev
236     cfcnsf(k) = 0.0
237     enddo
238     print *,'First year of emissions ', myyear
239     ! print *,'Emission will be fixed at year ',LYEAREM
240     print *,'Emission will be fixed at year ',1976+myyearlast
241     ! print *,'Emission will be fixed at year ',1891+myyearlast
242     !
243     #endif
244    
245     C
246     KDISK0=500+KDISK
247     ndaa=3
248     c LHORDIF=.false.
249     LHORDIF=.true.
250     if(.not.LHORDIF)print *,' NO HOR. DIFFUSION for Q'
251     if(LHORDIF)print *,' HOR. DIFFUSION for Q after COND'
252     print *,' RADIATION EVERY ',NRAD/NDYN,' HOURS'
253    
254     odifcarbon=.false.
255     #if ( defined CPL_OCEANCO2 && defined ML_2D )
256     odifcarbon=.true.
257     wind_amp=1.
258 jscott 1.4 dtco2=3600.*24.
259 jscott 1.1 c ncallgary=0
260     do j=1,jm
261     areaml(j)=dxyp(j)*(1-FDATA(1,J,2))
262     focean(j)=(1-FDATA(1,J,2))
263     DEPTHML(j)=ZOAV(j)
264     end do ! j
265     print *,' RCO2'
266 jscott 1.4 ! print 5001,((Rco2(j,k),j=1,jm),k=1,LMO)
267     print 5001,((Rco2(j,k)*1.e2,j=1,jm),k=1,LMO)
268 jscott 1.1 dzog(1)=10./SQRT(1.7010587)
269     dzg(2)=10.
270     do l=2,lmo-1
271     dzog(L)=dzog(L-1)*1.7010587
272     dzg(L+1)=dzg(L)*1.7010587
273     end do
274     zg(1)=50.
275     dzg(1)=100.
276     do l=2,lmo
277     zg(l)=zg(l-1)+0.5*(dzg(l-1)+dzg(l))
278     end do
279     do l=1,lmo
280     c edohd(l)=2.5e4/(zg(l)/zg(1))**1.0
281     C New coefficients for horizontal diffusion 11/16/00
282     edohd(l)=1.55e4-9.231e3*(atan((zg(l)-300)/50))
283     end do
284     #endif
285    
286    
287     #if ( defined CPL_CHEM )
288     !
289     ! --- Initialization of chemistry model:
290     !
291     c call cheminit(ISTRT1,T,q)
292     c 11/30/2000 ISTRTCHEM for restart of chemistry model
293     call cheminit(ISTRTCHEM,T,q)
294     c print *,'H2SO4 after cheminit ',h2so4(1,33,1)
295    
296    
297     !
298     ! --- tmp output
299     !
300     copen(124,file='OUTPUT/hro3',form='unformatted',
301     c & status='unknown')
302     hrcnt = 1
303     cnt3hr(:) = 0
304     sfc3hro3(:,:) = 0.0
305     !
306     #endif
307    
308     print *,' IRAND=',IRAND
309     print *,' NCNDS=',NCNDS
310     print *,' after INPUT MRCHT=',MRCHT
311     JDAY00=JDAY-1
312     AEXP4=AEXP
313     TAU4=TAU
314     print *,' ISTRT1=',ISTRT1
315     if(ISTART.eq.2.or.ISTRT1.eq.0)then
316     nwr=1
317     WRITE (546) AEXP4,TAU4,XLABEL
318     write(547)AEXP,nwr
319     elseif(ISTART.eq.10)then
320     read(547)AEXPX,nwr
321     print *,' NWR=',nwr
322     if(abs(AEXPX-AEXP).gt.0.05)then
323     print *,' DISAGREEMENT BETWEEN AEXPX AND AEXP FILE 47'
324     print *,' AEXPX=',AEXPX,' AEXP=',AEXP
325     stop
326     endif
327     READ (546) AEXP4
328     if(abs(AEXP4-AEXP).gt.0.05)then
329     print *,' DISAGREEMENT BETWEEN AEXP4 AND AEXP FILE 46'
330     print *,' AEXP4=',AEXP4,' AEXP=',AEXP
331     stop
332     endif
333     C***
334     do 245 nr=1,nwr-1
335     READ (546)
336     245 continue
337     endif
338     WRITE(503) OFFSSW 17.1
339     REWIND 503 17.2
340     c CALL FRTR0(IO) 18.
341     KBGN=KINC+1 18.5
342     KM2=KM*2-1 18.51
343     KM3=KM2 18.52
344     KB1=3 18.53
345     KB2=5 18.54
346     KB3=5 18.55
347     IS=IM 18.56
348     IF(KM.EQ.1) IS=1 18.57
349     FIO=IO 18.58
350     JMM2=JM-2 18.59
351     HR24=24. 18.6
352     HR12=12. 18.61
353     MSTART=MNOW+MDYN+MCNDS+MRAD+MSURF+MDIAG+MELSE 19.
354     LMT3P1=LM*3+1 20.
355     C**** INITIALIZE TIME PARAMETERS 21.
356     DTHR=DT/3600. 22.
357     IDTHR=INTFX(DTHR) 23.
358     DTFS=DT*2./3. 24.
359     DTLF=2.*DT 25.
360     NDYNO=MOD(NDYN,2) 26.
361     I24=INTFX(24.) 27.
362     NSTEP0=.5+TAUI/DTHR 28.
363     NSTEP=INT(.5+TAU/DTHR)-NSTEP0 29.
364     NSTEP1=NSTEP 29.5
365     NSTEP2=NSTEP 29.6
366     MRCHT=0. 29.7
367     ITAU=(NSTEP+NSTEP0)*IDTHR 30.
368     TAU=FLOAT(ITAU)/XINT 31.
369     IDAY=1+ITAU/I24 32.
370     TOFDAY=(ITAU-(IDAY-1)*I24)/XINT 33.
371 jscott 1.4 ! if(ISTART.eq.2.or.ISTRT1.eq.0.and..not.CONTRR)then
372     ! do 458 j=1,JM
373     ! TSURFD(j)=0.
374     ! TSURFT(j)=0.
375     ! 458 continue
376     ! endif
377 jscott 1.1 if(JDATE.eq.100)then
378     print *,JDATE,JMONTH,JYEAR
379     print *,' main before daily0'
380     print *,' T1 ocean'
381     print 5001,(ODATA(1,j,1),j=1,JM0)
382     print *,' T2 ocean'
383     print 5001,(ODATA(1,j,4),j=1,JM0)
384     print *,' T3 ocean'
385     print 5001,(ODATA(1,j,5),j=1,JM0)
386     print *,' sea ice'
387     print 5002,(ODATA(1,j,2),j=1,JM0)
388     endif
389     CALL DAILY_NEW0 34.
390     print *,' Main after DAILYNEW0 JYEAR=',JYEAR
391 jscott 1.4 print *,"DT2MGL"
392     print *,DT2MGL
393     print *,"DT2MLD"
394     print *,DT2MLD
395 jscott 1.1 #if( !defined OCEAN_3D&& !defined ML_2D )
396     CALL DAILY_OCEAN
397     print *,' AFTER DAILY_OCEAN IDAY=',IDAY,' IYEAR=',IYEAR
398     print *,' JYEAR=',JYEAR,' JDAY=',JDAY
399     print *,' JDATE=',JDATE,' JMONTH=',JMONTH
400     #endif
401     if(JDATE.eq.100)then
402     print *,JDATE,JMONTH,JYEAR
403     print *,' main after daily0'
404     print *,' T1 ocean'
405     print 5001,(ODATA(1,j,1),j=1,JM0)
406     print *,' T2 ocean'
407     print 5001,(ODATA(1,j,4),j=1,JM0)
408     print *,' T3 ocean'
409     print 5001,(ODATA(1,j,5),j=1,JM0)
410     print *,' sea ice'
411     print 5002,(ODATA(1,j,2),j=1,JM0)
412     endif
413     c print *,' Z12O'
414     c print *,(Z12O(1,j),j=1,JM0)
415     c99 CONTINUE 34.993
416     CALL CLOCKS (MLAST) 35.
417     MINC=MNOW-MLAST 36.
418     MELSE=MELSE+MINC 37.
419     PERCNT=100.*MELSE/(MSTART-MLAST+1.E-5) 38.
420     c WRITE (6,901) IDAY,TOFDAY,JDATE,JMONTH,MINC,MELSE,PERCNT,TAU 39.
421     DOPK=1. 40.
422     MODD5K=1000 41.
423     IF(TAU.GE.TAUE) GO TO 820 42.
424     HPRNT=.TRUE.
425     HPRNT=.FALSE.
426     JPR=7
427     JPR=1
428     JPR=14
429     LPR=1
430     print *,' MAIN MRCHT=',MRCHT
431     c
432     c print *,' TAUE=',TAUE
433     TAUEM=TAUE
434     FIRST=.FALSE.
435     C
436     c TAU for coupler
437     c
438     TAUATM=TAU
439     MONTHATM=JMONTH
440     JDATEATM=JDATE
441     JYEARATM=JYEAR
442     C
443 jscott 1.4 #if ( defined CPL_OCEANCO2 )
444     do j=1,jm
445     fluxco2(j)=0.0
446     enddo
447     #endif
448 jscott 1.1 #if ( defined CPL_CHEM) && ( defined CPL_TEM )
449     C For TEM
450     if(ISTRT1.eq.0) then
451     c New run
452     c Reading from flin_nep
453     read(537)adupt,temco2
454     else
455     c Restart of the run
456     c Reading from last_nep
457     read(367)adupt,temco2
458     ! & ,temch4,temn2o
459     rewind 367
460     endif
461     !
462     ! adupt= 1.459814341652516
463     ! adupt= 0.9078891180588442
464     ! adupt= 0.25
465     ! adupt= -0.1123070421398009
466     !
467 jscott 1.4 ! adupt= adupt+0.9 ! for vs23
468    
469 jscott 1.1 aduptd=adupt/(365.*JM)
470     temnepgl=0.0
471     do j=1,jm
472     temnepgl=temnepgl+temco2(j)
473     enddo
474     print *,'ADNEP=',adupt
475     print *,'Initial NEP=',adupt+temnepgl*1.e-3
476     temup0=0.0
477     #endif
478     #if (!defined OCEAN_3D && !defined ML_2D)
479     if(TRANSR)then
480     if(LMO.eq.11) then
481     call ODIFS
482     elseif(LMO.eq.12) then
483     call ODIFS12
484     else
485 jscott 1.4 Print *,' Wrong LMO',LMO
486 jscott 1.1 stop
487     endif
488     endif
489     #endif
490 jscott 1.4 !#if (defined PREDICTED_GASES)
491 jscott 1.1 #if (defined CPL_TEM || defined CPL_OCEANCO2 )
492     if(OBSFOR) then
493     call obsco2(iyear,imontha)
494     mnobco2=imonth
495     endif
496     #endif
497 jscott 1.4 !#endif
498     CJRS removed below from ocean_3d
499 jscott 1.1 #ifdef ML_2D
500     do j=1,jm
501     do i=1,io
502     CLAND4OCEAN(i,j)= C3LAND(I,J)
503     enddo
504     enddo
505     IDAYM=IDAY
506     JDAYM=JDAY
507     JDATEM=JDATE
508     JMONTHM=JMONTH
509     JYEARM=JYEAR
510     TAUML=TAU
511     TOFDAYML=TOFDAY
512     KOCEANM=KOCEAN
513     #endif
514     #if (defined OCEAN_3D || defined ML_2D)
515     do l=1,lm
516     sigfl(l)=sig(l)
517     enddo
518     print *,'SIGFL'
519     print *,sigfl
520     #endif
521     #if ( defined CPL_CHEM )
522     do j=1,jm
523     atm_co2(j)=zco2(1,j,1)
524     & *28.97296245/44.0*1.e-9
525     !ppb(m) to kg per volume base
526     & *1.e6
527     enddo
528     #else
529     if(.not.OBSFOR) then
530     CFF=1.
531     if (CO2.gt.0.0)CFF=CO2
532     do j=1,jm
533     atm_co2(j)=CFF*GHGBGR(1)
534     enddo
535     endif
536     #endif
537 jscott 1.4 #if (defined CPL_TEM || defined CPL_OCEANCO2 )
538 jscott 1.1 print *,'ATM_CO2'
539     print *,atm_co2
540 jscott 1.4 #endif
541 jscott 1.1 JDAYLAST=-1
542     ncallclm=0
543     NOCLM=.true.
544     #if ( defined CLM )
545     NOCLM=.false.
546     #endif
547     print *,' atmosphere DTATM=',DTATM
548 jscott 1.4 print *,' It is running'
549 jscott 1.1 print *,'End of atmospheric model initialization'
550     print *,' '
551     print *,' '
552     print *,' JMONTHM= ',JMONTHM
553     print *,' TOFDAYML= ',TOFDAYML
554     print *,' '
555     return
556     endif ! first
557     C**** 43.
558     C**** MAIN LOOP 44.
559     C**** 45.
560     C
561     cprint *,' atmosphere TAU=',tau
562     TAUE=TAU+DTATM
563     c HPRNT=TAU.ge.17520.0.and.TAU.lt.17545.0
564     c print *,' TAUE=',TAUE
565     #if ( defined OCEAN_3D || defined ML_2D)
566     do j=1,jm0
567     tauu(j)=0.
568     tauv(j)=0.
569     precip(j)=0.
570     evao(j)=0.
571     evai(j)=0.
572     hfluxo(j)=0.
573     hfluxi(j)=0.
574     dhfodtg(j)=0.
575     devodtg(j)=0.
576     dhfidtg(j)=0.
577     devidtg(j)=0.
578     dhfodtgeq(j)=0.
579     devodtgeq(j)=0.
580     dhfidtgeq(j)=0.
581     devidtgeq(j)=0.
582     tempr(j)=0.
583 jscott 1.4 cjrs change var name to arunoff
584 jscott 1.1 arunoff(j)=0.
585     solarinc_ice(j)=0.
586     solarnet_ice(j)=0.
587     solarinc_ocean(j)=0.
588     solarnet_ocean(j)=0.
589 jscott 1.5 Cjrs not used anymore (?) surfpr(j)=0.
590 jscott 1.1 naveo(j)=0.
591     navei(j)=0.
592     navrad(j)=0.
593     navrado(j)=0.
594     c
595     ps4ocean(j)=0.
596     do l=1,lm
597     qyz4ocean(j,l)=0.
598     tyz4ocean(j,l)=0.
599     enddo
600     c
601     enddo
602     #endif
603 jscott 1.5 #ifdef OCEAN_3D
604     C get data from atm-ocean common block
605     do j=1,jm0
606     ODATA(1,j,1)=mmsst(j)
607     ODATA(1,j,2)=mmfice(j)
608     GDATA(1,j,3)=mmtice(j)
609     GDATA(1,j,1)=mmsnowm(j)
610     ODATA(1,j,3)=mmicem(j)
611     GDATA(1,j,7)=0.5*(mmtice2(j)+mmtice1(j))
612     # ifdef CPL_OCEANCO2
613     fluxco2(j)=fluxco2(j) + dtatm*3600.*mmco2flux(j)
614     # endif
615 jscott 1.4 enddo
616 jscott 1.5 #endif
617 jscott 1.1 WLMMAX=0.0
618     C
619     100 IF(.NOT.EVENT(TAUT)) GO TO 200 46.
620     c HPRNT=TAU.ge.17520.00
621     NSTEP1=NSTEP 46.5
622     C**** WRITE RESTART INFORMATION ONTO DISK 47.
623     120 CALL RFINAL (IRAND) 48.
624     IF(NSTEP.EQ.NSTEP2) GO TO 116 48.3
625     DO 115 K=1,22 48.5
626     DO 115 J=1,JM 48.6
627     DO 115 I=2,IO 48.7
628     115 GDATA(I,J,K)=GDATA(1,J,K) 48.9
629     116 CONTINUE 48.91
630     if(wr25.and.ISTART.eq.2)then
631     print *,' main write'
632     print *,' T1 ocean'
633     print 5001,(ODATA(1,j,1),j=1,JM0)
634     print *,' T2 ocean'
635     print 5001,(ODATA(1,j,4),j=1,JM0)
636     print *,' T3 ocean'
637     print 5001,(ODATA(1,j,5),j=1,JM0)
638     REWIND KDISK0 49.
639     if(TRANSR)then
640     WRITE(KDISK0) AEXP,TAU,JC,C,RC,KEYNR,U,V,T,P,Q,ODATA,GDATA,BLDATA, 50.
641     * RQT,SRHR,TRHR,(AJ(K,1),K=1,KACC),TAU,TSSFC,CKS,CKN,WMGE,TPRIM2 51.
642 jscott 1.4 * ,MRCHT,TRSURF,SRSURF,TLANDW,TSURFW,DWAV0
643 jscott 1.1 * ,TG3M,RTGO,STG3,DTG3
644     print *,' STG3'
645     print 5001,(STG3(1,j),j=1,JM0)
646     print *,' DTG3/356'
647     print 5001,(DTG3(1,j)/365.,j=1,JM0)
648     print *,' RTGO'
649     print 5001,((RTGO(1,j,k),j=1,JM0),k=1,lmo)
650     5001 format(23f6.1)
651     5002 format(23f6.3)
652     else
653     WRITE(KDISK0) AEXP,TAU,JC,C,RC,KEYNR,U,V,T,P,Q,ODATA,GDATA,BLDATA, 50.
654     * RQT,SRHR,TRHR,(AJ(K,1),K=1,KACC),TAU,TSSFC,CKS,CKN,WMGE,TPRIM2 51.
655 jscott 1.4 * ,MRCHT,TRSURF,SRSURF,TLANDW,TSURFW,DWAV0
656 jscott 1.1 endif
657     REWIND KDISK0 52.
658     end if ! ISTART.eq.2
659     CALL CLOCKS (MNOW) 53.
660     MINC=MLAST-MNOW 54.
661     MELSE=MELSE+MINC 55.
662     PERCNT=100.*MELSE/(MSTART-MNOW+1.E-5) 56.
663     MLAST=MNOW 59.
664     C**** TEST FOR TERMINATION OF RUN 60.
665     200 READ (503,END=210) LABSSW 61.
666     c HPRNT=TAU.gt.45.0.and.TAU.lt.60.0
667     c HPRNT=TAU.gt.470.0.and.TAU.lt.550.0
668     NCOMP=0
669     210 REWIND 503 61.1
670     IF(LABSSW.EQ.LABEL1) KSS6=1 61.2
671 jscott 1.4 IF(KSS6.EQ.1) GO TO 800 62.
672 jscott 1.1 IF(TAU+.06125.GE.TAUE) GO TO 820 63.
673     JDAY00=JDAY
674     C**** IF TIME TO ZERO OUT DIAGNOSTIC ACCUMULATING ARRAYS, DO SO 64.
675     C**** (NORMALLY DONE AT THE BEGINNING OF A MONTH) 65.
676     IF (TAU.EQ.TAUI) GO TO 260 66.
677     IF(.NOT.EVENT(24.)) GO TO 290 67.
678     DO 250 K=1,13 68.
679     IF(JDAY.EQ.NDZERO(K)) GO TO 260 69.
680     250 CONTINUE 70.
681     GO TO 290 71.
682     260 CONTINUE
683     TAU0=TAU 72.
684     IDAY0=IDAY 73.
685     TOFDY0=TOFDAY 74.
686     JDATE0=JDATE 75.
687     JMNTH0=JMONTH 76.
688     JYEAR0=JYEAR 77.
689     DO 270 I=1,12 78.
690     270 IDACC(I)=0 79.
691     NODIFS=0
692     nsrps=0.
693     DO 280 K=1,KACC 80.
694     280 AJ(K,1)=0. 81.
695     do 5280 j=1,JM
696     tprmg(j)=0.
697     ntprmg(j)=0.
698     Z1OAV(j)=0.
699     NZ1OAV(j)=0
700     NFLIO(J)=0
701     FLIO(J)=0.
702     NCLR(J)=0
703     do 5282 k=1,lmo
704     RTGOAV(j,k)=0.
705     5282 continue
706     do 5281 n=1,12
707     AJCLR(J,n)=0.
708     BJCLR(J,n)=0.
709     CJCLR(J,n)=0.
710     5281 continue
711     5280 continue
712     do 5286 j=1,jm+3
713     srps(j)=0.
714     5286 continue
715    
716     c#if ( defined CPL_OCEANCO2 && defined ML_2D )
717     c call zerogary(ncallgary)
718     c call zerogary
719     c#endif
720    
721     CALL DIAG9A (1) 82.
722    
723     #if ( defined CPL_NEM )
724     C For Emission
725     c === 031097
726     DO 5655 MONTH=1,12
727     IF(JDAY.LE.JDOFM(MONTH+1)) GO TO 5656
728     5655 CONTINUE
729     5656 MNHTEM=MONTH-1
730     if(MNHTEM.eq.0)MNHTEM=12
731     do J=1,JM
732     c ECH4CHIEN(J)=temch4(J)/NDAYMN(MNHTEM)/1000.
733     c EN2OCHIEN(J)=temn2o(J)/NDAYMN(MNHTEM)/1000.
734    
735     ECH4CHIEN(J)=temch4(J)/NDAYMN(MNHTEM)
736     EN2OCHIEN(J)=temn2o(J)/NDAYMN(MNHTEM)
737    
738     enddo
739    
740     C For Emission
741     #endif
742    
743     290 CONTINUE 83.
744     C**** 84.
745     C**** INTEGRATE DYNAMIC TERMS 85.
746     if(HPRNT)then
747     print *,' main before comp1 1',' TAU=',TAU,JDATE,JMONTH
748     #include "PRNT.COM"
749     print *,' GDATA(1,7,5)=',GDATA(1,7,5),' GDATA(1,7,6)='
750     * ,GDATA(1,7,6)
751     endif
752     C**** 86.
753     MODD5D=MOD(NSTEP,NDA5D) 87.
754     IF(MODD5D.EQ.0) CALL DIAG5A (2,0) 88.
755     c IF(NDYNO.EQ.1) then
756     c print *,NDYNO
757     c print *,jm,im,LMT3P1
758     c endif
759     DO 310 J=1,JM 89.
760     DO 300 L=1,LMT3P1 90.
761     DO 300 I=1,IM 91.
762     UX(I,J,L)=U(I,J,L) 92.
763     300 UT(I,J,L)=U(I,J,L) 93.
764     DO 310 L=1,LM 94.
765     DO 310 I=1,IM 95.
766     310 QT(I,J,L)=Q(I,J,L) 96.
767     if(HPRNT)then
768     print *,' main before comp1 2',' TAU=',TAU
769     #include "PRNT.COM"
770     print *,' GDATA(1,7,5)=',GDATA(1,7,5),' GDATA(1,7,6)='
771     * ,GDATA(1,7,6)
772     endif
773     C**** INITIAL FORWARD STEP, QX = Q + .667*DT*F(Q) 97.
774     NS=0 98.
775     MRCH=0 99.
776     CALL COMP1 (UX,VX,TX,PX,Q,U,V,T,P,Q,DTFS,NS) 100.
777     if(HPRNT)then
778     print *,' main after comp1',' TAU=',TAU,' MRCH=',MRCH
779     #include "PRNT.COM"
780     endif
781     C
782     IF(NDYNO.EQ.1) GO TO 320 101.
783     C**** INITIAL BACKWARD STEP IS ODD, QT = QT + DT*F(QX) 102.
784     MRCH=-1 103.
785     CALL COMP1 (UT,VT,TT,PT,QT,UX,VX,TX,PX,Q,DT,NS) 104.
786     if(HPRNT)then
787     print *,' main after comp1',' TAU=',TAU,' MRCH=',MRCH
788     #include "PRNT.COM"
789     endif
790     C
791     GO TO 360 105.
792     C**** INITIAL BACKWARD STEP IS EVEN, Q = Q + DT*F(QX) 106.
793     320 NS=1 107.
794     MODD5K=MOD(NSTEP+NS-NDYN+NDA5K,NDA5K) 108.
795     MRCH=1 109.
796     CALL COMP1 (U,V,T,P,Q,UX,VX,TX,PX,QT,DT,NS) 110.
797     if(HPRNT)then
798     print *,' main after comp1',' TAU=',TAU,' MRCH=',MRCH
799     #include "PRNT.COM"
800     endif
801     C
802     CD DIAGA SHOULD BE CALLED HERE BUT THEN ARRAYS MUST BE CHANGED 111.
803     c
804     c
805     C**** ODD LEAP FROG STEP, QT = QT + 2*DT*F(Q) 112.
806     340 MRCH=-2 113.
807     CALL COMP1 (UT,VT,TT,PT,QT,U,V,T,P,Q,DTLF,NS) 114.
808     if(HPRNT)then
809     print *,' main after comp1',' TAU=',TAU,' MRCH=',MRCH
810     #include "PRNT.COM"
811     endif
812     C
813     C**** EVEN LEAP FROG STEP, Q = Q + 2*DT*F(QT) 115.
814     360 NS=NS+2 116.
815     MODD5K=MOD(NSTEP+NS-NDYN+NDA5K,NDA5K) 117.
816     MRCH=2 118.
817     CALL COMP1 (U,V,T,P,Q,UT,VT,TT,PT,QT,DTLF,NS) 119.
818     if(HPRNT)then
819     print *,' main after comp1',' TAU=',TAU,' MRCH=',MRCH
820     #include "PRNT.COM"
821     endif
822     C
823     IF(NS.LT.NDYN) GO TO 340 122.
824     c IF(MOD(NSTEP+NS-NDYN+NDAA,NDAA).LT.MRCH) THEN 120.
825     CALL DIAGA (UT,VT,TT,PT,QT,NOCLM) 121.
826     if(HPRNT)then
827     print *,' main after DIAGA',' TAU=',TAU,' MRCH=',MRCH
828     #include "PRNT.COM"
829     endif
830     c ENDIF
831     c IF(NS.LT.NDYN) GO TO 340 122.
832     DOPK=1. 123.
833     CALL CLOCKS (MNOW) 124.
834     MINC=MLAST-MNOW 125.
835     MDYN=MDYN+MINC 126.
836     MLAST=MNOW 127.
837     PERCNT=100.*MDYN/(MSTART-MNOW+1.E-5) 128.
838     C 130.
839     CALL DIAG9A (2) 131.
840     C**** 133.
841     C**** INTEGRATE SOURCE TERMS 134.
842     C**** 135.
843     #if (!defined PREDICTED_GASES)
844     #if (defined CPL_TEM || defined CPL_OCEANCO2 )
845     if(OBSFOR) then
846     if(JMONTH.ne.AMONTH(mnobco2)) then
847     mnobco2=mnobco2+1
848     if(mnobco2.eq.13)mnobco2=1
849     call obsco2(jyear,mnobco2)
850     endif
851     endif
852     #endif
853     #endif
854     C
855     C
856     MODRD=MOD(NSTEP,NRAD) 136.
857     MODD5S=MOD(NSTEP,NDA5S) 137.
858     IF(MODD5S.EQ.0) IDACC(8)=IDACC(8)+1 138.
859     C 139.
860     C**** CONDENSTATION, SUPER SATURATION AND MOIST CONVECTION 140.
861     if(HPRNT)then
862     print *,' main before conds',' TAU=',TAU
863     #include "PRNT.COM"
864     endif
865    
866    
867     #if ( defined CPL_CHEM )
868     !
869     ! ===== Calculate airmass at grid point
870     ! ===== Chien Wang, 092395
871     !
872    
873     call chemairmass(p4chem0)
874    
875     !
876     ! === Calculating total mass of tracers
877     ! === with long residence times:
878     !
879     call chemmass1(cfc11, cfc11mass)
880     call chemmass1(cfc12, cfc12mass)
881     call chemmass1(xn2o, xn2omass)
882     call chemmass1(zco2, zco2mass)
883     call chemmass1(ch4, ch4mass)
884    
885     ! === if hfc, pfc, and sf6 are included:
886     #ifdef INC_3GASES
887     ! === 032698
888     call chemmass1(hfc134a, hfc134amass)
889     call chemmass1(pfc, pfcmass)
890     call chemmass1(sf6, sf6mass)
891     ! ===
892     #endif
893    
894     !
895     ! === Calculating advection and eddy diffusion:
896     !
897     call chemadv0 (dt)
898    
899     !
900     ! === Calculate total n-s transport amount
901     ! === of cfc11 - temperary:
902     !
903     dth = 3600.0 ! for relexible setting of dt dt * 3.0
904     ! call chemtmp1 (dth,airmass0,p,pvv,cfc11)
905    
906     ! === Readjust mass of tracers, 1:
907     !
908     ! -----------------------------
909     ! Use tropospheric life time (yr) to calculate mass
910     ! loss of tracers and at the same time compensate
911     ! all the numerical loss back to tracer's mass
912     ! which equavalent to use adjcoe = 1.0 for
913     ! chemmass2.f
914     !
915     ! Chien Wang, September 12,1995
916     ! ------------------------------
917    
918     ! === 092595 update p
919    
920     call chemairmass(p4chem1)
921    
922     call chemmass6(46.0, 1.0, cfc11,cfc11mass)
923     call chemmass6(120.0,1.0, cfc12,cfc12mass)
924    
925     ! === 102596
926     ! === close tau type of ocean uptake of co2:
927    
928     call chemmass66(1.0, 1.0,zco2,zco2mass)
929    
930     call chemmass6(150.0,1.0,xn2o,xn2omass)
931     call chemmass2(1.0,ch4, ch4mass )
932    
933     ! === if hfc, pfc, and sf6 are included:
934     #ifdef INC_3GASES
935     ! === 032698
936     ! call chemmass6(14.6, 1.0,hfc134a, hfc134amass)
937     call chemmass2(1.0, hfc134a, hfc134amass)
938     call chemmass6(10000.0,1.0,pfc,pfcmass)
939     call chemmass6(3200.0, 1.0,sf6,sf6mass)
940     ! ===
941     #endif
942    
943     !
944     ! === Calculate tropospheric gaseous reactions
945     ! === every nhr_for_chem hours:
946     !
947     dt_chem_h = dth*float(nhr_for_chem)
948     if(ievenodd.eq.0) then
949     call chemtrop0(0, T, q, dt_chem_h, 1)
950    
951     c print *,'H2SO4 after chemtrop0 ',h2so4(1,33,1)
952     c print *,'SVIOD after chemtrop0 ',SVIOD(1,33,1)
953    
954     !
955     ! --- tmp output
956     !
957     cnt3hr(hrcnt) = cnt3hr(hrcnt) + 1
958     sfc3hro3(1:nlat,hrcnt) = sfc3hro3(1:nlat,hrcnt)
959     & + (o3(1,1:nlat,1))*29.0/48.0
960    
961     hrcnt = hrcnt + 1
962     if (hrcnt .gt. 8 ) hrcnt = 1
963     end if
964    
965     ievenodd = ievenodd + 1
966     if(ievenodd.ge.nhr_for_chem) ievenodd = 0
967    
968     !
969     ! === Calculating stratospheric processes:
970     ! === 092595
971     ! === adjust startospheric loss to whole global:
972     !
973     call chemmass1(cfc11, cfc11mass)
974     call chemmass1(cfc12, cfc12mass)
975     call chemmass1(xn2o, xn2omass)
976    
977     call chemstrat (dt)
978    
979     call chemmass2(1.0,cfc11,cfc11mass)
980     call chemmass2(1.0,cfc12,cfc12mass)
981     call chemmass2(1.0,xn2o, xn2omass)
982    
983     !
984     ! === Get total mass of chemically active species:
985     !
986     ! call chemairmass(p)
987     ! write(6,*)"P before 2nd", p
988     c print *,'H2SO4 after chemairmass ',h2so4(1,33,1)
989    
990     call chemmass1(cfc11, cfc11mass)
991     call chemmass1(cfc12, cfc12mass)
992     call chemmass1(xn2o, xn2omass)
993     call chemmass1(zco2, zco2mass)
994     call chemmass1(ch4, ch4mass)
995    
996     ! === if hfc, pfc, and sf6 are included:
997     #ifdef INC_3GASES
998     ! === 032698
999     call chemmass1(hfc134a, hfc134amass)
1000     call chemmass1(pfc, pfcmass)
1001     call chemmass1(sf6, sf6mass)
1002     ! ===
1003     #endif
1004    
1005     #endif
1006    
1007     c print *,'H2SO4 before CONDSE ',h2so4(1,33,1)
1008     CALL CONDSE(mndriver) 141.
1009     c print *,'H2SO4 after CONDSE ',h2so4(1,33,1)
1010    
1011     #if ( defined CPL_CHEM )
1012     !
1013     ! === Calculate emission once per hour:
1014     !
1015     ! timeinhr=1./(365.*24.) !hourly emission
1016     ! call chememission(timeinhr)
1017     !
1018     !
1019     ! === Print hourly
1020     !
1021     ! call chemprt
1022     !
1023     ! ==============================================
1024     #endif
1025    
1026     if(HPRNT)then
1027     print *,' main after conds',' TAU=',TAU
1028     #include "PRNT.COM"
1029     endif
1030     #if ( !defined CLM )
1031     CALL PRECIP_LAND(mndriver) 142.
1032     if(HPRNT)then
1033     print *,' main after preci',' TAU=',TAU
1034     #include "PRNT.COM"
1035     endif
1036     #endif
1037     CALL CLOCKS (MNOW) 143.
1038     MINC=MLAST-MNOW 144.
1039     MCNDS=MCNDS+MLAST-MNOW 145.
1040     MLAST=MNOW 146.
1041     C 147.
1042     CALL DIAG9A (3) 148.
1043     C**** RADIATION, SOLAR AND THERMAL 149.
1044     if(HPRNT)then
1045     print *,' main before radia',' TAU=',TAU
1046     #include "PRNT.COM"
1047     print *,' IRAND=',IRAND
1048     endif
1049     #if ( defined PREDICTED_GASES || defined PREDICTED_AEROSOL )
1050     call radia_chem
1051     #else
1052     if(OBSFOR) then
1053     CALL RADIAGSO
1054     else
1055     CALL RADIA
1056     endif
1057     #endif
1058     if(HPRNT)then
1059     print *,' main after radia',' TAU=',TAU
1060     #include "PRNT.COM"
1061     endif
1062     CALL CLOCKS (MNOW) 151.
1063     MINC=MINC+MLAST-MNOW 152.
1064     MRAD=MRAD+MLAST-MNOW 153.
1065     MLAST=MNOW 154.
1066     if(HPRNT)then
1067     print *,' main before diag9a',' TAU=',TAU
1068     #include "PRNT.COM"
1069     endif
1070     C 155.
1071     CALL DIAG9A (4) 156.
1072     if(HPRNT)then
1073     print *,' main after diag9a',' TAU=',TAU
1074     #include "PRNT.COM"
1075     endif
1076     C**** SURFACE INTERACTION AND GROUND CALCULATION 157.
1077    
1078     #if ( defined CLM )
1079     if(HPRNT)then
1080     print *,' main before surf4clm',' TAU=',TAU
1081     #include "PRNT.COM"
1082     endif
1083     CALL SUR4CLM
1084     if(HPRNT)then
1085     print *,' main after surf4clm',' TAU=',TAU
1086     #include "PRNT.COM"
1087     endif
1088 jscott 1.4 i=1
1089 jscott 1.1 do j=1,jm
1090 jscott 1.4 pcpl4clm(i,j)=pcpl4clm(i,j)*prlnd2total(j,mndriver)
1091     pcpc4clm(i,j)=pcpc4clm(i,j)*prlnd2total(j,mndriver)
1092 jscott 1.1 enddo
1093     ! print *,' main after surf4clm',' TAU=',TAU
1094     ! print ('2(12f7.2,/,11f7.2,/)'),ps4clm,pcpl4clm,
1095     ! & pcpc4clm,tpr4clm,
1096     ! & tsl4clm,
1097     ! & qs4clm,ws4clm
1098     ! & ,us4clm,vs4clm,
1099     ! & dsw4clm,
1100     ! & dlw4clm,pco24clm
1101     ! & ,swinr4clm,swvis4clm
1102     #if ( defined DATA4TEM )
1103     c if(JYEAR.gt.20)then
1104     write (935),ps4clm,pcpl4clm,
1105     & pcpc4clm,tpr4clm,
1106     & tsl4clm,
1107     & qs4clm,ws4clm
1108     & ,us4clm,vs4clm,
1109     & dsw4clm,
1110     & dlw4clm,pco24clm
1111     & ,swinr4clm,swvis4clm
1112     c endif
1113     #endif
1114    
1115     ncallclm=ncallclm+1
1116     c print *,'before clm4mit2d ncallclm=',ncallclm
1117     if(HPRNT)then
1118     print *,' main before clm4mit2d',' TAU=',TAU
1119     #include "PRNT.COM"
1120     endif
1121     call clm4mit2d
1122     if(HPRNT)then
1123     print *,' main after clm4mit2d',' TAU=',TAU
1124     #include "PRNT.COM"
1125     endif
1126    
1127     c if(JYEAR.gt.20)then
1128 jscott 1.4 ! write (934),tau,tsoiclm,snwdclm,snwcclm,
1129     ! & lwuclm,tref2mclm,tflxclm,tgndclm,
1130     ! & lhfclm,shfclm,tauxclm,tauyclm,
1131     ! & asdirclm,aldirclm,asdifclm,aldifclm,
1132     ! & sroclm,ssrclm,glrclm
1133     ! &,h2olclm,h2oiclm
1134 jscott 1.1 c endif
1135     ! print *,' main after clm4mit2d',' TAU=',TAU
1136     ! print ('2(12f7.2,/,11f7.2,/)'),tsoiclm,snwdclm,snwcclm,
1137     ! & lwuclm,tref2mclm,tflxclm,tgndclm,
1138     ! & lhfclm,shfclm,tauxclm,tauyclm,
1139     ! & asdirclm,aldirclm,asdifclm,aldifclm
1140    
1141     CALL SURF_CLM
1142     CALL SURF_OCEAN
1143     CALL GR_CLM
1144     if(HPRNT)then
1145     print *,' main after surfc',' TAU=',TAU
1146     #include "PRNT.COM"
1147     endif
1148     #else
1149     ! CALL SURF_LAND
1150     ! CALL SURF_OCEAN
1151     CALL SURFCE ! 07/14/2006
1152     if(HPRNT)then
1153     print *,' main after surfc',' TAU=',TAU
1154     #include "PRNT.COM"
1155     endif
1156     CALL GRLAND 159.
1157     #endif
1158    
1159     #if ( defined OCEAN_3D || defined ML_2D)
1160     CALL GRFOROCEAN
1161     #else
1162     CALL GROCEAN(mndriver) 159.
1163     #endif
1164     c
1165     if(HPRNT)then
1166     print *,' main after groun',' TAU=',TAU
1167     #include "PRNT.COM"
1168     endif
1169     c print *,'H2SO4 before DRYCNV ',h2so4(1,33,1)
1170     CALL DRYCNV 160.
1171     c print *,'H2SO4 after DRYCNV ',h2so4(1,33,1)
1172     if(HPRNT)then
1173     print *,' main after drycn',' TAU=',TAU
1174     #include "PRNT.COM"
1175     endif
1176     CALL CLOCKS (MNOW) 161.
1177     MINC=MINC+MLAST-MNOW 162.
1178     MSURF=MSURF+MLAST-MNOW 163.
1179     MLAST=MNOW 164.
1180     CALL DIAG9A (5) 165.
1181     C**** STRATOSPHERIC MOMENTUM DRAG 166.
1182     CALL SDRAG(WLMMAX,JWLMMAX) 167.
1183     if(HPRNT)then
1184     print *,' main after sdrag',' TAU=',TAU
1185     #include "PRNT.COM"
1186     endif
1187     CALL CLOCKS (MNOW) 168.
1188     MINC=MINC+MLAST-MNOW 169.
1189     MSURF=MSURF+MLAST-MNOW 170.
1190     MLAST=MNOW 171.
1191     C 172.
1192     CALL DIAG9A (6) 173.
1193     MSRCE=MCNDS+MRAD+MSURF 174.
1194     PERCNT=100.*MSRCE/(MSTART-MNOW+1.E-5) 175.
1195     C**** SEA LEVEL PRESSURE FILTER 177.
1196     NFILTR= NDYN 177.5
1197     IF(MFILTR.LE.0.OR.MOD(NSTEP,NFILTR).NE.0) GO TO 500 178.
1198     IDACC(10)=IDACC(10)+1 179.
1199     if(HPRNT)then
1200     print *,' main before filtr',' TAU=',TAU
1201     #include "PRNT.COM"
1202     endif
1203     C
1204     C 180.
1205     C ****************
1206     if(LHORDIF)then
1207     DTDIF=3600.
1208     if(JM.eq.24)then
1209     CALL HORDIFF(DTDIF)
1210     else if(JM.eq.46)then
1211     CALL HORDIFFALL(DTDIF)
1212     else
1213     print *,' Wromg JM=',JM
1214     stop
1215     endif
1216     end if
1217     C ****************
1218     c CALL FILTER 181.
1219     C
1220     #if ( defined CPL_CHEM )
1221     !
1222     ! === Readjust total mass of tracers:
1223     !
1224     call chemairmass(p)
1225    
1226     call chemmass2(1.00,cfc11,cfc11mass)
1227     call chemmass2(1.00,cfc12,cfc12mass)
1228     call chemmass2(1.00,xn2o ,xn2omass )
1229     call chemmass2(1.00,zco2,zco2mass)
1230     call chemmass2(1.00,ch4, ch4mass )
1231    
1232     ! === if hfc, pfc, and sf6 are included:
1233     #ifdef INC_3GASES
1234     ! === 032698
1235     call chemmass2(1.0,hfc134a, hfc134amass)
1236     call chemmass2(1.0,pfc, pfcmass)
1237     call chemmass2(1.0,sf6, sf6mass)
1238     ! ===
1239     #endif
1240    
1241     !
1242     ! === Accumulative calculation prepared for
1243     ! === carrying out monthly average:
1244     !
1245     call chemmonth1
1246     !
1247     #endif
1248    
1249     CALL CLOCKS (MNOW) 182.
1250     MDYN=MDYN+MLAST-MNOW 183.
1251     MLAST=MNOW 184.
1252     C 185.
1253     CALL DIAG9A (7) 186.
1254     C**** 187.
1255     C**** UPDATE MODEL TIME AND CALL DAILY IF REQUIRED 188.
1256     C**** 189.
1257     500 NSTEP=NSTEP+NDYN 190.
1258     ITAU=(NSTEP+NSTEP0)*IDTHR 191.
1259     TAU=FLOAT(ITAU)/XINT 192.
1260     IDAY=1+ITAU/I24 193.
1261     TOFDAYPR=TOFDAY+1.00
1262     TOFDAY=(ITAU-(IDAY-1)*I24)/XINT 194.
1263     IF(.NOT.EVENT(24.)) GO TO 550 195.
1264     C 196.
1265     do J=1,JM0
1266     TSURFW(J)=TSURFD(J)
1267 jscott 1.4 TLANDW(J)=TLANDD(J)
1268 jscott 1.1 enddo
1269    
1270     JDATECLM=JDATE
1271     JDATEPR=JDATE
1272     JMONTHPR=JMONTH
1273     JYEARPR=JYEAR
1274     CALL DAILY_NEW 197.
1275     c print *,' AFTER DAILY_NEW IDAY=',IDAY,' IYEAR=',IYEAR
1276     c print *,' JYEAR=',JYEAR,' JDAY=',JDAY
1277     #if( !defined OCEAN_3D && !defined ML_2D )
1278     CALL DAILY_OCEAN
1279     c print *,' AFTER DAILY_OCEAN IDAY=',IDAY,' IYEAR=',IYEAR
1280     c print *,' JYEAR=',JYEAR,' JDAY=',JDAY
1281     c print *,' JDATE=',JDATE,' JMONTH=',JMONTH
1282     #endif
1283     if(JDATE.eq.100)then
1284     print *,JDATE,JMONTH,JYEAR
1285     print *,' main after daily'
1286     print *,' T1 ocean'
1287     print 5001,(ODATA(1,j,1),j=1,JM0)
1288     print *,' T2 ocean'
1289     print 5001,(ODATA(1,j,4),j=1,JM0)
1290     print *,' T3 ocean'
1291     print 5001,(ODATA(1,j,5),j=1,JM0)
1292     print *,' sea ice'
1293     print 5002,(ODATA(1,j,2),j=1,JM0)
1294     endif
1295     CALL CLOCKS (MNOW) 198.
1296     MELSE=MELSE+(MLAST-MNOW) 199.
1297     MLAST=MNOW 200.
1298     NDAILY=SDAY/DT 201.
1299    
1300     #if ( defined CPL_CHEM )
1301     !
1302     ! === Calculate air mass, second step:
1303     !
1304     i=1
1305     call chemairmass(p)
1306    
1307     !
1308     ! === Calculate emission once per day:
1309     !
1310     c zxy=0.
1311     c do j=1,jm
1312     c zxy=zxy+zco2(1,j,1)
1313     c & *28.97296245/44.0*1.e-3
1314     c enddo
1315     c print *,' CO2 before emission ',zxy/jm
1316    
1317     timeinday=1./(365.) !daily emission
1318     call chememission(timeinday)
1319    
1320     c zxy=0.
1321     c do j=1,jm
1322     c zxy=zxy+zco2(1,j,1)
1323     c & *28.97296245/44.0*1.e-3
1324     c enddo
1325     c print *,' CO2 after emission ',zxy/jm
1326     !
1327     #endif
1328    
1329     C 202.
1330     CALL DIAG9A (8) 203.
1331     #if ( !defined OCEAN_3D && !defined ML_2D )
1332     IF(KOCEAN.EQ.1) THEN
1333     DO 540 J=1,JM 203.11
1334     DO 540 I=1,IM 203.12
1335     AIJ(I,J,59)=AIJ(I,J,59)+ODATA(I,J,4) 203.13
1336     540 AIJ(I,J,60)=AIJ(I,J,60)+ODATA(I,J,5) 203.14
1337     if(TRANSR)then
1338     if(LMO.eq.11) then
1339     CALL ODIFS
1340     elseif(LMO.eq.12) then
1341     CALL ODIFS12
1342     else
1343     Print *,' Wromng LMO',LMO
1344     stop
1345     endif
1346     NODIFS=NODIFS+1
1347     do 5283 j=1,JM0
1348     do 5283 k=1,lmo
1349     RTGOAV(j,k)=RTGOAV(j,k)+RTGO(1,j,k)
1350     5283 continue
1351     endif
1352     C**** RESTRUCTURE THE OCEAN LAYERS AND ELIMINATE SMALL ICE BERGS 203.15
1353     CALL OSTRUC 203.16
1354     if(JDATE.eq.41)then
1355     print *,JDATE,JMONTH,JYEAR
1356     print *,' main after ostruc'
1357     print *,'Z1O'
1358     print 5001,(Z1O(1,j),j=1,JM0)
1359     print *,' T1 ocean'
1360     print 5001,(ODATA(1,j,1),j=1,JM0)
1361     print *,' T2 top ocean'
1362     print 5001,((2.*ODATA(1,j,4)-ODATA(1,j,5)),j=1,JM0)
1363     print *,' T2 ocean'
1364     print 5001,(ODATA(1,j,4),j=1,JM0)
1365     print *,' T3 ocean'
1366     print 5001,(ODATA(1,j,5),j=1,JM0)
1367     endif
1368     ENDIF ! KOCEAN
1369     #endif
1370     #if ( defined ML_2D )
1371     if(TRANSR)then
1372     NODIFS=NODIFS+1
1373     do 5283 j=1,JM0
1374     do 5283 k=1,lmo
1375     RTGOAV(j,k)=RTGOAV(j,k)+RTGO(1,j,k)
1376     5283 continue
1377     endif
1378     #endif
1379    
1380     #if ( defined CPL_OCEANCO2 )
1381    
1382     C For OCM or 3D ocean model with carbon cycle
1383     #ifdef PREDICTED_GASES
1384     c -------
1385     c 102596
1386     c air co2 mixing ratio goes to ocean carbon model:
1387     c
1388     do j=1,jm
1389     pC_atm(j)=zco2(1,j,1)
1390     & *28.97296245/44.0*1.e-9
1391     !ppb(m) to kg per volume base
1392 jscott 1.4
1393     atm_co2(j)=pC_atm(j)*1.e6
1394    
1395 jscott 1.1 enddo ! j
1396     c
1397     c -------
1398     #else
1399    
1400     do j=1,jm
1401     pC_atm(j)=atm_co2(j)*1.e-6
1402     enddo ! j
1403    
1404     #endif
1405    
1406     #if ( defined ML_2D )
1407     CB Gary CO2 uptake by the ocean
1408    
1409     do j=1,jm
1410     if(NWMGEA(J).gt.0)then
1411     WSAV(J)=WSAV(J)/NWMGEA(J)
1412     NWMGEA(J)=0.
1413     else
1414     WSAV(J)=0.
1415     end if
1416     end do ! j
1417    
1418     do j=1,jm
1419     tggary(j)=ODATA(1,j,1)+273.16
1420     wsgary(j)=WSAV(J)
1421     WSAV(j)=0.
1422     arsrf(j)=areaml(j)*(1.-ODATA(1,j,2))
1423     DEPTHML(j)=ZOAV(J)
1424     co24ocnan(j)=co24ocnan(j)+pC_atm(j)*1.e6
1425     enddo ! j
1426     c print *,'CO2 for 2D ocean'
1427     c print *,JYEAR,JMONTH
1428     c print *,'co2=',pC_atm(27)*1.e6,' ws=',wsgary(27)
1429     c print *,'tem=',tggary(27)
1430     c print '12f7.1,/,2(11f7.1,/,),12f7.1',(pC_atm(j)*1.e6,j=1,jm)
1431     c print '12f7.1,/,2(11f7.1,/,),12f7.1',(rco2(j,1),j=1,jm)
1432     c ncallgary=ncallgary+1
1433 jscott 1.4 ! 10/28/06
1434     ! call carb_mxdlyr_chem(focean)
1435     ! call carb_airsea_flx
1436     !
1437     ! 3D ocean chemistry
1438     call carb_chem_ocmip(focean)
1439     call carb_airsea_flx(dtco2)
1440     ! 3D ocean chemistry
1441     ! 10/28/06
1442 jscott 1.1 c print *,'FCO2 ncallgary=',ncallgary
1443     c print '12f7.1,/,2(11f7.1,/,),12f7.1',
1444     c & (fluxco2(j)*12.e-15*365.,j=1,jm)
1445     #endif
1446    
1447     C For ocean carbon model
1448     c Annual oceanic CO2 uptake
1449     do j=1,jm
1450     OCUPT=OCUPT+fluxco2(j)
1451     enddo
1452 jscott 1.5 ! print *,' OCUPT=',OCUPT*12.e-15
1453 jscott 1.1
1454     #if ( defined CPL_CHEM )
1455     !
1456     ! === Calculate ocean uptake of CO2
1457     ! === once per day:
1458     !
1459     i=1
1460     call chemairmass(p) !update airmass
1461    
1462     call chemoceanco2(fluxco2)
1463     !
1464     #endif
1465    
1466 jscott 1.4 do j=1,jm
1467     fluxco2(j)=0.0
1468     enddo
1469 jscott 1.1 #endif
1470    
1471     #if ( defined CPL_TEM )
1472     !#if ( defined CLM )
1473     c print *,'JDATE for TEM=',JDATECLM
1474     do j=1,jm
1475     if(npred4tem(j).gt.0)then
1476     c pred4tem(j)=pred4tem(j)/npred4tem(j)
1477     ewvd4tem(j)=ewvd4tem(j)/npred4tem(j)
1478     pre4tem(J)=pre4tem(J)+pred4tem(j)
1479     endif
1480     c
1481     if(nt2md4tem(j).gt.0)then
1482     t2md4tem(j)=t2md4tem(j)/nt2md4tem(j)
1483     endif
1484     temp4tem(j)=temp4tem(j)+t2md4tem(j)
1485     dtem4tem(JDATECLM,j)=t2md4tem(j)
1486     c
1487     if(nradd4tem(j).gt.0)then
1488     cldd4tem(j)=cldd4tem(j)/ncldd4tem(j)
1489     swtd4tem(j)=swtd4tem(j)/nradd4tem(j)
1490     swsd4tem(j)=swsd4tem(j)/nradd4tem(j)
1491     endif
1492     sws4tem(j)=sws4tem(j)+swsd4tem(j)
1493     c
1494     enddo
1495     c
1496     do j=1,jm
1497     pred4tem(j)=0.0
1498     ewvd4tem(j)=0.0
1499     t2md4tem(j)=0.0
1500     cldd4tem(j)=0.0
1501     swtd4tem(j)=0.0
1502     swsd4tem(j)=0.0
1503     npred4tem(j)=0
1504     ncldd4tem(j)=0
1505     nradd4tem(j)=0
1506     nt2md4tem(j)=0.
1507     enddo
1508     #endif
1509    
1510     #if ( defined CPL_OCEANCO2 && defined ML_2D )
1511     C For OCM
1512    
1513 jscott 1.4 ! dtco2=3600.*24.
1514     ! cfkvct=1.0
1515     ! if (JYEAR.ge.1991) then begin
1516     ! if (JYEAR.le.2100) then begin
1517     ! cfkvct=(1.0*(2100-JYEAR)+0.25*(JYEAR-1990))/110.
1518     ! esle
1519     ! cfkvct=0.25
1520     ! endif
1521     ! endif
1522     ! do j=1,jm
1523     ! edzcatr(j)=cfkvct*edzcar(j)
1524     ! enddo
1525     call diffusco2(lmo,jm,dtco2,0.5,edzcart,depthml,focean,
1526 jscott 1.1 & dzg,dzog,rco2)
1527     call hdocean(rco2,focean,dxv,dyv,DXYP,depthml,edohd,dtco2)
1528     call avegary
1529     CB Gary CO2 uptake by the ocean
1530     #endif
1531    
1532     #if ( defined CPL_CHEM) && ( defined CPL_TEM )
1533     C take into accout land uptake form TEM for previous month
1534     do j=1,jm
1535     fluxnep(j)=aduptd+1.e-3*temco2(j)/NDAYMN(mndriver)
1536     c Annual TEM CO2 uptake
1537     TEMUPTANN=TEMUPTANN+fluxnep(j)
1538     enddo
1539     if(jdate.eq.1)then
1540     print *,'Monthly TEM uptake'
1541     c print *,mndriver,adupt,temuptann-temup0
1542     temup0=temuptann
1543     endif
1544     C
1545     c
1546     i=1
1547     call chemairmass(p) !update airmass
1548    
1549     call chemtemco2(fluxnep)
1550     C
1551     c
1552     !
1553     #endif
1554    
1555     c End of month
1556     if(JDATE.eq.1)then
1557    
1558     #if ( defined CPL_CHEM )
1559     !
1560     ! === Calculating monthly averaged mixing ratios:
1561     !
1562     call chemmonth2
1563     c print *,'Atmosphere after chemmonth2'
1564    
1565     ! === Calculate and print monthly n-s transport
1566     ! === of cfc11:
1567     !
1568     ! call chemtmp2
1569     !
1570     do nhr = 1,8
1571     sfc3hro3(1:nlat,nhr) = sfc3hro3(1:nlat,nhr)
1572     & /float(cnt3hr(nhr))
1573     #if ( defined CPL_TEM )
1574     o34tem(nhr,1:nlat)=sfc3hro3(1:nlat,nhr)
1575     #endif
1576     end do
1577     cwrite(124)sfc3hro3
1578     cnt3hr(1:8) = 0
1579     sfc3hro3(1:nlat,1:8) = 0.0
1580    
1581     ! === Writing rawdata every month:
1582     !
1583     ! call chemprt !closed 032697
1584    
1585     call chembudget (p)
1586     print *,' Atmosphre after chembudget mymonth=',mymonth
1587     ! === 09/26/94
1588     ! === Reset year and month index:
1589     !
1590     mymonth = mymonth + 1
1591     if(mymonth.gt.12)then
1592     myyear = myyear +1
1593     ! myyear = min(myyear,nchemyr)
1594     myyear = min(myyear,myyearlast)
1595     mymonth = 1
1596     ! endif ! 27/8/2005
1597    
1598     ! === 092295
1599     ! === write rawdata for possible renew run
1600     ! === at end of each month:
1601     ! === at end of each year: 27/8/2005
1602     !
1603     rewind 178
1604     print *,'For chem restart ',myyear,mymonth
1605     write(178)myyear,mymonth,airmass,
1606     & cfc11,cfc110,
1607     & cfc11m,
1608     & cfc11sd,
1609     & cfc12,cfc12m,
1610     & cfc12sd,
1611     & xn2o ,xn2om ,
1612     & xn2osd,
1613     & hfc134a,hfc134am,
1614     & pfc ,pfcm,
1615     & sf6 ,sf6m,
1616     & bcarbon,bcm,
1617     & ocarbon,ocm,
1618     & atomo ,
1619     & o1d ,
1620     & o3 ,o3m ,
1621     & co ,com ,
1622     & zco2 ,zco2m,
1623     & atomh ,
1624     & ho ,
1625     & ho2 ,hoxm ,
1626     & h2o2 ,
1627     & xno ,
1628     & xno2 ,xnoxm,
1629     & xno3 ,
1630     & xn2o5 ,xnoym,
1631     & hno3 ,
1632     & ch4 ,ch4m ,
1633     & ch3 ,
1634     & cho ,
1635     & ch2o ,
1636     & ch3o ,
1637     & ch3o2 ,
1638     & ch3o2h,
1639     & so2 ,so2m ,
1640     & hoso2 ,
1641     & so3 ,
1642     & h2so4 ,h2so4m,
1643     & sviod ,sviodm
1644     rewind 178
1645    
1646     endif
1647     #endif
1648    
1649     #if ( defined CPL_TEM )
1650     !#if ( defined CLM )
1651     do J=1,JM
1652     c
1653     #ifdef PREDICTED_GASES
1654     co24tem(j)=zco2(1,j,1)
1655     & *28.97296245/44.0*1.e-3
1656     !ppm(m) to kg per volume base
1657     #else
1658     co24tem(j)=atm_co2(j)
1659    
1660     #endif
1661     c
1662     enddo
1663     #endif
1664    
1665     endif ! end of month
1666     550 continue
1667     C END of DAilY
1668    
1669    
1670     C CALL CHECKT (11) 203.17
1671     CALL CLOCKS (MNOW) 203.18
1672     MSURF=MSURF+(MLAST-MNOW) 203.19
1673     MLAST=MNOW 203.2
1674     C**** 204.
1675     C**** WRITE INFORMATION ONTO A TAPE EVERY USET HOURS 205.
1676     C**** 206.
1677     IF(USET.LE.0.) GO TO 600 207.
1678     IF(.NOT.EVENT(USET)) GO TO 600 208.
1679     C COMPUTATIONS FOR XXXXXX 209.
1680     WRITE (520) TAU,XXXXXX 210.
1681     CALL CLOCKS (MNOW) 211.
1682     MINC=MLAST-MNOW 212.
1683     MELSE=MELSE+MINC 213.
1684     PERCNT=100.*MELSE/(MSTART-MNOW+1.E-5) 214.
1685     c WRITE (6,910) MINC,MELSE,PERCNT,TAU 215.
1686     C**** 216.
1687     C**** CALL DIAGNOSTIC ROUTINES 217.
1688     C**** 218.
1689     600 IF(MOD(NSTEP-NDYN,NDA4).EQ.0) CALL DIAG4A 219.
1690     C 220.
1691     IF(NDPRNT(1).GE.0) GO TO 610 221.
1692     C**** PRINT CURRENT DIAGNOSTICS (INCLUDING THE INITIAL CONDITIONS) 222.
1693     IF(KDIAG(1).LT.9) CALL DIAG1(NOCLM) 223.
1694     IF(KDIAG(2).LT.9) CALL DIAG2 224.
1695     IF(KDIAG(7).LT.9) CALL DIAG7P 225.
1696     IF(KDIAG(3).LT.9) CALL DIAG3 226.
1697     C 227.
1698     C 228.
1699     IF(KDIAG(4).LT.9) CALL DIAG4 229.
1700     IF(TAU.LE.TAUI+DTHR*(NDYN+.5)) CALL DIAGKN 230.
1701     NDPRNT(1)=NDPRNT(1)+1 231.
1702     C 690 changed to 691 02/21/2003
1703     610 IF(.NOT.EVENT(24.)) GO TO 691 232.
1704     C**** PRINT DIAGNOSTIC TIME AVERAGED QUANTITIES ON NDPRNT-TH DAY OF RUN 233.
1705     DO 620 K=1,13 234.
1706     IF (JDAY.EQ.NDPRNT(K)) GO TO 630 235.
1707     620 CONTINUE 236.
1708     GO TO 640 237.
1709     c 630 WRITE (6,920) 238.
1710     630 continue
1711     IF(KDIAG(1).LT.9) CALL DIAG1(NOCLM) 239.
1712     IF(KDIAG(2).LT.9) CALL DIAG2 240.
1713     IF(KDIAG(7).LT.9) CALL DIAG7P 241.
1714     IF(KDIAG(3).LT.9) CALL DIAG3 242.
1715     C 243.
1716     C 244.
1717     C IF(KDIAG(6).LT.9) CALL DIAG6 245.
1718     IF(KDIAG(4).LT.9) CALL DIAG4 246.
1719     IF(KDIAG(8).LT.9) CALL DIAG8 (IPFLAG) 247.
1720     C**** THINGS TO DO BEFORE ZEROING OUT THE ACCUMULATING ARRAYS 248.
1721     C**** (NORMALLY DONE AT THE END OF A MONTH) 249.
1722     640 DO 650 K=1,13 250.
1723     IF(JDAY.EQ.NDZERO(K)) GO TO 660 251.
1724     650 CONTINUE 252.
1725     GO TO 690 253.
1726     C**** PRINT THE KEY DIAGNOSTICS 254.
1727     660 CONTINUE 255.
1728     CALL DIAGKN 255.
1729     C**** PRINT AND ZERO OUT THE TIMING NUMBERS 256.
1730     CALL CLOCKS (MNOW) 257.
1731     MDIAG=MDIAG+(MLAST-MNOW) 258.
1732     MLAST=MNOW 259.
1733     TOTALT=.01*(MSTART-MNOW) 260.
1734     PDYN=MDYN/TOTALT 261.
1735     PCDNS=MCNDS/TOTALT 262.
1736     PRAD=MRAD/TOTALT 263.
1737     PSURF=MSURF/TOTALT 264.
1738     PDIAG=MDIAG/TOTALT 265.
1739     PELSE=MELSE/TOTALT 266.
1740     DTIME=24.*TOTALT/(60.*(TAU-TAU0)) 267.
1741     c WRITE (6,909) DTIME,PDYN,PCDNS,PRAD,PSURF,PDIAG,PELSE 268.
1742     MDYN=0 269.
1743     MCNDS=0 270.
1744     MRAD=0 271.
1745     MSURF=0 272.
1746     MDIAG=0 273.
1747     MELSE=0 274.
1748     MSTART=MNOW 275.
1749     if(TRANSR)then
1750     do 5284 j=1,JM0
1751     do 5284 k=1,lmo
1752     RTGOAV(j,k)=RTGOAV(j,k)/NODIFS
1753     5284 continue
1754     c print *,'ATM RTGOAV monthly'
1755     c print *,(RTGOAV(J,1),j=1,jm)
1756     endif ! TRANSR
1757     SPGAV=srps(jm+3)/nsrps+PTOP
1758     do 5285 j=1,JM+3
1759     GBUDG(j,38,1)=(srps(j)/nsrps+PTOP)*1013./SPGAV
1760     5285 continue
1761     c do 5287 j=1,JM+3
1762     c GBUDG(j,38,1)=GBUDG(j,37,1)*1013./GBUDG(jm+3,37,1)
1763     c5287 continue
1764 jscott 1.4 ! print *,'FRMDICE'
1765     ! print '6(1PE12.4)',FRMDICE
1766 jscott 1.1 ENKE=0.0
1767     ENPT=0.0
1768     do ii=1,4
1769     ENKE=ENKE+SPECA(1,19,ii)
1770     ENPT=ENPT+SPECA(1,20,ii)
1771     enddo
1772     c print *,'ENKE=',ENKE,' ENTP=',ENPT,' ENTT=',ENKE+ENPT
1773     print *,'ENKE=',QTABLE(JMP3,19,10),' ENTP=',QTABLE(JMP3,20,10),
1774     & ' ENTT=',QTABLE(JMP3,19,10)+QTABLE(JMP3,20,10)
1775     print *,'WLMMAX=',WLMMAX,' JWLMMAX=',JWLMMAX
1776     c print *,'AJ(*,37)'
1777     c print *,(AJ(J,37),j=1,jm)
1778     c print *,'AJ(*,28)'
1779     c print *,(AJ(J,28),j=1,jm)
1780     c IF(USEP.LE.0.) GO TO 680 276.
1781     C**** WRITE SELECTED DIAGNOSTICS ONTO A DISK DATA SET FOR PLOTTING 277.
1782     c IF(TAU.LE.TAUI+1080.) GO TO 675 278.
1783     c 670 READ (16) TAUX 279.
1784     c IF(TAU.GT.TAUX+1080.) GO TO 670 280.
1785     675 WRITE (546) AEXP4,JDATE,JMONTH,JYEAR,JDATE0,JMNTH0,JYEAR0, 281.
1786     * GBUDG,QMAPS,QTABLE,INQTAB,J1QT,INQMAP,RTGOAV
1787     print *,'From atm write(546) ',JMNTH0,JYEAR0,' ',JMONTH,JYEAR
1788     c print *,0.1*GBUDG(1,26,2),0.1*GBUDG(1,35,2)
1789     nwr=nwr+1
1790     690 continue
1791     if(JDAY.eq.1)then
1792     rewind 547
1793     write(547)AEXP,nwr
1794     rewind 547
1795     print *,'From atm write(547) ',AEXP,nwr
1796     endif
1797     if(wr25.and.JDAY.eq.1)then
1798     c Write a restart file once a year
1799     print *,'Write a restart file once a year.'
1800    
1801     #if ( defined CPL_OCEANCO2 && defined ML_2D )
1802     C Data for possible restart for OCM
1803     write(369)jyear-1,vdfocm
1804     write(369)Hg
1805     write(369)Rco2
1806     rewind 369
1807     #endif
1808    
1809     print *,' KDISK0=',KDISK0
1810     CALL RFINAL (IRAND)
1811     REWIND KDISK0
1812     if(TRANSR)then
1813     WRITE(KDISK0) AEXP,TAU,JC,C,RC,KEYNR,U,V,T,P,Q,ODATA,GDATA,
1814     * BLDATA,
1815     * RQT,SRHR,TRHR,(AJ(K,1),K=1,KACC),TAU,TSSFC,CKS,CKN,WMGE,TPRIM2
1816 jscott 1.4 * ,MRCHT,TRSURF,SRSURF,TLANDW,TSURFW,DWAV0
1817 jscott 1.1 * ,TG3M,RTGO,STG3,DTG3
1818     else
1819     WRITE(KDISK0) AEXP,TAU,JC,C,RC,KEYNR,U,V,T,P,Q,ODATA,GDATA,
1820     * BLDATA,
1821     * RQT,SRHR,TRHR,(AJ(K,1),K=1,KACC),TAU,TSSFC,CKS,CKN,WMGE,TPRIM2
1822 jscott 1.4 * ,MRCHT,TRSURF,SRSURF,TLANDW,TSURFW,DWAV0
1823 jscott 1.1 c print *,' TSURFT'
1824     c print 5001,TSURFT
1825     c print *,' TSURFW'
1826     c print 5001,TSURFW
1827     endif
1828     REWIND KDISK0
1829     KDISK=3-KDISK
1830     KDISK0=500+KDISK
1831     end if
1832     C 690 changed to 691 02/21/2003
1833     680 IF(KCOPY.LE.0) GO TO 691 284.
1834     C**** WRITE A COPY OF THE FINAL RESTART DATA SET ONTO DISK 285.
1835     CALL RFINAL (IRAND) 285.5
1836     print *,' after 680'
1837     print *,' TAU=',TAU,' IRAND=',IRAND
1838     IF(KCOPY.GT.99) GO TO 687 286.
1839     683 READ (KCOPY) TAUX 286.5
1840     IF(TAU.GT.TAUX+3240.) GO TO 683 287.
1841     685 WRITE (KCOPY) TAU,JC,C,RC,KEYNR,U,V,T,P,Q,ODATA,GDATA,BLDATA, 287.5
1842     * RQT,SRHR,TRHR,(AJ(K,1),K=1,KACC),TAU,TSSFC,CKS,CKN 288.
1843 jscott 1.4 * ,WMGE,TPRIM2,MRCHT,TRSURF,SRSURF,TLANDW,TSURFW,DWAV0
1844 jscott 1.1 * ,TG3M,RTGO,STG3,DTG3
1845     REWIND KCOPY 288.5
1846     C 690 changed to 691 02/21/2003
1847     GO TO 691 289.
1848     687 KCOPY=KCOPY-100 289.5
1849     GO TO 685 289.6
1850     C**** TIME FOR CALLING DIAGNOSTICS 290.
1851     C 690 changed to 691 02/21/2003
1852     691 CALL CLOCKS (MNOW) 291.
1853     MDIAG=MDIAG+(MLAST-MNOW) 292.
1854     MLAST=MNOW 293.
1855     780 IF(TAU.LE.TAUI+DTHR*(NDYN+.5)) GO TO 120 294.
1856     GO TO 100 295.
1857     C**** 296.
1858     C**** END OF MAIN LOOP 297.
1859     C**** 298.
1860     C**** RUN TERMINATED BECAUSE SENSE SWITCH 6 WAS TURNED ON 299.
1861     800 WRITE (6,904) 300.
1862     IF(EVENT(TAUT)) GO TO 820 301.
1863     CALL RFINAL (IRAND) 302.
1864     print *,' after 800'
1865     print *,' TAU=',TAU,' IRAND=',IRAND
1866     if(wr25) then
1867     REWIND KDISK0 303.
1868     WRITE(KDISK0) AEXP,TAU,JC,C,RC,KEYNR,U,V,T,P,Q,ODATA,GDATA,BLDATA, 304.
1869     * RQT,SRHR,TRHR,(AJ(K,1),K=1,KACC),TAU,TSSFC,CKS,CKN,WMGE,TPRIM2 305.
1870 jscott 1.4 * ,MRCHT,TRSURF,SRSURF,TLANDW,TSURFW,DWAV0
1871 jscott 1.1 end if
1872     C WRITE (6,908) 306.
1873     C**** RUN TERMINATED BECAUSE IT REACHED TAUE (OR SS6 WAS TURNED ON) 307.
1874     c 820 WRITE (6,905) TAU,IDAY,TOFDAY 308.
1875     820 continue
1876    
1877     #if ( defined OCEAN_3D || defined ML_2D )
1878     C DTATM time step of atm model in hours
1879     C precip and evap in mm/day or kg/m**2/day
1880     do j=1,jm0
1881     #if ( defined OCEAN_3D && defined CPL_OCEANCO2 )
1882     ncallatm=ncallatm+1
1883 jscott 1.4 ! 020107
1884     ! co24ocean(j)=pC_atm(j)*1.e6
1885     co24ocean(j)=atm_co2(j)
1886     ! 020107
1887 jscott 1.1 co24ocnan(j)=co24ocnan(j)+co24ocean(j)
1888     #endif
1889     #ifdef ML_2D
1890 jscott 1.4 cjrs block only MD_2D
1891 jscott 1.1 rseaice(j)=ODATA(1,J,2)
1892     #endif
1893     tauu(j)=tauu(j)/(NSURF*DTATM)
1894     tauv(j)=tauv(j)/(NSURF*DTATM)
1895     tempr(j)=tempr(j)/DTATM
1896     precip(j)=precip(j)/(DTATM/24.)
1897     fland=FDATA(1,J,2)
1898     if (fland.lt.1.0)then
1899     precip(j)=precip(j)*(1.-fland*prlnd2total(j,mndriver))
1900     & /(1.-fland)
1901     endif
1902 jscott 1.4 Cjrs surfpr(j)=surfpr(j)/(DTATM/24.)
1903 jscott 1.1 c
1904     if(naveo(j).gt.0)then
1905     hfluxo(j)=NSURF*hfluxo(j)/(NDYN*DT*naveo(j))
1906     dhfodtg(j)=dhfodtg(j)/naveo(j)
1907     dhfodtgeq(j)=dhfodtgeq(j)/naveo(j)
1908     evao(j)=NSURF*evao(j)/(NDYN*DT*naveo(j))
1909     devodtg(j)=devodtg(j)/naveo(j)
1910     devodtgeq(j)=devodtgeq(j)/naveo(j)
1911     C From mm/sec to mm/day
1912     evao(j)=24.*3600.*evao(j)
1913     devodtg(j)=24.*3600.*devodtg(j)
1914     devodtgeq(j)=24.*3600.*devodtgeq(j)
1915     endif
1916     C
1917     if(navei(j).gt.0)then
1918     hfluxi(j)=NSURF*hfluxi(j)/(NDYN*DT*navei(j))
1919     dhfidtg(j)=dhfidtg(j)/navei(j)
1920     dhfidtgeq(j)=dhfidtgeq(j)/navei(j)
1921     evai(j)=NSURF*evai(j)/(NDYN*DT*navei(j))
1922     devidtg(j)=devidtg(j)/navei(j)
1923     devidtgeq(j)=devidtgeq(j)/navei(j)
1924     C From mm/sec to mm/day
1925     evai(j)=24.*3600.*evai(j)
1926     devidtg(j)=24.*3600.*devidtg(j)
1927     devidtgeq(j)=24.*3600.*devidtgeq(j)
1928     endif
1929     C
1930     if(navrad(j).gt.0)then
1931     solarinc_ice(j)=solarinc_ice(j)/navrad(j)
1932     solarnet_ice(j)=solarnet_ice(j)/navrad(j)
1933     endif
1934     if(navrado(j).gt.0)then
1935     solarinc_ocean(j)=solarinc_ocean(j)/navrado(j)
1936     solarnet_ocean(j)=solarnet_ocean(j)/navrado(j)
1937     endif
1938     c Runoff is a flux of water from land in mm/day
1939     c not for m**2
1940 jscott 1.4 cjrs change runoff to new name arunoff
1941 jscott 1.1 arunoff(j)=arunoff(j)/(DTATM/24.)*FDATA(1,j,2)
1942     & *DXYP(J)
1943     if(NWMGEA(J).gt.0)then
1944     wsocean(J)=WSAV(J)/NWMGEA(J)
1945     NWMGEA(J)=0.
1946     else
1947     wsocean(J)=0.
1948     end if
1949     WSAV(J)=0.
1950     c dhdtav(j)=dhdtav(j)+dhfdtg(j)
1951     c devdtav(j)=devdtav(j)+devdtg(j)
1952     c
1953     c
1954     ps4ocean(j)=ps4ocean(j)/DTATM
1955     do l=1,lm
1956     qyz4ocean(j,l)=qyz4ocean(j,l)/DTATM
1957     tyz4ocean(j,l)=tyz4ocean(j,l)/DTATM
1958     enddo
1959     c
1960     c
1961     end do ! j
1962     rungl=0.0
1963     runn=0.0
1964     runt=0.0
1965     runs=0.0
1966     SLAND=0.0
1967     CLAT=20.*TWOPI/360.
1968     do j=1,jm
1969     SLAND=SLAND+FDATA(1,j,2)*DXYP(J)
1970 jscott 1.4 c jrs runoff->arunoff
1971 jscott 1.1 rungl=rungl+arunoff(j)
1972     if(LAT(J).lt.-CLAT)then
1973     runs=runs+arunoff(j)
1974     else if(LAT(J).lt.CLAT)then
1975     runt=runt+arunoff(j)
1976     else
1977     runn=runn+arunoff(j)
1978     endif
1979     enddo
1980     c print *,'RUNOFF TOFDAY=',TOFDAY
1981     c print *,rungl/SLAND,rungl,runs,runt,runn
1982     c nmonth=JMNTH0
1983     #ifdef ML_2D
1984 jscott 1.4 c jrs only ML_2D
1985 jscott 1.1 nmonth=AMONTH(mndriver)
1986     #endif
1987     jdatefl=jdate-1
1988     c if(JMONTH.ne.JMNTH0)jdatefl=NDAYMN(mndriver)
1989     c if(JMONTH.ne.JMNTH0)then
1990     c print *,'OCEAN FLUXS'
1991     c print *,JMONTH,JMNTH0,JDATE
1992     c print *,'Month=',AMONTH(mndriver),' day=',jdatefl
1993     c print *,' TAIR_OCEAN'
1994     c print *,(tairo(j),j=1,jm0)
1995     c print *,' TAIR_ICE'
1996     c print *,(tairi(j),j=1,jm0)
1997     c print *,' TAUU'
1998     c print *,(tauu(j),j=1,jm0)
1999     c print *,' TAUV'
2000     c print *,(tauv(j),j=1,jm0)
2001     c print *,' EVA-E OCEAN'
2002     c print *,(evao(j),j=1,jm0)
2003     c print *,' EVA-I OCEAN'
2004     c print *,(evai(j),j=1,jm0)
2005     c print *,' P-E ICE'
2006     c print *,(pmei(j),j=1,jm0)
2007     c print *,' HEAT FLUX OCEAN'
2008     c print *,(hfluxo(j),j=1,jm0)
2009     c print *,' DHEATFLUX_OCEAN/DTG'
2010     c print *,(dhfodtg(j),j=1,jm0)
2011     c print *,' EVA_OCEAN/DTG'
2012     c print *,(devodtg(j),j=1,jm0)
2013     c print *,' EVA_ICE/DTG'
2014     c print *,(devidtg(j),j=1,jm0)
2015     c print *,' HEAT FLUX ICE'
2016     c print *,(hfluxi(j),j=1,jm0)
2017     c print *,' DHEATFLUX_ICE/DTG'
2018     c print *,(dhfidtg(j),j=1,jm0)
2019     c print *,' DLH_OCEAN/DTG'
2020     c print *,(devidtg(j),j=1,jm0)
2021     c print *,'PS4OCEAN'
2022     c print *,(ps4ocean(j),j=1,jm0)
2023     c print *,'QYZ4OCEAN'
2024     c do l=1,lm
2025     c print *,(qyz4ocean(j,l),j=1,jm0)
2026     c enddo
2027     c endif
2028     c go to 587
2029     c write(893),nmonth,jdatefl,tempr,tauu,tauv,precip,evao,
2030     c & evai,hfluxo,dhfodtg,devodtg,hfluxi,dhfidtg,devidtg,
2031     c & solarinc_ice,solarnet_ice,rseaice
2032     #ifdef ML_2D
2033 jscott 1.4 c jrs only ML_2D
2034 jscott 1.1 do j=1,jm
2035     osst(j)=ODATA(1,j,1)
2036     aoice(j)=ODATA(1,j,3)
2037     foice(j)=ODATA(1,j,2)
2038     snowice(j)=GDATA(1,j,1)
2039     tice1(j)=GDATA(1,j,3)
2040     tice2(j)=GDATA(1,j,7)
2041     enddo
2042     #endif
2043     c write (894),nmonth,jdatefl,osst,aoice,foice,snowice,tice1,tice2
2044     587 continue
2045     #endif
2046    
2047     C
2048     c TAU for coupler
2049     c
2050     TAUATM=TAU
2051     MONTHATM=JMONTH
2052     JDATEATM=JDATE
2053     JYEARATM=JYEAR
2054     C
2055     #ifdef ML_2D
2056 jscott 1.4 Cjrs change this block to only ML_2D
2057 jscott 1.1 IDAYM=IDAY
2058     JDAYM=JDAY
2059     JDATEM=JDATE
2060     JMONTHM=JMONTH
2061     JYEARM=JYEAR
2062     TAUML=TAU
2063     TOFDAYML=TOFDAY
2064     #endif
2065     C
2066     if(JDAY.ne.JDAYLAST)then
2067     c print *,'co24ocean=',co24ocean(jm/2)
2068     if(JDATEPR.ne.0)then
2069     c WRITE (6,905) TOFDAY,JDATE,JMONTH,JYEAR
2070     WRITE (6,905) TOFDAYPR,JDATEPR,JMONTHPR,JYEARPR
2071     endif
2072 jscott 1.4 cjrs print *,'ncallclm=',ncallclm
2073 jscott 1.1 JDAYLAST=JDAY
2074     c if(ncallclm.gt.6) stop
2075     c stop
2076     endif
2077     return
2078     C CALL ENQJOB 309.
2079     C CALL ENQJOB 310.
2080     IF(KSS6.EQ.1) STOP 12 310.1
2081     IF(IPFLAG.EQ.0) STOP 13 311.
2082     STOP 1 312.
2083     C**** 313.
2084     901 FORMAT ('0CLIMATE MODEL STARTED UP',14X,'DAY',I5,', HR',F6.2,I6, 314.
2085     * A5,I27,I7,F7.1,' TAU',F9.2) 315.
2086     902 FORMAT (' DYNAMIC TERMS INTEGRATED, MRCH=',I1,4X,'DAY',I5, 316.
2087     * ', HR',F6.2,I6,A5,2I7,F7.1,23X,'TAU',F9.2) 317.
2088     903 FORMAT (' SOURCE TERMS INTEGRATED',64X,2I7,F7.1) 318.
2089     904 FORMAT ('0SENSE SWITCH 6 HAS BEEN TURNED ON.') 319.
2090     905 FORMAT (/1(1X,33('****')/)/ 320.
2091     c * ' PROGRAM HAS TERMINATED NORMALLY. TAU,IDAY,TOFDAY=',F9.2, 321.
2092     * ' ATM HAS TERMINATED NORMALLY. AT ',F9.2,I3,A5,i5,
2093     * /1(1X,33('****')/)) 322.
2094     906 FORMAT (' OUTPUT RECORD WRITTEN ON UNIT',I3,55X,2I7,F7.1, 323.
2095     * ' TAU',F9.2,' ON ',A4) 324.
2096     908 FORMAT (' OUTPUT RECORD WRITTEN ON UNIT',I3,79X,'TAU',F9.2, 325.
2097     * ' ON ',A4) 326.
2098     909 FORMAT (/'0TIME',F7.2,'(MINUTES) DYNAMICS',F5.1, 327.
2099     * ' CONDENSATION',F5.1,' RADIATION',F5.1,' SURFACE', 328.
2100     * F5.1,' DIAGNOSTICS',F5.1,' OTHER',F5.1//) 329.
2101     910 FORMAT (' INFORMATION WRITTEN ON UNIT 20',57X,2I7,F7.1, 330.
2102     * ' TAU',F9.2,' ON TAPE') 331.
2103     920 FORMAT ('1'/64(1X/)) 332.
2104     END 333.

  ViewVC Help
Powered by ViewVC 1.1.22