/[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.10 - (hide annotations) (download)
Wed Oct 17 21:07:21 2007 UTC (17 years, 9 months ago) by jscott
Branch: MAIN
Changes since 1.9: +3 -3 lines
do pov_nep write in atmos library code to avoid endian problem

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

  ViewVC Help
Powered by ViewVC 1.1.22