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

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

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


Revision 1.3 - (show annotations) (download)
Wed Sep 13 15:23:04 2006 UTC (18 years, 10 months ago) by jscott
Branch: MAIN
Changes since 1.2: +1 -1 lines
comment debugging print statement

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

  ViewVC Help
Powered by ViewVC 1.1.22