| blob | 6505dde7ec23a8eca2540c658c5ebcb6357afb4c |
1 MODULE module_radm
2 USE module_data_radm2
3 USE module_data_sorgam
4 integer numchem
5 parameter (numchem=numchem_radm)
6 ! ..
7 CONTAINS
8 subroutine radm_driver(id,ktau,dtstep,config_flags, &
9 gmt,julday,t_phy,moist,p8w,t8w, &
10 p_phy,chem,rho_phy,dz8w,z,z_at_w,vdrog3, &
11 ph_o31d,ph_o33p,ph_no2,ph_no3o2,ph_no3o,ph_hno2, &
12 ph_hno3,ph_hno4,ph_h2o2,ph_ch2or,ph_ch2om,ph_ch3cho, &
13 ph_ch3coch3,ph_ch3coc2h5,ph_hcocho,ph_ch3cocho, &
14 ph_hcochest,ph_ch3o2h,ph_ch3coo2h,ph_ch3ono2,ph_hcochob,&
15 ids,ide, jds,jde, kds,kde, &
16 ims,ime, jms,jme, kms,kme, &
17 its,ite, jts,jte, kts,kte )
18 USE module_configure
19 USE module_state_description
20 USE module_model_constants
21 IMPLICIT NONE
23 INTEGER, INTENT(IN ) :: id,julday, &
24 ids,ide, jds,jde, kds,kde, &
25 ims,ime, jms,jme, kms,kme, &
26 its,ite, jts,jte, kts,kte
27 INTEGER, INTENT(IN ) :: &
28 ktau
29 REAL, INTENT(IN ) :: &
30 dtstep,gmt
31 !
32 ! advected moisture variables
33 !
34 REAL, DIMENSION( ims:ime, kms:kme, jms:jme, num_moist ), &
35 INTENT(IN ) :: moist
36 !
37 ! advected chemical tracers
38 !
39 REAL, DIMENSION( ims:ime, kms:kme, jms:jme, num_chem ), &
40 INTENT(INOUT ) :: chem
42 !
43 ! arrays that hold photolysis rates
44 !
45 REAL, DIMENSION( ims:ime, kms:kme, jms:jme ), &
46 INTENT(INOUT ) :: &
47 ph_o31d,ph_o33p,ph_no2,ph_no3o2,ph_no3o,ph_hno2, &
48 ph_hno3,ph_hno4,ph_h2o2,ph_ch2or,ph_ch2om,ph_ch3cho, &
49 ph_ch3coch3,ph_ch3coc2h5,ph_hcocho,ph_ch3cocho, &
50 ph_hcochest,ph_ch3o2h,ph_ch3coo2h,ph_ch3ono2,ph_hcochob
51 !
52 ! on input from met model
53 !
54 REAL, DIMENSION( ims:ime , kms:kme , jms:jme ) , &
55 INTENT(IN ) :: &
56 t_phy, &
57 p_phy, &
58 dz8w, &
59 z , &
60 t8w,p8w,z_at_w , &
61 rho_phy
62 !
63 ! for interaction with aerosols (really is output)
64 !
65 real , INTENT(INOUT) :: &
66 vdrog3(ims:ime,kms:kme-0,jms:jme,ldrog)
67 TYPE(grid_config_rec_type), INTENT(IN ) :: config_flags
69 ! ..
70 ! .. Local Scalars ..
71 REAL :: clwchem, dt60, dtcmax, dtcmin, xtime
72 INTEGER :: i,j,k,iprt, jce, jcs, n, nr, ipr,jpr,nvr
73 ! ..
74 ! .. Local Arrays ..
75 REAL :: p(kts:kte-1), rh(kts:kte-1), rj(kts:kte-1,nreacj), &
76 t(kts:kte-1), vcinp(kts:kte-1,numchem),wlc(kts:kte-1)
77 real :: vdrog1(kts:kte-1,ldrog)
79 !
81 INTEGER :: ixhour,iaerosol_sorgam
82 real :: xhour,xmin,xtimin
83 xtime=ktau*dtstep/60.
84 ixhour=ifix(gmt+.01)+ifix(xtime/60.)
85 xhour=float(ixhour)
86 xmin=60.*gmt+(xtime-xhour*60.)
87 ipr=-10
88 jpr=-10
89 nvr=5
90 !
91 ! following is for combination radm/sorgam only, p_nu0 must be defined
92 ! in that case
93 !
94 iaerosol_sorgam=0
95 if(p_nu0.gt.1)iaerosol_sorgam=1
97 !
98 chem=max(chem,epsilc)
99 do 100 j=jts,jte
100 do 100 i=its,ite
101 vcinp=epsilc
102 vdrog1=0.
103 iprt=0
104 ! if(xtime/60.ge.2.)then
105 ! if((i.eq.12.and.j.eq.17).or. &
106 ! (i.eq.12.and.j.eq.7).or. &
107 ! (i.eq.1.and.j.eq.17))iprt=2
108 ! endif
110 ! reorder
111 !
112 ! if(iprt.eq.2)print *,'k,chem(i,k,j,p_sulf),vcinp(k,lsulf)'
113 do k=kts,kte-1
114 vcinp(k,lso2) = max(chem(i,k,j,p_so2),epsilc)
115 vcinp(k,Lsulf) = max(chem(i,k,j,p_sulf),epsilc)
116 vcinp(k,Lno2) = max(chem(i,k,j,p_no2),epsilc)
117 vcinp(k,Lno) = max(chem(i,k,j,p_no),1.e-6)
118 ! vcinp(k,Lno) = max(chem(i,k,j,p_no),epsilc)
119 vcinp(k,Lo3) = max(chem(i,k,j,p_o3),epsilc)
120 vcinp(k,Lhno3) = max(chem(i,k,j,p_hno3),epsilc)
121 vcinp(k,Lh2o2) = max(chem(i,k,j,p_h2o2),epsilc)
122 vcinp(k,Lald) = max(chem(i,k,j,p_ald),epsilc)
123 vcinp(k,Lhcho) = max(chem(i,k,j,p_hcho),epsilc)
124 vcinp(k,Lop1) = max(chem(i,k,j,p_op1),epsilc)
125 vcinp(k,Lop2) = max(chem(i,k,j,p_op2),epsilc)
126 vcinp(k,Lpaa) = max(chem(i,k,j,p_paa),epsilc)
127 vcinp(k,Lora1) = max(chem(i,k,j,p_ora1),epsilc)
128 vcinp(k,Lora2) = max(chem(i,k,j,p_ora2),epsilc)
129 vcinp(k,Lnh3) = max(chem(i,k,j,p_nh3),epsilc)
130 vcinp(k,Ln2o5) = max(chem(i,k,j,p_n2o5),epsilc)
131 vcinp(k,Lno3) = max(chem(i,k,j,p_no3),epsilc)
132 vcinp(k,Lpan) = max(chem(i,k,j,p_pan),epsilc)
133 vcinp(k,Lhc3) = max(chem(i,k,j,p_hc3),epsilc)
134 vcinp(k,Lhc5) = max(chem(i,k,j,p_hc5),epsilc)
135 vcinp(k,Lhc8) = max(chem(i,k,j,p_hc8),epsilc)
136 vcinp(k,Leth) = max(chem(i,k,j,p_eth),epsilc)
137 vcinp(k,Lco) = max(chem(i,k,j,p_co),epsilc)
138 vcinp(k,Lol2) = max(chem(i,k,j,p_ol2),epsilc)
139 vcinp(k,Lolt) = max(chem(i,k,j,p_olt),epsilc)
140 vcinp(k,Loli) = max(chem(i,k,j,p_oli),epsilc)
141 vcinp(k,Ltol) = max(chem(i,k,j,p_tol),epsilc)
142 vcinp(k,Lxyl) = max(chem(i,k,j,p_xyl),epsilc)
143 vcinp(k,Laco3) = max(chem(i,k,j,p_aco3),epsilc)
144 vcinp(k,Ltpan) = max(chem(i,k,j,p_tpan),epsilc)
145 vcinp(k,Lhono) = max(chem(i,k,j,p_hono),epsilc)
146 vcinp(k,Lhno4) = max(chem(i,k,j,p_hno4),epsilc)
147 vcinp(k,Lket) = max(chem(i,k,j,p_ket),epsilc)
148 vcinp(k,Lgly) = max(chem(i,k,j,p_gly),epsilc)
149 vcinp(k,Lmgly) = max(chem(i,k,j,p_mgly),epsilc)
150 vcinp(k,Ldcb) = max(chem(i,k,j,p_dcb),epsilc)
151 vcinp(k,Lonit) = max(chem(i,k,j,p_onit),epsilc)
152 vcinp(k,Lcsl) = max(chem(i,k,j,p_csl),epsilc)
153 vcinp(k,Lxyl) = max(chem(i,k,j,p_xyl),epsilc)
154 vcinp(k,Liso) = max(chem(i,k,j,p_iso),epsilc)
155 vcinp(k,Lho) = max(chem(i,k,j,p_ho),epsilc)
156 vcinp(k,Lho2) = max(chem(i,k,j,p_ho2),epsilc)
157 ! if(iprt.eq.2)then
158 ! print *,k,chem(i,k,j,p_sulf),vcinp(k,lsulf)
159 ! endif
160 enddo
161 !--- now do chemistry, need some input here
163 do k=kts,kte-1
164 t(k) = t_phy(i,k,j)
165 p(k) = .001*p_phy(i,k,j)
166 rh(k) = .95
167 rh(k) = MIN( .95, moist(i,k,j,p_qv) / &
168 (3.80*exp(17.27*(t_phy(i,k,j)-273.)/ &
169 (t_phy(i,k,j)-36.))/(.01*p_phy(i,k,j))))
170 rh(k)=max(.1,rh(k))
171 ! wlc(k) = moist(i,k,j,p_qc)
172 wlc(k) = 0.
173 END DO
174 dt60 = dtstep/60.
175 xtimin = max(0.,xtime-dt60)
176 dtcmin = min(.05,xtime-xtimin)
177 dtcmax = min(5.,dt60)
178 dtcmax = min(dtcmax,xtime-xtimin)
179 !
180 ! radm here is called with a vertical stack
181 !
182 jcs = kts
183 jce = kte-1
186 !
187 ! fill photolysis rates for use in radm module
188 !
189 do k=kts,kte-1
190 rj(k,1) = ph_no2(i,k,j)
191 rj(k,2) = ph_o31d(i,k,j)
192 rj(k,3) = ph_o33p(i,k,j)
193 rj(k,4) = ph_hno2(i,k,j)
194 rj(k,5) = ph_hno3(i,k,j)
195 rj(k,6) = ph_hno4(i,k,j)
196 rj(k,7) = ph_no3o2(i,k,j)
197 rj(k,8) = ph_no3o(i,k,j)
198 rj(k,9) = ph_h2o2(i,k,j)
199 rj(k,10) = ph_ch2om(i,k,j)
200 rj(k,11) = ph_ch2or(i,k,j)
201 rj(k,12) = ph_ch3cho(i,k,j)
202 rj(k,13) = ph_ch3o2h(i,k,j)
203 rj(k,14) = ph_ch3coch3(i,k,j)
204 rj(k,15) = ph_ch3coo2h(i,k,j)
205 rj(k,16) = ph_ch3coc2h5(i,k,j)
206 rj(k,17) = ph_hcocho(i,k,j)
207 rj(k,18) = ph_hcochob(i,k,j)
208 rj(k,19) = ph_ch3cocho(i,k,j)
209 rj(k,20) = ph_hcochest(i,k,j)
210 rj(k,21) = ph_ch3ono2(i,k,j)
211 END DO
212 ! print *,'before radm, i,j = ',i,j
213 ! iprt=0
214 ! if((i.eq.87.and.j.eq.15).or.(i.eq.87.and.j.eq.4))then
215 ! iprt=1
216 ! endif
217 CALL radm(rj,wlc,vcinp,t,p,rh,xtime,xtimin,kts,kte-1, &
218 iprt,dt60,dtcmax,dtcmin,vdrog1,iaerosol_sorgam)
219 ! print *,'after radm, i,j = ',i,j
220 do k=kts,kte-1
221 chem(i,k,j,p_so2) = max(vcinp(k,lso2),epsilc)
222 chem(i,k,j,p_sulf) = max(vcinp(k,Lsulf),epsilc)
223 chem(i,k,j,p_no2) = max(vcinp(k,Lno2),epsilc)
224 chem(i,k,j,p_no) = max(vcinp(k,Lno),1.e-6)
225 chem(i,k,j,p_o3) = max(vcinp(k,Lo3),epsilc)
226 chem(i,k,j,p_hno3) = max(vcinp(k,Lhno3),epsilc)
227 chem(i,k,j,p_h2o2) = max(vcinp(k,Lh2o2),epsilc)
228 chem(i,k,j,p_ald) = max(vcinp(k,Lald),epsilc)
229 chem(i,k,j,p_hcho) = max(vcinp(k,Lhcho),epsilc)
230 chem(i,k,j,p_op1) = max(vcinp(k,Lop1),epsilc)
231 chem(i,k,j,p_op2) = max(vcinp(k,Lop2),epsilc)
232 chem(i,k,j,p_paa) = max(vcinp(k,Lpaa),epsilc)
233 chem(i,k,j,p_ora1) = max(vcinp(k,Lora1),epsilc)
234 chem(i,k,j,p_ora2) = max(vcinp(k,Lora2),epsilc)
235 chem(i,k,j,p_nh3) = max(vcinp(k,Lnh3),epsilc)
236 chem(i,k,j,p_n2o5) = max(vcinp(k,Ln2o5),epsilc)
237 chem(i,k,j,p_no3) = max(vcinp(k,Lno3),epsilc)
238 chem(i,k,j,p_pan) = max(vcinp(k,Lpan),epsilc)
239 chem(i,k,j,p_hc3) = max(vcinp(k,Lhc3),epsilc)
240 chem(i,k,j,p_hc5) = max(vcinp(k,Lhc5),epsilc)
241 chem(i,k,j,p_hc8) = max(vcinp(k,Lhc8),epsilc)
242 chem(i,k,j,p_eth) = max(vcinp(k,Leth),epsilc)
243 chem(i,k,j,p_co) = max(vcinp(k,Lco),epsilc)
244 chem(i,k,j,p_ol2) = max(vcinp(k,Lol2),epsilc)
245 chem(i,k,j,p_olt) = max(vcinp(k,Lolt),epsilc)
246 chem(i,k,j,p_oli) = max(vcinp(k,Loli),epsilc)
247 chem(i,k,j,p_tol) = max(vcinp(k,Ltol),epsilc)
248 chem(i,k,j,p_xyl) = max(vcinp(k,Lxyl),epsilc)
249 chem(i,k,j,p_aco3) = max(vcinp(k,Laco3),epsilc)
250 chem(i,k,j,p_tpan) = max(vcinp(k,Ltpan),epsilc)
251 chem(i,k,j,p_hono) = max(vcinp(k,Lhono),epsilc)
252 chem(i,k,j,p_hno4) = max(vcinp(k,Lhno4),epsilc)
253 chem(i,k,j,p_ket) = max(vcinp(k,Lket),epsilc)
254 chem(i,k,j,p_gly) = max(vcinp(k,Lgly),epsilc)
255 chem(i,k,j,p_mgly) = max(vcinp(k,Lmgly),epsilc)
256 chem(i,k,j,p_dcb) = max(vcinp(k,Ldcb),epsilc)
257 chem(i,k,j,p_onit) = max(vcinp(k,Lonit),epsilc)
258 chem(i,k,j,p_csl) = max(vcinp(k,Lcsl),epsilc)
259 chem(i,k,j,p_iso) = max(vcinp(k,Liso),epsilc)
260 chem(i,k,j,p_ho) = max(vcinp(k,Lho),epsilc)
261 chem(i,k,j,p_ho2) = max(vcinp(k,Lho2),epsilc)
262 if(p_nu0.gt.1)then
263 VDROG3(i,k,j,PXYL ) = VDROG1(k,PXYL )
264 VDROG3(i,k,j,PTOL ) = VDROG1(k,PTOL )
265 VDROG3(i,k,j,PCSL1) = VDROG1(k,PCSL1)
266 VDROG3(i,k,j,PCSL2) = VDROG1(k,PCSL2)
267 VDROG3(i,k,j,PHC8 ) = VDROG1(k,PHC8 )
268 VDROG3(i,k,j,POLI1) = VDROG1(k,POLI1)
269 VDROG3(i,k,j,POLI2) = VDROG1(k,POLI2)
270 VDROG3(i,k,j,POLI3) = VDROG1(k,POLI3)
271 VDROG3(i,k,j,POLT1) = VDROG1(k,POLT1)
272 VDROG3(i,k,j,POLT2) = VDROG1(k,POLT2)
273 VDROG3(i,k,j,POLT3) = VDROG1(k,POLT3)
274 endif
275 END DO
276 ! if(iprt.eq.2)then
277 ! print *,'after radm, k,chem(i,k,j,p_sulf),vcinp(k,lsulf)'
278 ! do k=kts,kte
279 ! print *,k,chem(i,k,j,p_sulf),vcinp(k,lsulf)
280 ! enddo
281 ! endif
282 100 continue
285 END SUBROUTINE radm_driver
288 SUBROUTINE radm(rjj,wlcc,vcinp,tinp,pinp,rhinp,tstart,timemx, &
289 jcs,jce,iprt,dt60,dtcmax,dtcmin,vdrog,iaerosol_sorgam)
290 implicit none
291 ! .. Parameters ..
292 REAL, PARAMETER :: c302 = 5417.4, c303 = 19.83
293 ! ..
294 ! .. Scalar Arguments ..
295 REAL,INTENT(IN) :: dt60, dtcmax, dtcmin, timemx, tstart
296 INTEGER, INTENT(IN) :: iprt, jce, jcs
297 !
298 !
301 integer, intent (in) :: iaerosol_sorgam
302 REAL,INTENT(IN) :: rjj(jcs:jce,nreacj), &
303 wlcc(jcs:jce), tinp(jcs:jce),pinp(jcs:jce),rhinp(jcs:jce)
304 ! ..
305 real,intent (INOUT) :: vdrog(jcs:jce,ldrog),vcinp(jcs:jce,lspec)
306 ! ..
307 ! .. Local Scalars ..
308 REAL :: dtc, r, timenow, tsqrd, xk0, xk2, xk3
309 INTEGER :: i, ir, irdum, j, k, kdum, l, nr
310 ! ..
311 ! .. Local Arrays ..
312 REAL :: prdrog(jcs:jce,ldrog)
313 REAL :: aquad(jcs:jce), bquad(jcs:jce), &
314 crj(jcs:jce,nreacj), crk(jcs:jce,nreack), &
315 dum(jcs:jce), dvc(jcs:jce,ldiag), dvca(jcs:jce,ldiag), &
316 dvcg(jcs:jce,ldiag), h2o(jcs:jce,1), &
317 loss(jcs:jce,lpred), lossl(jcs:jce,lump), &
318 p(jcs:jce,1), patmot1(jcs:jce), &
319 patmot2(jcs:jce), patmot3(jcs:jce), &
320 pot(jcs:jce), prod(jcs:jce,lpred), &
321 prodl(jcs:jce,lump), rh(jcs:jce), &
322 rj(jcs:jce,nreacj), rk(jcs:jce,nreack), &
323 t(jcs:jce,1), tin(jcs:jce), to300(jcs:jce), &
324 vc(jcs:jce,1,lspec), vca(jcs:jce,1,lspec), &
325 vcg(jcs:jce,1,lspec), vcl(jcs:jce,lump), wlc(jcs:jce)
326 ! ..
327 ! .. Intrinsic Functions ..
328 INTRINSIC amax1, amin1, exp, log10
329 ! ..
330 IF (iprt==1) PRINT *, 'in radm ', jcs, jce, vcinp(jcs:jce,3), &
331 vcinp(jcs:jce,lho2)
332 IF (iprt==1) PRINT *, 'in radm ', lspec, lho2
333 IF (iprt==2) PRINT *, 'in radm ', lsulf,vcinp(jcs:jcs+5,lsulf)
334 r = 0.0820578
335 do nr=1,ldrog
336 do j=jcs,jce
337 VDROG(j,nr)=0.
338 enddo
339 enddo
340 DO nr = 1, nreacj
341 DO j = jcs, jce
342 rj(j,nr) = rjj(j,nr)
343 END DO
344 END DO
345 DO j = jcs, jce
346 wlc(j) = wlcc(j)
347 t(j,1) = tinp(j)
348 p(j,1) = pinp(j)
349 rh(j) = rhinp(j)
350 END DO
351 DO l = 1, lspec
352 DO j = jcs, jce
353 vca(j,1,l) = epsilc
354 vcg(j,1,l) = amax1(epsilc,vcinp(j,l))
355 vc(j,1,l) = amax1(epsilc,vcinp(j,l))
356 END DO
357 END DO
358 IF (iprt==1) PRINT *, ' radm', lho2, vc(jcs:jce,1,3), vc(jcs:jce,1,7), &
359 vc(jcs:jce,1,lho2)
360 DO l = 1, lpred
361 DO j = jcs, jce
362 prod(j,l) = 0.
363 loss(j,l) = epsilc
364 END DO
365 END DO
366 DO l = 1, nreacj
367 DO j = jcs, jce
368 crj(j,l) = 0.
369 END DO
370 END DO
371 DO l = 1, ldiag
372 DO j = jcs, jce
373 dvca(j,l) = epsilc
374 dvcg(j,l) = epsilc
375 dvc(j,l) = epsilc
376 END DO
377 END DO
378 DO l = 1, nreack
379 DO j = jcs, jce
380 rk(j,l) = 0.
381 crk(j,l) = epsilc
382 END DO
383 END DO
384 DO l = 1, lump
385 DO j = jcs, jce
386 vcl(j,l) = 1.e-9
387 lossl(j,l) = epsilc
388 prodl(j,l) = 0.
389 END DO
390 END DO
392 dtc = dtcmin
394 DO j = jcs, jce
395 h2o(j,1) = .611E6*rh(j)*exp(c303-c302/t(j,1))/p(j,1)
396 END DO
398 k = 1
399 i = 1
400 kdum = k
401 DO j = jcs, jce
402 tin(j) = 1./t(j,1) !RADM2.0 I --> IMRCHEM
403 !RADM2.0 I --> IMRCHEM
404 pot(j) = p(j,1)*tin(j)/101.3
405 !RADM2.0 I --> IMRCHEM
406 to300(j) = t(j,1)/300.
407 patmot1(j) = const(1)*pot(j)
408 patmot2(j) = const(2)*pot(j)
409 patmot3(j) = const(3)*pot(j)*pot(j)
410 END DO
411 DO ir = 1, nreack
412 DO j = jcs, jce
413 rk(j,ir) = thafac(ir)*exp(-eor(ir)*tin(j))*patmot2(j)
414 END DO
415 END DO
416 DO j = jcs, jce
417 !3RD ORDER
418 rk(j,16) = rk(j,16)*patmot3(j)/patmot2(j)*1.E-20
419 !1ST ORDER
420 rk(j,54) = rk(j,54)/patmot2(j)*60.
421 !1ST ORDER
422 rk(j,56) = rk(j,56)/patmot2(j)*60.
423 END DO
424 DO ir = 1, ntroe
425 irdum = itroe(ir)
426 DO j = jcs, jce
427 aquad(j) = xk0300(ir)*to300(j)**(-xntroe(ir))
428 aquad(j) = aquad(j)*patmot1(j)
429 bquad(j) = xkf300(ir)*to300(j)**(-xmtroe(ir))
430 bquad(j) = aquad(j)/bquad(j)
431 END DO
432 DO j = jcs, jce
433 rk(j,irdum) = aquad(j)/(1.+bquad(j))*0.6**(1./(1.+(log10(bquad(j)) &
434 )**2))
435 END DO
436 IF (ir>2) THEN
437 DO j = jcs, jce
438 rk(j,irdum) = rk(j,irdum)*patmot2(j)
439 END DO
440 ELSE
441 DO j = jcs, jce
442 !changed RADM2.0 IMRCHEM
443 rk(j,irdum) = rk(j,irdum)/(afac(ir)*exp(bfac(ir)/t(j,1)))*60.
444 END DO
445 END IF
446 !END DO 90 LOOP
447 END DO
448 DO j = jcs, jce
449 tsqrd = t(j,1)*t(j,1) !was Imrchem 3d I --> IMRCHEM
450 rk(j,30) = rk(j,30)*tsqrd
451 rk(j,31) = rk(j,31)*tsqrd
452 rk(j,50) = rk(j,50)*tsqrd
453 END DO
454 DO j = jcs, jce
455 rk(j,1) = patmot1(j)*6.E-34*to300(j)**(-2.3)*patmot2(j)
456 rk(j,12) = (2.2E-13*exp(620.*tin(j))+1.9E-33*patmot1(j)*exp(980.*tin &
457 (j)))*patmot2(j)
458 ! IF (iprt==1 .AND. j==jce) THEN
459 ! PRINT *, j, tin(j), patmot1(j), patmot2(j), &
460 ! 1.9E-33*patmot1(j)*exp(980.*tin(j))
461 ! PRINT *, rk(j,12), 2.2E-13*exp(620.*tin(j)), const(3), p(j,1)
462 ! END IF
463 xk0 = 7.2E-15*exp(785.*tin(j))
464 xk2 = 4.1E-16*exp(1440.*tin(j))
465 xk3 = 1.9E-33*exp(725.*tin(j))*patmot1(j)
466 rk(j,25) = (xk0+xk3/(1.+xk3/xk2))*patmot2(j)
467 rk(j,29) = (1.5E-13*(1.+2.439E-20*patmot1(j)))*patmot2(j)
468 rk(j,13) = (3.08E-34*exp(2820.*tin(j))+2.66E-34*patmot1(j)*1.E-20* &
469 exp(3180.*tin(j)))*patmot3(j)
470 END DO
471 DO j = jcs, jce
472 dum(j) = amin1(rh(j),1.)
473 dum(j) = amax1(dum(j),0.)
474 ! RK(J,137) = 1./(600.*EXP(-(DUM(J)/.28)**2.8)+5.)
475 ! RK(J,137)= CVMGP(0.2,0.0,DUM(J) - .70) ! HETEROGENOUS N2O5
476 rk(j,23) = 0.0
477 ! HOMOGENEOUS N2O5
478 IF (dum(j)-.7>=0.) THEN
479 rk(j,137) = 0.2
480 ELSE
481 rk(j,137) = 0.
482 END IF
483 END DO
484 !**
485 DO j = jcs, jce
486 vcl(j,lnox) = vc(j,1,lno) + vc(j,1,lno2)
487 vcl(j,lhox) = max(1.e-9,vc(j,1,lho) + vc(j,1,lho2))
488 vcl(j,lpao3) = vc(j,1,lpan) + vc(j,1,laco3)
489 vcl(j,ln2n3) = vc(j,1,lno3) + vc(j,1,ln2o5)
490 END DO
491 !**********************************************************************
492 ! C H E M I C A L S O L V E R
493 !**********************************************************************
494 timenow = 0.
495 10 CONTINUE
497 ! Chemical solver
498 CALL predraten(jcs,jce,iprt,crj,crk,rj,rk,vc,dvc,vca, &
499 wlc,dvca,p,h2o,dvcg,t,r)
501 CALL producn(jcs,jce,iprt,crj,crk,loss,prod,prodl,lossl, &
502 prdrog,iaerosol_sorgam)
504 CALL setdtc(jcs,jce,dtc,dtcmax,dtcmin,dt60,prod,loss,vc,timenow)
506 CALL integ1n(jcs,jce,iprt,dtc,vc,loss,prod,vcl,lossl,prodl, &
507 rk,dvc,h2o,rj,vdrog,prdrog,iaerosol_sorgam)
509 timenow = timenow + dtc
510 IF (iprt==2) PRINT *, 'end radm', timenow,vc(jcs:jce,1,lsulf)
511 IF ((timenow+0.001)<dt60) GO TO 10
512 !********************END K - LEVEL LOOP
513 IF (iprt==1) PRINT *, 'end radm', vc(jcs:jce,1,3), vc(jcs:jce,1,7)
514 DO l = 1, lspec
515 DO j = jcs, jce
516 vcinp(j,l) = amax1(epsilc,vc(j,1,l))
517 END DO
518 END DO
522 RETURN
523 END SUBROUTINE radm
525 SUBROUTINE integ1n(jcs,jce,iprt,dtc,vc,loss,prod,vcl,lossl, &
526 prodl,rk,dvc,h2o,rj,vdrog,prdrog,iaerosol_sorgam)
527 implicit none
528 !***********************************************************************
529 !* DESCRIPTION:
530 !* This is the chemical integration routine. The routine employs the
531 !* FOLLOWING APPROACH TO THE SOLUTION OF THE DIFFERENTIAL
532 !* equations by the following algorithm:
533 !*
534 !* dt * Pn + Yn (1 - dt/2 * (Ln/Yn))
535 !* Yn+1= ---------------------------------
536 !* (1 + dt/2 * (Ln/Yn))
537 !*
538 !***********************************************************************
539 !** Revision history:
540 !** NO. DATE WHO WHAT
541 !** __ ____ ___ _______________________________________
542 !** 03 89279 MCB INTEGRATED VECTORIZATION CONTROL FROM
543 !** JEFF BROOKS OF CRAY RESEARCH. THIS
544 !** INCLUDED PUSHING J LOOPS INSIDE THE
545 !** IF BLOCKS AND INSIDE REST OF CODE FOR
546 !** VECTORIZATION. (THIS CLAIMS TO GIVE
547 !** APPROXIMATELY A FACTOR OF 2.3 INCREASE
548 !** IN EXECUTION TIME ALONE!)
549 !** 02 89153 MCB REMOVED CALL TO PREDRT2 - EMBEDDED CODE
550 !** 01 89149 MCB INTEGRATIED FINAL STAND ALONE VERSION
551 !**
552 !** Input files: None
553 !**
554 !** Output files:None
555 !**
556 !** Called by: CHEM
557 !**
558 !** CALLS THE FOLLOWING SUBROUTINES: NONE
559 !**********************************************************************
560 !*****THE FOLLOWING TEMPORARY ARRAYS WERE CREATED BY JEFF BROOKS FOR
561 !*****VECTORIZABLE CODE.
562 ! .. Scalar Arguments ..
563 REAL, INTENT(IN) :: dtc
564 INTEGER, intent(in) :: iprt, jce, jcs
565 integer, intent (in) :: iaerosol_sorgam
566 real,intent(in)::prdrog(jcs:jce,ldrog)
567 real,intent(inout)::vdrog(jcs:jce,ldrog)
569 ! ..
570 ! .. Array Arguments ..
571 REAL :: dvc(jcs:jce,ldiag), h2o(jcs:jce,1), &
572 loss(jcs:jce,lpred), lossl(jcs:jce,lump), prod(jcs:jce,lpred), &
573 prodl(jcs:jce,lump), rj(jcs:jce,nreacj), rk(jcs:jce,nreack), &
574 vc(jcs:jce,1,lspec), vcl(jcs:jce,lump)
575 ! ..
576 ! .. Local Scalars ..
577 REAL :: eqno,alow
578 INTEGER :: inumho, inumhox, iter, j, l, niter
579 ! ..
580 ! .. Local Arrays ..
581 REAL :: crji(nreacj), crki(nreack), h2o2sv(jcs:jce), &
582 hoxsv(jcs:jce),jb(jcs:jce), jbb(jcs:jce), jc(jcs:jce), &
583 jloss(jcs:jce), jprod(jcs:jce), tauinv(jcs:jce), &
584 taulinv(jcs:jce)
585 ! ..
586 ! .. Intrinsic Functions ..
587 INTRINSIC exp, max
588 ! ..
589 alow=1.e-17
590 DO j = jcs, jce
591 jprod(j) = 0.
592 jb(j) = 0.
593 jc(j) = 0.
594 jloss(j) = 0.
595 jbb(j) = 0.
596 END DO
597 ! K = 1
598 ! if(iprt.eq.1)print *,'in integ1'
599 ! if(iprt.eq.2)print *,'in integ1',lsulf,intgrt(lsulf),qdtc(lsulf)
600 inumho = 3
601 inumhox = 3
602 DO l = 1, lpred
603 IF (intgrt(l)==1) THEN
605 DO j = jcs, jce
606 ! if(iprt.eq.1.and.(l.eq.3.or.l.eq.7))then
607 ! print *,'i1',l,j,vc(j,1,l),loss(j,l),cmin(l)
608 ! endif
609 tauinv(j) = (loss(j,l)/vc(j,1,l))*0.5
610 vc(j,1,l) = (vc(j,1,l)*(1.-dtc*tauinv(j))+prod(j,l)*dtc)/ &
611 (1+dtc*tauinv(j))
612 vc(j,1,l) = max(cmin(l),vc(j,1,l))
613 ! if(iprt.eq.1.and.(l.eq.3.or.l.eq.7))then
614 ! print *,'i1',l,j,vc(j,1,l),dtc
615 ! endif
616 ! if(iprt.eq.2.and.(l.eq.lsulf))then
617 ! print *,'i1',l,j,vc(j,1,l),dtc,loss(j,l),prod(j,l)
618 ! endif
619 END DO
620 END IF
621 END DO
623 !... do H2O2 exponentially and store old conc in H2O2SV for iteration
625 DO j = jcs, jce
626 h2o2sv(j) = vc(j,1,lh2o2)
627 tauinv(j) = loss(j,lh2o2)/vc(j,1,lh2o2)
628 vc(j,1,lh2o2) = prod(j,lh2o2)/tauinv(j) + (h2o2sv(j)-prod(j,lh2o2)/ &
629 tauinv(j))*exp(-dtc*tauinv(j))
630 vc(j,1,lh2o2) = max(cmin(lh2o2),vc(j,1,lh2o2))
631 END DO
634 !*****LUMPED SPECIES INTEGRATION:
635 DO l = 1, lump
636 IF (l==lhox) THEN
637 DO j = jcs, jce
638 hoxsv(j) = vcl(j,lhox)
639 taulinv(j) = lossl(j,l)/vcl(j,l)
641 vcl(j,l) = prodl(j,l)/taulinv(j) + (hoxsv(j)-prodl(j,l)/taulinv( &
642 j))*exp(-dtc*taulinv(j))
643 ! vcl(j,l)=max(cmin(l),vcl(j,l))
644 vcl(j,l)=max(1.e-9,vcl(j,l))
645 END DO
646 ELSE
647 DO j = jcs, jce
648 taulinv(j) = (lossl(j,l)/vcl(j,l))*0.5
649 vcl(j,l) = (vcl(j,l)*(1.-dtc*taulinv(j))+prodl(j,l)*dtc)/ &
650 (1.+dtc*taulinv(j))
651 ! vcl(j,l)=max(cmin(l),vcl(j,l))
652 vcl(j,l)=max(1.e-9,vcl(j,l))
653 END DO
654 END IF
655 END DO
657 DO j = jcs, jce
659 !*******[N2O5]:
660 vc(j,1,ln2o5) = vcl(j,ln2n3) - vc(j,1,lno3)
661 vc(j,1,ln2o5) = max(cmin(ln2o5),vc(j,1,ln2o5))
663 !*******[NO]:
664 ! VC(J,1,LNO) = VCL(J,LNOX) - VC(J,1,LNO2)
665 ! VC(J,1,LNO) = MAX(CMIN(LNO),VC(J,1,LNO))
666 ! 8/31/00 Calculate steady state NO from NO2,O3 ...etc., renormalize NO,
667 !23456789 123456789 123456789 123456789 123456789 123456789 123456789 12
668 jprod(j) = rj(j,1)*vc(j,1,lno2) + rj(j,4)*vc(j,1,lhono) + &
669 rj(j,8)*vc(j,1,lno3)
670 jloss(j) = rk(j,6)*vc(j,1,lo3) + rk(j,9)*vc(j,1,lho2) + &
671 rk(j,15)*vc(j,1,lho) + rk(j,16)*vc(j,1,lno) + &
672 rk(j,18)*vc(j,1,lno3) + rk(j,57)*dvc(j,lmo2) + &
673 rk(j,58)*dvc(j,lhc3p) + rk(j,60)*dvc(j,lhc5p) + &
674 rk(j,62)*dvc(j,lhc8p) + rk(j,64)*dvc(j,lol2p) + &
675 rk(j,65)*dvc(j,loltp) + rk(j,66)*dvc(j,lolip) + &
676 rk(j,67)*dvc(j,ltco3) + rk(j,68)*dvc(j,ltco3) + &
677 rk(j,69)*dvc(j,ltolp) + rk(j,70)*dvc(j,lxylp) + &
678 rk(j,71)*dvc(j,lethp) + rk(j,72)*dvc(j,lketp) + &
679 rk(j,73)*dvc(j,loln) + rk(j,131)*dvc(j,lxo2)
680 eqno = max(cmin(lno),jprod(j)/max(epsilc,jloss(j)))
681 EQNO=MAX(1.e-6,JPROD(J)/MAX(EPSILC,JLOSS(J)))
682 !*******[NO]:
683 vc(j,1,lno) = max(1.e-6,eqno*vcl(j,lnox)/(eqno+vc(j,1,lno2)))
684 vc(j,1,lno2) = vc(j,1,lno2)*vcl(j,lnox)/(eqno+vc(j,1,lno2))
686 !*******[PAN]:
687 vc(j,1,lpan) = vcl(j,lpao3) - vc(j,1,laco3)
688 vc(j,1,lpan) = max(cmin(lpan),vc(j,1,lpan))
689 ! if(iprt.ge. 1.and.j.eq.jcs)then
690 ! print *,'i1,no',j,vc(j,1,lno),vc(j,1,lno2),EQNO,JPROD(J),JLOSS(J),rj(j,1), &
691 ! vc(j,1,lno2),rj(j,8),vc(j,1,lno3),rj(j,4),vc(j,1,lhono)
692 ! endif
693 END DO
695 DO j = jcs, jce
696 jprod(j) = rj(j,4)*vc(j,1,lhono) + rj(j,5)*vc(j,1,lhno3) + &
697 2.*rj(j,9)*vc(j,1,lh2o2) + rj(j,13)*vc(j,1,lop1) + &
698 rj(j,14)*vc(j,1,lop2) + rj(j,15)*vc(j,1,lpaa) + &
699 2.E0*rk(j,5)*dvc(j,lo1d)*h2o(j,1) + .1E0*(rk(j,85)*vc(j,1,lolt)*vc &
700 (j,1,lo3)+rk(j,87)*vc(j,1,liso)*vc(j,1,lo3)) + &
701 .14E0*rk(j,86)*vc(j,1,loli)*vc(j,1,lo3)
703 jloss(j) = rk(j,7)*vc(j,1,lo3) + rk(j,14)*vc(j,1,lh2o2) + &
704 rk(j,15)*vc(j,1,lno) + rk(j,24)*vc(j,1,lno2) + &
705 rk(j,25)*vc(j,1,lhno3) + rk(j,26)*vc(j,1,lhno4) + &
706 rk(j,28)*vc(j,1,lso2) + rk(j,29)*vc(j,1,lco) + rk(j,30)*ch4 + &
707 rk(j,31)*vc(j,1,leth) + rk(j,32)*vc(j,1,lhc3) + &
708 rk(j,33)*vc(j,1,lhc5) + rk(j,34)*vc(j,1,lhc8) + &
709 rk(j,35)*vc(j,1,lol2) + rk(j,36)*vc(j,1,lolt) + &
710 rk(j,37)*vc(j,1,loli) + rk(j,38)*vc(j,1,ltol) + &
711 rk(j,39)*vc(j,1,lxyl) + rk(j,40)*vc(j,1,lcsl) + &
712 rk(j,41)*vc(j,1,lhcho) + rk(j,42)*vc(j,1,lald) + &
713 rk(j,43)*vc(j,1,lket) + rk(j,44)*vc(j,1,lgly) + &
714 rk(j,45)*vc(j,1,lmgly) + rk(j,46)*vc(j,1,ldcb) + &
715 .5E0*rk(j,47)*vc(j,1,lop1) + .5E0*rk(j,48)*vc(j,1,lop2) + &
716 rk(j,49)*vc(j,1,lpaa) + rk(j,50)*vc(j,1,lpan) + &
717 rk(j,51)*vc(j,1,lonit) + rk(j,52)*vc(j,1,liso)
718 jbb(j) = rk(j,8)*vc(j,1,lo3) + rk(j,9)*vc(j,1,lno) + jloss(j)
719 END DO
720 ! if(iprt.eq.1)print *,'jloss '
722 DO niter = 1, inumhox
723 DO j = jcs, jce
724 crji(9) = rj(j,9)*vc(j,1,lh2o2)
725 crki(10) = rk(j,10)*vc(j,1,lho2)*vc(j,1,lno2)
726 crki(12) = rk(j,12)*vc(j,1,lho2)*vc(j,1,lho2)
727 crki(13) = rk(j,13)*vc(j,1,lho2)*vc(j,1,lho2)*h2o(j,1)
728 crki(14) = rk(j,14)*vc(j,1,lh2o2)*vc(j,1,lho)
729 crki(15) = rk(j,15)*vc(j,1,lno)*vc(j,1,lho)
730 crki(20) = rk(j,20)*vc(j,1,lno3)*vc(j,1,lho2)
731 crki(24) = rk(j,24)*vc(j,1,lho)*vc(j,1,lno2)
732 crki(25) = rk(j,25)*vc(j,1,lho)*vc(j,1,lhno3)
733 crki(26) = rk(j,26)*vc(j,1,lho)*vc(j,1,lhno4)
734 crki(27) = rk(j,27)*vc(j,1,lho)*vc(j,1,lho2)
735 crki(30) = rk(j,30)*ch4*vc(j,1,lho)
736 crki(31) = rk(j,31)*vc(j,1,leth)*vc(j,1,lho)
737 crki(32) = rk(j,32)*vc(j,1,lhc3)*vc(j,1,lho)
738 crki(33) = rk(j,33)*vc(j,1,lhc5)*vc(j,1,lho)
739 crki(34) = rk(j,34)*vc(j,1,lhc8)*vc(j,1,lho)
740 crki(35) = rk(j,35)*vc(j,1,lol2)*vc(j,1,lho)
741 crki(36) = rk(j,36)*vc(j,1,lolt)*vc(j,1,lho)
742 crki(37) = rk(j,37)*vc(j,1,loli)*vc(j,1,lho)
743 crki(38) = rk(j,38)*vc(j,1,ltol)*vc(j,1,lho)
744 crki(39) = rk(j,39)*vc(j,1,lxyl)*vc(j,1,lho)
745 crki(40) = rk(j,40)*vc(j,1,lcsl)*vc(j,1,lho)
746 crki(42) = rk(j,42)*vc(j,1,lald)*vc(j,1,lho)
747 crki(43) = rk(j,43)*vc(j,1,lket)*vc(j,1,lho)
748 crki(45) = rk(j,45)*vc(j,1,lmgly)*vc(j,1,lho)
749 crki(46) = rk(j,46)*vc(j,1,ldcb)*vc(j,1,lho)
750 crki(47) = rk(j,47)*vc(j,1,lop1)*vc(j,1,lho)
751 crki(48) = rk(j,48)*vc(j,1,lop2)*vc(j,1,lho)
752 crki(49) = rk(j,49)*vc(j,1,lpaa)*vc(j,1,lho)
753 crki(50) = rk(j,50)*vc(j,1,lpan)*vc(j,1,lho)
754 crki(51) = rk(j,51)*vc(j,1,lonit)*vc(j,1,lho)
755 crki(52) = rk(j,52)*vc(j,1,liso)*vc(j,1,lho)
756 crki(88) = rk(j,88)*vc(j,1,lho2)*dvc(j,lmo2)
757 crki(89) = rk(j,89)*vc(j,1,lho2)*dvc(j,lethp)
758 crki(90) = rk(j,90)*vc(j,1,lho2)*dvc(j,lhc3p)
759 crki(91) = rk(j,91)*vc(j,1,lho2)*dvc(j,lhc5p)
760 crki(92) = rk(j,92)*vc(j,1,lho2)*dvc(j,lhc8p)
761 crki(93) = rk(j,93)*vc(j,1,lho2)*dvc(j,lol2p)
762 crki(94) = rk(j,94)*vc(j,1,lho2)*dvc(j,loltp)
763 crki(95) = rk(j,95)*vc(j,1,lho2)*dvc(j,lolip)
764 crki(96) = rk(j,96)*vc(j,1,lho2)*dvc(j,lketp)
765 crki(97) = rk(j,97)*vc(j,1,lho2)*vc(j,1,laco3)
766 crki(98) = rk(j,98)*vc(j,1,lho2)*dvc(j,ltolp)
767 crki(99) = rk(j,99)*vc(j,1,lho2)*dvc(j,lxylp)
768 crki(100) = rk(j,100)*vc(j,1,lho2)*dvc(j,ltco3)
769 crki(101) = rk(j,101)*vc(j,1,lho2)*dvc(j,loln)
770 crki(127) = rk(j,127)*dvc(j,lxo2)*vc(j,1,lho2)
771 crki(133) = rk(j,133)*dvc(j,lxno2)*vc(j,1,lho2)
773 lossl(j,lhox) = crki(10) + 2.*crki(12) + 2.*crki(13) + crki(15) + &
774 crki(20) + crki(24) + crki(25) + crki(26) + 2.*crki(27) + &
775 crki(30) + crki(31) + .83*crki(32) + crki(33) + crki(34) + &
776 crki(35) + crki(36) + crki(37) + .75*crki(38) + .83*crki(39) + &
777 1.8*crki(40) + crki(42) + crki(43) + crki(45) + crki(46) + &
778 .5*crki(47) + .5*crki(48) + crki(49) + crki(50) + crki(51) + &
779 crki(52) + crki(88) + crki(89) + crki(90) + crki(91) + &
780 crki(92) + crki(93) + crki(94) + crki(95) + crki(96) + &
781 crki(97) + crki(98) + crki(99) + crki(100) + crki(101) + &
782 crki(127) + crki(133)
783 lossl(j,lhox) = max(alow,lossl(j,lhox))
784 prod(j,lh2o2) = crki(12) + crki(13)
785 loss(j,lh2o2) = crji(9) + crki(14)
786 END DO
788 ! if(iprt.eq.1)print *,'i1, hox,niter ',niter
789 DO j = jcs, jce
790 !********[HOX]:
791 taulinv(j) = (lossl(j,lhox)/hoxsv(j))
792 vcl(j,lhox) = prodl(j,lhox)/taulinv(j) + &
793 (hoxsv(j)-prodl(j,lhox)/taulinv(j))*exp(-dtc*taulinv(j))
794 ! vcl(j,lhox) = max(cmin(lhox),vcl(j,lhox))
795 vcl(j,lhox)=max(1.e-9,vcl(j,lhox))
796 ! if(iprt.eq.2.and.j.eq.jcs)print *,lossl(j,lhox),hoxsv(j),prodl(j,lhox),cmin(lhox),vcl(j,lhox)
797 jc(j) = jprod(j) + (rk(j,9)*vc(j,1,lno)+rk(j,8)*vc(j,1,lo3))*vcl(j &
798 ,lhox)
799 jb(j) = jbb(j) + (rk(j,27)*vcl(j,lhox))
801 !********[H2O2]:
802 tauinv(j) = loss(j,lh2o2)/h2o2sv(j)
803 vc(j,1,lh2o2) = prod(j,lh2o2)/tauinv(j) + &
804 (h2o2sv(j)-prod(j,lh2o2)/tauinv(j))*exp(-dtc*tauinv(j))
805 vc(j,1,lh2o2) = max(cmin(lh2o2),vc(j,1,lh2o2))
807 END DO
809 !********[HO]:
810 ! if(iprt.eq.1)print *,'i1, ho and ho2 '
811 DO iter = 1, inumho
812 ! if(iprt.eq.2)print *,'i1, ho and ho2 iter = ',iter
813 DO j = jcs, jce
814 vc(j,1,lho) = vc(j,1,lho) - (rk(j,27)*(vc(j,1, &
815 lho)**2)-jb(j)*vc(j,1,lho)+jc(j))/(2*rk(j,27)*vc(j,1,lho)-jb(j &
816 ))
817 vc(j,1,lho) = max(cmin(lho),vc(j,1,lho))
818 ! if(iprt.eq.2.and.j.eq.jcs)print *,vc(j,1,lho),cmin(lho),&
819 ! rk(j,27),jb(j),jc(j)
820 END DO
821 END DO
823 !********[HO2]:
824 DO j = jcs, jce
825 vc(j,1,lho2) = vcl(j,lhox) - vc(j,1,lho)
826 vc(j,1,lho2) = max(cmin(lho2),vc(j,1,lho2))
827 ! if(iprt.eq.2.and.j.eq.jcs)print *,vc(j,1,lho2),vcl(j,lhox),vc(j,1,lho),cmin(lho2)
828 END DO
830 END DO
831 !
832 ! if(iprt.eq.1)print *,'i1, aerosols '
833 if(iaerosol_sorgam==1)then
834 ! if(iprt.eq.2)print *,'i1, aerosols '
835 DO L = 1, LDROG
836 DO J = JCS, JCE
837 VDROG( J, L ) = VDROG( J , L) + PRDROG( J, L ) * DTC
838 VDROG( J, L ) = MAX( 0., VDROG( J, L ) )
839 ENDDO
840 ENDDO
841 endif
842 !
843 RETURN
844 END SUBROUTINE integ1n
845 SUBROUTINE predraten(jcs,jce,iprt,crj,crk,rj,rk,vc,dvc,vca, &
846 wlc,dvca,p,h2o,dvcg,t,r)
847 implicit none
848 ! .. Scalar Arguments ..
849 REAL, INTENT(IN) :: r
850 INTEGER, INTENT(IN) :: iprt, jce, jcs
852 ! ..
853 ! .. Array Arguments ..
854 REAL :: crj(jcs:jce,nreacj), crk(jcs:jce,nreack), &
855 dvc(jcs:jce,ldiag), dvca(jcs:jce,ldiag), &
856 dvcg(jcs:jce,ldiag), h2o(jcs:jce,1), p(jcs:jce,1), &
857 rj(jcs:jce,nreacj), rk(jcs:jce,nreack), t(jcs:jce,1), &
858 vc(jcs:jce,1,lspec), vca(jcs:jce,1,lspec), wlc(jcs:jce)
859 ! ..
860 ! .. Local Scalars ..
861 INTEGER :: i, j, k, l
862 ! ..
863 k = 1
864 i = 1
865 DO j = jcs, jce
867 crj(j,1) = rj(j,1)*vc(j,1,lno2)
868 crj(j,2) = rj(j,2)*vc(j,1,lo3)
869 crj(j,3) = rj(j,3)*vc(j,1,lo3)
870 crk(j,15) = rk(j,15)*vc(j,1,lno)*vc(j,1,lho)
871 crj(j,4) = rj(j,4)*vc(j,1,lhono)
872 crj(j,5) = rj(j,5)*vc(j,1,lhno3)
873 crk(j,10) = rk(j,10)*vc(j,1,lho2)*vc(j,1,lno2)
874 crj(j,6) = rj(j,6)*vc(j,1,lhno4)
875 crj(j,7) = rj(j,7)*vc(j,1,lno3)
876 crj(j,8) = rj(j,8)*vc(j,1,lno3)
877 crj(j,9) = rj(j,9)*vc(j,1,lh2o2)
878 crj(j,10) = rj(j,10)*vc(j,1,lhcho)
879 crj(j,11) = rj(j,11)*vc(j,1,lhcho)
880 crj(j,12) = rj(j,12)*vc(j,1,lald)
881 crj(j,13) = rj(j,13)*vc(j,1,lop1)
882 crj(j,14) = rj(j,14)*vc(j,1,lop2)
883 crj(j,15) = rj(j,15)*vc(j,1,lpaa)
884 crj(j,16) = rj(j,16)*vc(j,1,lket)
885 crj(j,17) = rj(j,17)*vc(j,1,lgly)
886 crj(j,18) = rj(j,18)*vc(j,k,lgly)
887 crj(j,19) = rj(j,19)*vc(j,k,lmgly)
888 crj(j,20) = rj(j,20)*vc(j,k,ldcb)
889 crj(j,21) = rj(j,21)*vc(j,k,lonit)
890 END DO
891 !*
892 !* end of photolysis equilimaxium species
893 DO j = jcs, jce
894 dvc(j,lo1d) = crj(j,2)/(rk(j,3)*n2+rk(j,4)*o2+rk(j,5)*h2o(j,1))
895 END DO
896 DO j = jcs, jce
897 dvc(j,lo3p) = (crj(j,1)+crj(j,3)+crj(j,8)+rk(j,3)*dvc(j,lo1d)*n2+rk( &
898 j,4)*dvc(j,lo1d)*o2)/(rk(j,2)*vc(j,1,lno2)+rk(j,1)*o2)
899 END DO
900 DO j = jcs, jce
901 crk(j,1) = rk(j,1)*dvc(j,lo3p)*o2
902 crk(j,2) = rk(j,2)*dvc(j,lo3p)*vc(j,1,lno2)
903 crk(j,3) = rk(j,3)*dvc(j,lo1d)*n2
904 crk(j,4) = rk(j,4)*dvc(j,lo1d)*o2
905 crk(j,5) = rk(j,5)*dvc(j,lo1d)*h2o(j,1)
906 crk(j,6) = rk(j,6)*vc(j,k,lo3)*vc(j,1,lno)
907 crk(j,7) = rk(j,7)*vc(j,k,lo3)*vc(j,1,lho)
908 crk(j,8) = rk(j,8)*vc(j,k,lo3)*vc(j,1,lho2)
909 crk(j,9) = rk(j,9)*vc(j,k,lho2)*vc(j,1,lno)
910 crk(j,11) = rk(j,11)*vc(j,1,lhno4)
911 crk(j,12) = rk(j,12)*vc(j,1,lho2)*vc(j,1,lho2)
912 crk(j,13) = rk(j,13)*vc(j,1,lho2)*vc(j,1,lho2)*h2o(j,1)
913 ! if(iprt.eq.1.and.j.eq.jce)then
914 ! print *,'in predate ',CRK(J, 12),CRK(J, 13),RK(J, 12),RK(J,13)
915 ! print *,VC(J,1, LHO2),H2O(J,1)
916 ! endif
917 crk(j,14) = rk(j,14)*vc(j,1,lh2o2)*vc(j,1,lho)
918 crk(j,16) = rk(j,16)*vc(j,1,lno)*vc(j,1,lno)*o2
919 crk(j,17) = rk(j,17)*vc(j,1,lo3)*vc(j,1,lno2)
920 crk(j,22) = rk(j,22)*vc(j,1,ln2o5)
921 crk(j,23) = rk(j,23)*vc(j,1,ln2o5)*h2o(j,1)
922 crk(j,24) = rk(j,24)*vc(j,1,lho)*vc(j,1,lno2)
923 crk(j,25) = rk(j,25)*vc(j,1,lho)*vc(j,1,lhno3)
924 crk(j,26) = rk(j,26)*vc(j,1,lho)*vc(j,1,lhno4)
925 crk(j,27) = rk(j,27)*vc(j,1,lho)*vc(j,1,lho2)
926 crk(j,28) = rk(j,28)*vc(j,1,lho)*vc(j,1,lso2)
927 ! if(iprt.eq.2.and.j.eq.jcs)then
928 ! print *,'in predate ',crk(j,28),rk(j,28),vc(j,1,lho),vc(j,1,lso2)
929 ! endif
930 crk(j,29) = rk(j,29)*vc(j,1,lco)*vc(j,1,lho)
931 crk(j,30) = rk(j,30)*ch4*vc(j,1,lho)
932 crk(j,31) = rk(j,31)*vc(j,1,leth)*vc(j,1,lho)
933 crk(j,32) = rk(j,32)*vc(j,1,lhc3)*vc(j,1,lho)
934 crk(j,33) = rk(j,33)*vc(j,1,lhc5)*vc(j,1,lho)
935 crk(j,34) = rk(j,34)*vc(j,1,lhc8)*vc(j,1,lho)
936 crk(j,35) = rk(j,35)*vc(j,1,lol2)*vc(j,1,lho)
937 crk(j,36) = rk(j,36)*vc(j,1,lolt)*vc(j,1,lho)
938 crk(j,37) = rk(j,37)*vc(j,1,loli)*vc(j,1,lho)
939 crk(j,38) = rk(j,38)*vc(j,1,ltol)*vc(j,1,lho)
940 crk(j,39) = rk(j,39)*vc(j,1,lxyl)*vc(j,1,lho)
941 crk(j,40) = rk(j,40)*vc(j,1,lcsl)*vc(j,1,lho)
942 crk(j,41) = rk(j,41)*vc(j,1,lhcho)*vc(j,1,lho)
943 crk(j,42) = rk(j,42)*vc(j,1,lald)*vc(j,1,lho)
944 crk(j,43) = rk(j,43)*vc(j,1,lket)*vc(j,1,lho)
945 crk(j,44) = rk(j,44)*vc(j,1,lgly)*vc(j,1,lho)
946 crk(j,45) = rk(j,45)*vc(j,1,lmgly)*vc(j,1,lho)
947 crk(j,46) = rk(j,46)*vc(j,1,ldcb)*vc(j,1,lho)
948 crk(j,47) = rk(j,47)*vc(j,1,lop1)*vc(j,1,lho)
949 crk(j,48) = rk(j,48)*vc(j,1,lop2)*vc(j,1,lho)
950 crk(j,49) = rk(j,49)*vc(j,1,lpaa)*vc(j,1,lho)
951 crk(j,50) = rk(j,50)*vc(j,1,lpan)*vc(j,1,lho)
952 crk(j,18) = rk(j,18)*vc(j,1,lno3)*vc(j,1,lno)
953 crk(j,19) = rk(j,19)*vc(j,1,lno3)*vc(j,1,lno2)
954 crk(j,20) = rk(j,20)*vc(j,1,lno3)*vc(j,1,lho2)
955 crk(j,21) = rk(j,21)*vc(j,1,lno3)*vc(j,1,lno2)
956 crk(j,51) = rk(j,51)*vc(j,1,lonit)*vc(j,1,lho)
957 crk(j,52) = rk(j,52)*vc(j,1,liso)*vc(j,1,lho)
958 crk(j,53) = rk(j,53)*vc(j,1,laco3)*vc(j,1,lno2)
959 crk(j,54) = rk(j,54)*vc(j,1,lpan)
960 crk(j,56) = rk(j,56)*vc(j,1,ltpan)
961 crk(j,78) = rk(j,78)*vc(j,1,ldcb)*vc(j,1,lno3)
962 END DO
963 DO j = jcs, jce
964 dvc(j,ltco3) = (crj(j,20)+crk(j,56)+0.9*crk(j,40)+crk(j,46)+crk(j,78 &
965 ))/(rk(j,55)*vc(j,k,lno2)+rk(j,68)*vc(j,k,lno)+rk(j,100)*vc(j,k, &
966 lho2)+rk(j,126)*vc(j,k,laco3)+rk(j,114)*dvc(j,lmo2))
967 END DO
968 DO j = jcs, jce
969 crk(j,55) = rk(j,55)*dvc(j,ltco3)*vc(j,k,lno2)
971 dvc(j,lhc3p) = (.83*crk(j,32)+0.5*crk(j,48)+crk(j,51))/ &
972 (rk(j,58)*vc(j,k,lno)+rk(j,90)*vc(j,k,lho2)+ &
973 rk(j,116)*vc(j,k,laco3)+rk(j,104)*dvc(j,lmo2))
974 crk(j,58) = rk(j,58)*dvc(j,lhc3p)*vc(j,k,lno)
975 END DO
976 DO j = jcs, jce
977 dvc(j,lhc5p) = crk(j,33)/(rk(j,60)*vc(j,k,lno)+rk(j,91)*vc(j,k,lho2) &
978 +rk(j,117)*vc(j,k,laco3)+rk(j,105)*dvc(j,lmo2))
979 crk(j,60) = rk(j,60)*dvc(j,lhc5p)*vc(j,k,lno)
980 dvc(j,lhc8p) = crk(j,34)/(rk(j,62)*vc(j,k,lno)+rk(j,92)*vc(j,k,lho2) &
981 +rk(j,118)*vc(j,k,laco3)+rk(j,106)*dvc(j,lmo2))
982 crk(j,62) = rk(j,62)*dvc(j,lhc8p)*vc(j,k,lno)
983 dvc(j,lol2p) = crk(j,35)/(rk(j,64)*vc(j,k,lno)+rk(j,93)*vc(j,k,lho2) &
984 +rk(j,119)*vc(j,k,laco3)+rk(j,107)*dvc(j,lmo2))
985 crk(j,64) = rk(j,64)*dvc(j,lol2p)*vc(j,k,lno)
986 dvc(j,loltp) = (crk(j,36)+crk(j,52))/(rk(j,65)*vc(j,k,lno)+rk(j,94)* &
987 vc(j,k,lho2)+rk(j,120)*vc(j,k,laco3)+rk(j,108)*dvc(j,lmo2))
988 crk(j,65) = rk(j,65)*dvc(j,loltp)*vc(j,k,lno)
989 dvc(j,lolip) = crk(j,37)/(rk(j,66)*vc(j,k,lno)+rk(j,95)*vc(j,k,lho2) &
990 +rk(j,121)*vc(j,k,laco3)+rk(j,109)*dvc(j,lmo2))
991 crk(j,66) = rk(j,66)*dvc(j,lolip)*vc(j,k,lno)
992 crk(j,67) = rk(j,67)*vc(j,k,laco3)*vc(j,k,lno)
993 crk(j,68) = rk(j,68)*dvc(j,ltco3)*vc(j,k,lno)
994 dvc(j,ltolp) = 0.75*crk(j,38)/(rk(j,69)*vc(j,k,lno)+rk(j,98)*vc(j,k, &
995 lho2)+rk(j,124)*vc(j,k,laco3)+rk(j,112)*dvc(j,lmo2))
996 crk(j,69) = rk(j,69)*dvc(j,ltolp)*vc(j,k,lno)
997 dvc(j,lxylp) = 0.83*crk(j,39)/(rk(j,70)*vc(j,k,lno)+rk(j,99)*vc(j,k, &
998 lho2)+rk(j,125)*vc(j,k,laco3)+rk(j,113)*dvc(j,lmo2))
1000 crk(j,70) = rk(j,70)*dvc(j,lxylp)*vc(j,k,lno)
1001 dvc(j,lethp) = (crj(j,16)+crk(j,31))/(rk(j,71)*vc(j,k,lno)+rk(j,89)* &
1002 vc(j,k,lho2)+rk(j,115)*vc(j,k,laco3)+rk(j,103)*dvc(j,lmo2))
1004 crk(j,71) = rk(j,71)*dvc(j,lethp)*vc(j,k,lno)
1005 dvc(j,lketp) = crk(j,43)/(rk(j,72)*vc(j,k,lno)+rk(j,96)*vc(j,k,lho2) &
1006 +rk(j,122)*vc(j,k,laco3)+rk(j,110)*dvc(j,lmo2))
1008 END DO
1009 DO j = jcs, jce
1010 crk(j,72) = rk(j,72)*dvc(j,lketp)*vc(j,k,lno)
1011 crk(j,80) = rk(j,80)*vc(j,k,lol2)*vc(j,k,lno3)
1012 crk(j,81) = rk(j,81)*vc(j,k,lolt)*vc(j,k,lno3)
1013 crk(j,82) = rk(j,82)*vc(j,k,loli)*vc(j,k,lno3)
1014 crk(j,83) = rk(j,83)*vc(j,k,liso)*vc(j,k,lno3)
1015 dvc(j,loln) = (crk(j,80)+crk(j,81)+crk(j,82)+crk(j,83))/ &
1016 (rk(j,73)*vc(j,k,lno)+rk(j,101)*vc(j,k,lho2)+ &
1017 rk(j,138)*dvc(j,lmo2)+rk(j,139)*vc(j,k,laco3)+ &
1018 rk(j,140)*dvc(j,loln))
1019 ! IF (iprt==1 .AND. j==jce-1) THEN
1020 ! PRINT *, dvc(j,loln), crk(j,80), crk(j,81), crk(j,82), crk(j,83), &
1021 ! rk(j,73), vc(j,1,lno), rk(j,101), vc(j,1,lho2), rk(j,138), &
1022 ! dvc(j,lmo2), rk(j,139), vc(j,1,laco3), rk(j,140), dvc(j,loln)
1023 ! PRINT *, dvc(j,loln), crk(j,80), crk(j,81), crk(j,82), crk(j,83), &
1024 ! rk(j,73), vc(j,k,lno), rk(j,101), vc(j,k,lho2), rk(j,138), &
1025 ! dvc(j,lmo2), rk(j,139), vc(j,k,laco3), rk(j,140), dvc(j,loln)
1026 ! END IF
1028 END DO
1029 DO j = jcs, jce
1030 crk(j,73) = rk(j,73)*dvc(j,loln)*vc(j,k,lno)
1031 crk(j,74) = rk(j,74)*vc(j,k,lhcho)*vc(j,k,lno3)
1032 crk(j,75) = rk(j,75)*vc(j,k,lald)*vc(j,k,lno3)
1033 crk(j,76) = rk(j,76)*vc(j,k,lgly)*vc(j,k,lno3)
1034 crk(j,77) = rk(j,77)*vc(j,k,lmgly)*vc(j,k,lno3)
1035 crk(j,79) = rk(j,79)*vc(j,k,lcsl)*vc(j,k,lno3)
1036 crk(j,84) = rk(j,84)*vc(j,k,lol2)*vc(j,k,lo3)
1037 crk(j,85) = rk(j,85)*vc(j,k,lolt)*vc(j,k,lo3)
1038 crk(j,86) = rk(j,86)*vc(j,k,loli)*vc(j,k,lo3)
1039 crk(j,87) = rk(j,87)*vc(j,k,liso)*vc(j,k,lo3)
1040 crk(j,89) = rk(j,89)*vc(j,k,lho2)*dvc(j,lethp)
1041 crk(j,90) = rk(j,90)*vc(j,k,lho2)*dvc(j,lhc3p)
1042 crk(j,91) = rk(j,91)*vc(j,k,lho2)*dvc(j,lhc5p)
1043 crk(j,92) = rk(j,92)*vc(j,k,lho2)*dvc(j,lhc8p)
1044 crk(j,93) = rk(j,93)*vc(j,k,lho2)*dvc(j,lol2p)
1045 crk(j,94) = rk(j,94)*vc(j,k,lho2)*dvc(j,loltp)
1046 crk(j,95) = rk(j,95)*vc(j,k,lho2)*dvc(j,lolip)
1047 crk(j,96) = rk(j,96)*vc(j,k,lho2)*dvc(j,lketp)
1048 crk(j,97) = rk(j,97)*vc(j,k,lho2)*vc(j,k,laco3)
1049 crk(j,98) = rk(j,98)*vc(j,k,lho2)*dvc(j,ltolp)
1050 crk(j,99) = rk(j,99)*vc(j,k,lho2)*dvc(j,lxylp)
1051 crk(j,100) = rk(j,100)*vc(j,k,lho2)*dvc(j,ltco3)
1052 crk(j,101) = rk(j,101)*vc(j,k,lho2)*dvc(j,loln)
1053 crk(j,115) = rk(j,115)*dvc(j,lethp)*vc(j,k,laco3)
1054 crk(j,116) = rk(j,116)*dvc(j,lhc3p)*vc(j,k,laco3)
1055 crk(j,117) = rk(j,117)*dvc(j,lhc5p)*vc(j,k,laco3)
1056 crk(j,118) = rk(j,118)*dvc(j,lhc8p)*vc(j,k,laco3)
1057 crk(j,119) = rk(j,119)*dvc(j,lol2p)*vc(j,k,laco3)
1058 crk(j,120) = rk(j,120)*dvc(j,loltp)*vc(j,k,laco3)
1059 crk(j,121) = rk(j,121)*dvc(j,lolip)*vc(j,k,laco3)
1060 crk(j,122) = rk(j,122)*dvc(j,lketp)*vc(j,k,laco3)
1061 crk(j,123) = rk(j,123)*vc(j,k,laco3)*vc(j,k,laco3)
1062 crk(j,124) = rk(j,124)*vc(j,k,laco3)*dvc(j,ltolp)
1063 crk(j,125) = rk(j,125)*vc(j,k,laco3)*dvc(j,lxylp)
1064 crk(j,126) = rk(j,126)*vc(j,k,laco3)*dvc(j,ltco3)
1065 END DO
1066 DO j = jcs, jce
1067 dvc(j,lxo2) = (0.25*crk(j,33)+0.75*crk(j,34)+0.9*crk(j,40)+crk(j,50) &
1068 +2.0*crk(j,68)+2.0*crk(j,126)+crk(j,114))/ &
1069 (rk(j,131)*vc(j,k,lno)+rk(j,127)*vc(j,k,lho2)+ &
1070 rk(j,128)*dvc(j,lmo2)+rk(j,129)*vc(j,k,laco3))
1072 crk(j,127) = rk(j,127)*dvc(j,lxo2)*vc(j,k,lho2)
1073 crk(j,129) = rk(j,129)*dvc(j,lxo2)*vc(j,k,laco3)
1074 crk(j,130) = 0.
1075 crk(j,131) = rk(j,131)*dvc(j,lxo2)*vc(j,k,lno)
1076 dvc(j,lxno2) = crk(j,79)/(rk(j,132)*vc(j,k,lno2)+rk(j,133)*vc(j,k, &
1077 lho2)+rk(j,134)*dvc(j,lmo2)+rk(j,135)*vc(j,k,laco3))
1078 END DO
1079 DO j = jcs, jce
1080 crk(j,132) = rk(j,132)*dvc(j,lxno2)*vc(j,k,lno2)
1081 crk(j,133) = rk(j,133)*dvc(j,lxno2)*vc(j,k,lho2)
1082 crk(j,135) = rk(j,135)*dvc(j,lxno2)*vc(j,k,laco3)
1083 crk(j,136) = 0.
1084 crk(j,137) = rk(j,137)*vc(j,k,ln2o5)
1085 crk(j,138) = rk(j,138)*dvc(j,lmo2)*dvc(j,loln)
1086 crk(j,139) = rk(j,139)*vc(j,k,laco3)*dvc(j,loln)
1087 crk(j,140) = rk(j,140)*dvc(j,loln)*dvc(j,loln)
1088 ! IF (iprt==1 .AND. j==jce-1) THEN
1089 ! PRINT *, 'crk(j,140),RK(J, 140),DVC(J,LOLN)'
1090 ! PRINT *, crk(j,140), rk(j,140), dvc(j,loln)
1091 ! END IF
1093 dvc(j,lmo2) = (crj(j,12)+crj(j,15)+crk(j,30)+.5*crk(j,47)+crk(j,67)+ &
1094 .22*crk(j,85)+.31*crk(j,86)+.22*crk(j,87)+.5*(crk(j,111)+crk(j, &
1095 115)+crk(j,116)+crk(j,117)+crk(j,118)+crk(j,119)+crk(j,120)+crk(j, &
1096 121)+crk(j,122))+2.0*crk(j,123)+crk(j,124)+crk(j,125)+crk(j,126)+ &
1097 crk(j,129)+crk(j,135)+.5*crk(j,139))/(rk(j,57)*vc(j,k,lno)+rk(j,88 &
1098 )*vc(j,k,lho2)+rk(j,102)*dvc(j,lmo2)+rk(j,103)*dvc(j,lethp)+ &
1099 rk(j,104)*dvc(j,lhc3p)+rk(j,105)*dvc(j,lhc5p)+ &
1100 rk(j,106)*dvc(j,lhc8p)+rk(j,107)*dvc(j,lol2p)+ &
1101 rk(j,108)*dvc(j,loltp)+rk(j,109)*dvc(j,lolip)+ &
1102 rk(j,110)*dvc(j,lketp)+rk(j,111)*vc(j,k,laco3)+ &
1103 rk(j,112)*dvc(j,ltolp)+rk(j,113)*dvc(j,lxylp)+ &
1104 rk(j,114)*dvc(j,ltco3)+rk(j,128)*dvc(j,lxo2)+ &
1105 rk(j,134)*dvc(j,lxno2)+rk(j,138)*dvc(j,loln))
1106 END DO
1107 DO j = jcs, jce
1108 crk(j,57) = rk(j,57)*dvc(j,lmo2)*vc(j,k,lno)
1109 crk(j,88) = rk(j,88)*vc(j,k,lho2)*dvc(j,lmo2)
1110 crk(j,102) = rk(j,102)*dvc(j,lmo2)*dvc(j,lmo2)
1111 crk(j,103) = rk(j,103)*dvc(j,lmo2)*dvc(j,lethp)
1112 crk(j,104) = rk(j,104)*dvc(j,lmo2)*dvc(j,lhc3p)
1113 crk(j,105) = rk(j,105)*dvc(j,lmo2)*dvc(j,lhc5p)
1114 crk(j,106) = rk(j,106)*dvc(j,lmo2)*dvc(j,lhc8p)
1115 crk(j,107) = rk(j,107)*dvc(j,lmo2)*dvc(j,lol2p)
1116 crk(j,108) = rk(j,108)*dvc(j,lmo2)*dvc(j,loltp)
1117 crk(j,109) = rk(j,109)*dvc(j,lmo2)*dvc(j,lolip)
1118 crk(j,110) = rk(j,110)*dvc(j,lmo2)*dvc(j,lketp)
1119 crk(j,111) = rk(j,111)*dvc(j,lmo2)*vc(j,k,laco3)
1120 crk(j,112) = rk(j,112)*dvc(j,lmo2)*dvc(j,ltolp)
1121 crk(j,113) = rk(j,113)*dvc(j,lmo2)*dvc(j,lxylp)
1122 crk(j,114) = rk(j,114)*dvc(j,lmo2)*dvc(j,ltco3)
1123 crk(j,128) = rk(j,128)*dvc(j,lxo2)*dvc(j,lmo2)
1124 crk(j,134) = rk(j,134)*dvc(j,lxno2)*dvc(j,lmo2)
1125 END DO
1126 ! if(iprt.eq.2)then
1127 ! print *,'list of crks ',vc(jcs,1,lho2),vc(jcs,1,lho)
1128 ! do l=1,140
1129 ! print *,'crk at l ',l,crk(1,l)
1130 ! enddo
1131 ! print *,'list of crjs '
1132 ! do l=1,21
1133 ! print *,'crj at l ',l,crj(1,l)
1134 ! enddo
1135 ! endif
1137 RETURN
1138 END SUBROUTINE predraten
1139 SUBROUTINE producn(jcs,jce,iprt,crj,crk,loss,prod,prodl,lossl, &
1140 prdrog,iaerosol_sorgam)
1141 !***********************************************************************
1142 !** DESCRIPTION:
1143 !** Subroutine PRODUC derives the differential equations from the
1144 !** reaction rates computed in the PREDRATE subroutine.
1145 !***********************************************************************
1146 ! .. Scalar Arguments ..
1147 implicit none
1148 INTEGER,INTENT(IN) :: iprt, jce, jcs
1149 real , intent (out) :: prdrog(jcs:jce,ldrog)
1150 integer, intent (in) :: iaerosol_sorgam
1152 ! ..
1153 ! .. Array Arguments ..
1154 REAL :: crj(jcs:jce,nreacj), crk(jcs:jce,nreack), &
1155 loss(jcs:jce,lpred), lossl(jcs:jce,lump), &
1156 prod(jcs:jce,lpred), prodl(jcs:jce,lump)
1157 ! ..
1158 ! .. Local Scalars ..
1159 REAL :: alow
1160 INTEGER :: j, l
1161 ! ..
1162 ! .. Intrinsic Functions ..
1163 INTRINSIC max
1164 ! ..
1165 alow = 1.e-17
1166 PRDROG = 0.
1167 DO l = 1, lump
1168 DO j = jcs, jce
1169 prodl(j,l) = 0.0
1170 END DO
1171 END DO
1172 DO l = 1, lpred
1173 DO j = jcs, jce
1174 prod(j,l) = 0.0
1175 END DO
1176 END DO
1178 DO j = jcs, jce
1179 prodl(j,lpao3) = crj(j,16) + crj(j,19) + .02*crj(j,20) + crk(j,42) + &
1180 crk(j,45) + crk(j,49) + .05*crk(j,68) + crk(j,75) + crk(j,77) + &
1181 .03*crk(j,114)
1182 END DO
1184 DO j = jcs, jce
1185 lossl(j,lpao3) = crk(j,67) + crk(j,97) + crk(j,111) + crk(j,115) + &
1186 crk(j,116) + crk(j,117) + crk(j,118) + crk(j,119) + crk(j,120) + &
1187 crk(j,121) + crk(j,122) + 2.*crk(j,123) + crk(j,124) + &
1188 crk(j,125) + crk(j,129) + crk(j,135) + crk(j,50) + &
1189 .95*crk(j,126) + crk(j,139)
1190 lossl(j,lpao3) = max(alow,lossl(j,lpao3))
1191 END DO
1193 DO j = jcs, jce
1194 loss(j,lo3) = crk(j,2) + crk(j,5) + crk(j,6) + crk(j,7) + crk(j,8) + &
1195 crk(j,17) + crk(j,84) + crk(j,85) + crk(j,86) + crk(j,87)
1196 END DO
1198 DO j = jcs, jce
1199 prod(j,lo3) = crj(j,1) + crj(j,8)
1200 ! PROD(J,LO3) = CRK(J,1)
1201 END DO
1203 DO j = jcs, jce
1204 prod(j,lno2) = crj(j,5) + crj(j,6) + crj(j,7) + crj(j,21) + &
1205 crk(j,6) + crk(j,9) + crk(j,11) + 2.*crk(j,16) + 2.*crk(j,18) + &
1206 crk(j,19) + crk(j,22) + crk(j,26) + crk(j,51) + crk(j,54) + &
1207 crk(j,56) + crk(j,57) + .964*crk(j,58) + .92*crk(j,60) + &
1208 .76*crk(j,62) + crk(j,64) + crk(j,65) + crk(j,66) + crk(j,67) + &
1209 crk(j,68) + crk(j,69) + crk(j,70) + crk(j,71) + crk(j,72) + &
1210 crk(j,131) + crk(j,138) + crk(j,139) + 2.0*crk(j,140)
1211 ! IF (iprt==1 .AND. j==jce) THEN
1212 !! print *,'CRK(J,58),CRK(J,60),CRK(J,62),CRK(J,64),CRK(J,66),'
1213 !! 1 ,'CRK(J,71),CRK(J,138),CRK(J,139),CRK(J,140)'
1214 ! PRINT *, 'CRK(J,58),CRK(J,60),CRK(J,62),CRK(J,64),..'
1215 ! PRINT *, crk(j,58), crk(j,60), crk(j,62), crk(j,64), crk(j,65), &
1216 ! crk(j,66), crk(j,138), crk(j,139), crk(j,140)
1217 ! END IF
1218 END DO
1220 DO j = jcs, jce
1221 loss(j,lno2) = crj(j,1) + crk(j,2) + crk(j,10) + crk(j,17) + &
1222 crk(j,19) + crk(j,21) + crk(j,24) + crk(j,53) + crk(j,55) + &
1223 crk(j,132)
1224 ! if(iprt.eq.1.and.j.eq.jce)then
1225 ! print *,LNO2,LOSS(J,LNO2),CRJ(J, 9),CRK(J, 14)
1226 ! endif
1227 END DO
1229 DO j = jcs, jce
1230 prodl(j,lnox) = crj(j,5) + crj(j,6) + crj(j,8) + crj(j,21) + &
1231 crk(j,11) + crk(j,18) + crk(j,22) + crk(j,26) + crk(j,54) + &
1232 crk(j,56) + crk(j,73) + crj(j,4) + crj(j,7) + crk(j,19) + &
1233 crk(j,51)
1234 END DO
1236 DO j = jcs, jce
1237 lossl(j,lnox) = crk(j,10) + crk(j,17) + crk(j,21) + crk(j,24) + &
1238 crk(j,53) + crk(j,55) + crk(j,132) + crk(j,15) + .036*crk(j,58) + &
1239 .08*crk(j,60) + .24*crk(j,62)
1240 END DO
1242 !....... NTOTAL=NOX+HNO3+NO3+2*N2O5+HONO+HNO4+PAN+TPAN+ONIT
1243 !....... ----> TAKE THE COMMENTED STUFF
1244 !........ UNCOMMENTED IS: NTOTAL=NOX+HNO3+NO3+2*N2O5+HONO+HNO4
1245 !-------> DON'T FORGET TO CHANGE INTEG AND CHEM
1248 ! LOSSL(J,LNTOTAL)= CRK(J,80)+CRK(J,81)+
1249 ! 1 CRK(J,82)+CRK(J,83)
1251 ! LOSSL(J,LNTOTAL)= CRK(J,53)+CRK(J,55)+
1252 ! 1 CRK(J,132)+
1253 ! 1 CRK(J,80)+CRK(J,81)+
1254 ! 1 CRK(J,82)+CRK(J,83)
1255 ! PRODL(J,LNTOTAL)=
1256 ! 1 CRK(J,20)+
1257 ! 1 CRK(J,50)+CRK(J,54)+CRK(J,56)+CRK(J,73)+
1258 ! 1 CRK(J,51)
1260 ! PRODL(J,LNTOTAL)= CRK(J,101) + CRK(J,73)
1262 DO j = jcs, jce
1263 prodl(j,ln2n3) = crk(j,17) + crk(j,25) + crk(j,50)
1264 END DO
1266 DO j = jcs, jce
1267 lossl(j,ln2n3) = crk(j,23) + crj(j,7) + crj(j,8) + crk(j,18) + &
1268 crk(j,19) + crk(j,20) + crk(j,74) + crk(j,75) + crk(j,76) + &
1269 crk(j,77) + crk(j,78) + crk(j,79) + crk(j,80) + crk(j,81) + &
1270 crk(j,82) + crk(j,83) + crk(j,137)
1271 END DO
1273 DO j = jcs, jce
1274 loss(j,lpan) = crk(j,50) + crk(j,54)
1275 END DO
1277 DO j = jcs, jce
1278 prod(j,lpan) = crk(j,53)
1279 END DO
1281 DO j = jcs, jce
1282 loss(j,lhno3) = crj(j,5) + crk(j,25)
1283 END DO
1285 DO j = jcs, jce
1286 prod(j,lhno3) = crk(j,20) + 2.D0*crk(j,23) + crk(j,24) + crk(j,74) + &
1287 crk(j,75) + crk(j,76) + crk(j,77) + crk(j,78) + crk(j,79) + &
1288 2.*crk(j,137)
1289 END DO
1291 DO j = jcs, jce
1292 loss(j,lh2o2) = max(alow,crj(j,9) + crk(j,14) )
1293 ! if(iprt.eq.1.and.j.eq.jce)then
1294 ! print *,LH2O2,LOSS(J,LH2O2),CRJ(J, 9),CRK(J, 14)
1295 ! endif
1296 END DO
1298 DO j = jcs, jce
1299 prod(j,lh2o2) = crk(j,12) + crk(j,13)
1300 ! if(iprt.eq.1.and.j.eq.jce)then
1301 ! print *,LH2O2,prod(J,LH2O2),CRK(J, 12),CRK(J, 13)
1302 ! endif
1303 END DO
1305 DO j = jcs, jce
1306 loss(j,lhcho) = crj(j,10) + crj(j,11) + crk(j,41) + crk(j,74)
1307 END DO
1309 DO j = jcs, jce
1310 prod(j,lhcho) = crj(j,13) + .13*crj(j,17) + .45*crj(j,18) + &
1311 .009*crk(j,32) + .5*crk(j,47) + crk(j,50) + crk(j,57) + &
1312 .09*crk(j,58) + .04*crk(j,62) + 1.6*crk(j,64) + crk(j,65) + &
1313 .28*crk(j,66) + crk(j,73) + crk(j,84) + .53*crk(j,85) + &
1314 .18*crk(j,86) + .53*crk(j,87) + 1.5*crk(j,102) + .75*crk(j,103) + &
1315 .75*crk(j,104) + .77*crk(j,105) + .80*crk(j,106) + &
1316 1.55*crk(j,107) + 1.25*crk(j,108) + .89*crk(j,109) + &
1317 .75*crk(j,110) + crk(j,111) + crk(j,112) + crk(j,113) + &
1318 .5*crk(j,114) + .8*crk(j,119) + .5*crk(j,120) + .14*crk(j,121) + &
1319 crk(j,128) + crk(j,134) + 1.75*crk(j,138) + crk(j,139) + &
1320 2.0*crk(j,140)
1321 END DO
1323 DO j = jcs, jce
1324 prod(j,lhono) = crk(j,15)
1325 END DO
1327 DO j = jcs, jce
1328 loss(j,lhono) = crj(j,4)
1329 END DO
1331 DO j = jcs, jce
1332 prod(j,lhno4) = crk(j,10)
1333 END DO
1335 DO j = jcs, jce
1336 loss(j,lhno4) = crj(j,6) + crk(j,11) + crk(j,26)
1337 END DO
1339 DO j = jcs, jce
1340 prod(j,ln2o5) = crk(j,21)
1341 END DO
1343 DO j = jcs, jce
1344 loss(j,ln2o5) = crk(j,22) + crk(j,23) + crk(j,137)
1345 END DO
1347 DO j = jcs, jce
1348 prod(j,lno3) = crk(j,17) + crk(j,22) + crk(j,25) + crk(j,50)
1349 END DO
1351 DO j = jcs, jce
1352 loss(j,lno3) = crj(j,7) + crj(j,8) + crk(j,18) + crk(j,19) + &
1353 crk(j,20) + crk(j,21) + crk(j,74) + crk(j,75) + crk(j,76) + &
1354 crk(j,77) + crk(j,78) + crk(j,79) + crk(j,80) + crk(j,81) + &
1355 crk(j,82) + crk(j,83)
1356 END DO
1358 DO j = jcs, jce
1359 loss(j,lco) = crk(j,29)
1360 END DO
1362 DO j = jcs, jce
1363 prod(j,lco) = crj(j,10) + crj(j,11) + crj(j,12) + 1.87*crj(j,17) + &
1364 1.55*crj(j,18) + crj(j,19) + crk(j,41) + 2.*crk(j,44) + &
1365 crk(j,45) + .95*crk(j,68) + crk(j,74) + 2.*crk(j,76) + crk(j,77) + &
1366 .42*crk(j,84) + .33*crk(j,85) + .23*crk(j,86) + .33*crk(j,87) + &
1367 .475*crk(j,114) + .95*crk(j,126)
1368 END DO
1370 DO j = jcs, jce
1371 loss(j,lald) = crj(j,12) + crk(j,42) + crk(j,75)
1372 END DO
1374 DO j = jcs, jce
1375 prod(j,lald) = crj(j,14) + .075*crk(j,32) + .2*crj(j,21) + &
1376 .5*crk(j,48) + .75*crk(j,58) + .38*crk(j,60) + .35*crk(j,62) + &
1377 .2*crk(j,64) + crk(j,65) + 1.45*crk(j,66) + crk(j,73) + &
1378 crk(j,71) + .5*crk(j,85) + .72*crk(j,86) + .5*crk(j,87) + &
1379 .75*crk(j,103) + .15*crk(j,104) + .41*crk(j,105) + &
1380 .46*crk(j,106) + .35*crk(j,107) + .75*crk(j,108) + &
1381 .725*crk(j,109) + crk(j,115) + .2*crk(j,116) + .14*crk(j,117) + &
1382 .1*crk(j,118) + .6*crk(j,119) + crk(j,120) + .725*crk(j,121) + &
1383 crk(j,138) + crk(j,139) + 2.0*crk(j,140)
1384 END DO
1386 DO j = jcs, jce
1387 loss(j,lop1) = crj(j,13) + crk(j,47)
1388 END DO
1390 DO j = jcs, jce
1391 prod(j,lop1) = crk(j,88)
1392 END DO
1394 DO j = jcs, jce
1395 loss(j,lop2) = crj(j,14) + crk(j,48)
1396 END DO
1398 DO j = jcs, jce
1399 prod(j,lop2) = crk(j,89) + crk(j,90) + crk(j,91) + crk(j,92) + &
1400 crk(j,93) + crk(j,94) + crk(j,95) + crk(j,96) + crk(j,98) + &
1401 crk(j,99) + crk(j,100) + crk(j,127) + crk(j,133)
1402 END DO
1404 DO j = jcs, jce
1405 loss(j,lpaa) = crj(j,15) + crk(j,49)
1406 END DO
1408 DO j = jcs, jce
1409 prod(j,lpaa) = crk(j,97)
1410 END DO
1412 DO j = jcs, jce
1413 loss(j,lket) = crj(j,16) + crk(j,43)
1414 END DO
1416 DO j = jcs, jce
1417 prod(j,lket) = .8*crj(j,21) + .025*crk(j,32) + .25*crk(j,58) + &
1418 .69*crk(j,60) + 1.06*crk(j,62) + .10*crk(j,66) + .10*crk(j,86) + &
1419 .6*crk(j,104) + .75*crk(j,105) + 1.39*crk(j,106) + &
1420 .55*crk(j,109) + .8*crk(j,116) + .86*crk(j,117) + .9*crk(j,118) + &
1421 .55*crk(j,121)
1422 END DO
1424 DO j = jcs, jce
1425 loss(j,lgly) = crj(j,17) + crj(j,18) + crk(j,44) + crk(j,76)
1426 END DO
1428 DO j = jcs, jce
1429 prod(j,lgly) = .89*crk(j,68) + .16*crk(j,69) + .16*crk(j,112) + &
1430 .44*crk(j,114) + .2*crk(j,124) + .89*crk(j,126)
1431 END DO
1433 DO j = jcs, jce
1434 loss(j,lmgly) = crj(j,19) + crk(j,45) + crk(j,77)
1435 END DO
1437 DO j = jcs, jce
1438 prod(j,lmgly) = .11*crk(j,68) + .17*crk(j,69) + .450*crk(j,70) + &
1439 crk(j,72) + .75*crk(j,110) + .17*crk(j,112) + .45*crk(j,113) + &
1440 .05*crk(j,114) + crk(j,122) + .8*crk(j,124) + crk(j,125) + &
1441 .11*crk(j,126)
1442 END DO
1444 DO j = jcs, jce
1445 loss(j,ldcb) = crj(j,20) + crk(j,46) + crk(j,78)
1446 END DO
1448 DO j = jcs, jce
1449 loss(j,ldcb) = max(alow,loss(j,ldcb))
1450 END DO
1452 DO j = jcs, jce
1453 prod(j,ldcb) = .70*crk(j,69) + .806*crk(j,70) + .7*crk(j,112) + &
1454 .806*crk(j,113) + crk(j,124) + crk(j,125)
1455 END DO
1457 DO j = jcs, jce
1458 loss(j,lonit) = crj(j,21) + crk(j,51)
1459 END DO
1461 DO j = jcs, jce
1462 prod(j,lonit) = .036*crk(j,58) + .08*crk(j,60) + .24*crk(j,62) + &
1463 crk(j,101) + crk(j,132)
1464 END DO
1466 DO j = jcs, jce
1467 loss(j,lso2) = crk(j,28)
1468 END DO
1470 DO j = jcs, jce
1471 loss(j,lsulf) = 0.
1472 END DO
1474 DO j = jcs, jce
1475 prod(j,lsulf) = crk(j,28)
1476 ! if(iprt==2)print *,' j,prod = ',j,prod(j,lsulf)
1477 END DO
1479 DO j = jcs, jce
1480 loss(j,leth) = crk(j,31)
1481 END DO
1483 DO j = jcs, jce
1484 loss(j,lhc3) = crk(j,32)
1485 END DO
1487 DO j = jcs, jce
1488 loss(j,lhc5) = crk(j,33)
1489 END DO
1491 DO j = jcs, jce
1492 loss(j,lhc8) = crk(j,34)
1493 END DO
1495 DO j = jcs, jce
1496 loss(j,lol2) = crk(j,35) + crk(j,80) + crk(j,84)
1497 END DO
1499 DO j = jcs, jce
1500 loss(j,lolt) = crk(j,36) + crk(j,81) + crk(j,85)
1501 END DO
1503 DO j = jcs, jce
1504 loss(j,loli) = crk(j,37) + crk(j,82) + crk(j,86)
1505 END DO
1507 DO j = jcs, jce
1508 loss(j,ltol) = crk(j,38)
1509 END DO
1511 DO j = jcs, jce
1512 loss(j,lcsl) = crk(j,40) + .5*crk(j,79)
1513 END DO
1515 DO j = jcs, jce
1516 prod(j,lcsl) = .25*crk(j,38) + .17*crk(j,39)
1517 END DO
1519 DO j = jcs, jce
1520 loss(j,lxyl) = crk(j,39)
1521 END DO
1523 DO j = jcs, jce
1524 loss(j,laco3) = crk(j,53) + crk(j,67) + crk(j,97) + crk(j,111) + &
1525 crk(j,115) + crk(j,116) + crk(j,117) + crk(j,118) + crk(j,119) + &
1526 crk(j,120) + crk(j,121) + crk(j,122) + 2.*crk(j,123) + &
1527 crk(j,124) + crk(j,125) + .95*crk(j,126) + crk(j,129) + &
1528 crk(j,135) + crk(j,139)
1529 END DO
1531 DO j = jcs, jce
1532 prod(j,laco3) = crj(j,16) + crj(j,19) + .02*crj(j,20) + crk(j,42) + &
1533 crk(j,45) + crk(j,49) + crk(j,54) + .05*crk(j,68) + crk(j,75) + &
1534 crk(j,77) + .03*crk(j,114)
1535 END DO
1537 DO j = jcs, jce
1538 loss(j,liso) = crk(j,52) + crk(j,83) + crk(j,87)
1539 END DO
1541 DO j = jcs, jce
1542 loss(j,ltpan) = crk(j,56)
1543 END DO
1545 DO j = jcs, jce
1546 prod(j,ltpan) = crk(j,55)
1547 END DO
1549 DO j = jcs, jce
1550 loss(j,lora1) = 1.E-27
1551 END DO
1553 DO j = jcs, jce
1554 prod(j,lora1) = .4*crk(j,84) + .06*crk(j,86) + .2*crk(j,85) + &
1555 .2*crk(j,87)
1556 END DO
1558 DO j = jcs, jce
1559 loss(j,lora2) = 1.E-27
1560 END DO
1562 DO j = jcs, jce
1563 prod(j,lora2) = .2*crk(j,85) + .29*crk(j,86) + .2*crk(j,87) + &
1564 .5*crk(j,111) + .5*crk(j,114) + .5*crk(j,115) + .5*crk(j,116) + &
1565 .5*crk(j,117) + .5*crk(j,118) + .5*crk(j,119) + .5*crk(j,120) + &
1566 .5*crk(j,121) + .5*crk(j,122) + .5*crk(j,139)
1567 END DO
1569 DO j = jcs, jce
1570 lossl(j,lhox) = crk(j,15) + crk(j,24) + crk(j,25) + crk(j,26) + &
1571 crk(j,27) + crk(j,30) + crk(j,31) + .83*crk(j,32) + crk(j,33) + &
1572 crk(j,34) + crk(j,35) + crk(j,36) + crk(j,37) + .75*crk(j,38) + &
1573 .83*crk(j,39) + 1.8*crk(j,40) + crk(j,42) + crk(j,43) + &
1574 crk(j,45) + crk(j,46) + crk(j,49) + crk(j,50) + crk(j,51) + &
1575 crk(j,52) + crk(j,10) + 2.*crk(j,12) + 2.*crk(j,13) + crk(j,20) + &
1576 crk(j,27) + crk(j,88) + crk(j,89) + crk(j,90) + crk(j,91) + &
1577 .5*crk(j,47) + .5*crk(j,48) + crk(j,92) + crk(j,93) + crk(j,94) + &
1578 crk(j,95) + crk(j,96) + crk(j,97) + crk(j,98) + crk(j,99) + &
1579 crk(j,100) + crk(j,101) + crk(j,127) + crk(j,133)
1580 lossl(j,lhox) = max(alow,lossl(j,lhox))
1581 END DO
1583 DO j = jcs, jce
1584 prodl(j,lhox) = crj(j,4) + crj(j,5) + crj(j,6) + 2.*crj(j,9) + &
1585 crj(j,13) + crj(j,14) + crj(j,15) + 2.*crk(j,5) + 2.*crj(j,11) + &
1586 crj(j,12) + crj(j,13) + crj(j,14) + .8*crj(j,18) + crj(j,19) + &
1587 .98*crj(j,20) + crj(j,21) + crk(j,11) + crk(j,57) + &
1588 .964*crk(j,58) + .92*crk(j,60) + .76*crk(j,62) + crk(j,64) + &
1589 crk(j,65) + crk(j,66) + .92*crk(j,68) + crk(j,69) + crk(j,70) + &
1590 crk(j,71) + crk(j,72) + crk(j,74) + crk(j,76) + .12*crk(j,84) + &
1591 .33*crk(j,85) + .40*crk(j,86) + .33*crk(j,87) + crk(j,102) + &
1592 crk(j,103) + crk(j,104) + crk(j,105) + crk(j,106) + crk(j,107) + &
1593 crk(j,108) + crk(j,109) + crk(j,110) + .5*crk(j,111) + &
1594 2.0*crk(j,112) + 2.*crk(j,113) + .46*crk(j,114) + .5*crk(j,115) + &
1595 .5*crk(j,116) + .5*crk(j,117) + .5*crk(j,118) + .5*crk(j,119) + &
1596 .5*crk(j,120) + .5*crk(j,121) + .5*crk(j,122) + crk(j,124) + &
1597 crk(j,125) + .92*crk(j,126) + crk(j,128) + crk(j,134) + &
1598 .5*crk(j,138)
1599 END DO
1601 ! DO 850 L=1 ,LPRED
1602 ! DO 850 J=JCS,JCE
1603 ! PROD(J,L)= PROD(J,L) + PRODS(J,L)
1604 !850 CONTINUE
1606 ! DO 900 J=JCS,JCE
1607 ! PRODL(J,LNOX) = PRODL(J,LNOX)+ PRODS(J,LNO) + PRODS(J,LNO2)
1608 ! PRODL(J,LNTOTAL)= PRODL(J,LNTOTAL)+PRODS(J,LNO)+PRODS(J,LNO2)
1609 !900 CONTINUE
1610 DO J = JCS, JCE
1611 PRDROG(J,PXYL) = CRK(J, 39)
1612 PRDROG(J,PTOL) = CRK(J, 38)
1613 PRDROG(J,PCSL1) = CRK(J, 40)
1614 PRDROG(J,PCSL2) = 0.50 * CRK(J, 79)
1615 PRDROG(J,PHC8) = CRK(J, 34)
1616 PRDROG(J,POLI1) = CRK(J, 37)
1617 PRDROG(J,POLI2) = CRK(J, 82)
1618 PRDROG(J,POLI3) = CRK(J, 86)
1619 PRDROG(J,POLT1) = CRK(J, 36)
1620 PRDROG(J,POLT2) = CRK(J, 81)
1621 PRDROG(J,POLT3) = CRK(J, 85)
1622 !
1623 ! next lines for radm only, RACM would be different
1624 !
1625 PRDROG(J,PAPI1) = 0.
1626 PRDROG(J,PAPI2) = 0.
1627 PRDROG(J,PAPI3) = 0.
1628 PRDROG(J,PLIM1) = 0.
1629 PRDROG(J,PLIM2) = 0.
1630 PRDROG(J,PLIM3) = 0.
1631 ENDDO
1633 RETURN
1634 END SUBROUTINE producn
1635 SUBROUTINE setdtc(jcs,jce,dtc,dtcmax,dtcmin,dt60,prod,loss,vc, &
1636 timenow )
1637 implicit none
1638 REAL, PARAMETER :: huge=1.e10
1639 ! .. Scalar Arguments ..
1640 REAL, intent(in) :: dt60, dtcmax, dtcmin, timenow
1641 INTEGER, intent(in) :: jce, jcs
1642 REAL, intent(in) :: loss(jcs:jce,lpred), &
1643 prod(jcs:jce,lpred), vc(jcs:jce,1,lspec)
1644 real, intent(inout) :: dtc
1646 ! ..
1647 ! ..
1648 ! .. Local Scalars ..
1649 INTEGER :: j, k, l
1650 ! ..
1651 ! .. Local Arrays ..
1652 REAL :: dtlsp(lspec), dum(jcs:jce)
1653 ! ..
1654 ! .. Intrinsic Functions ..
1655 INTRINSIC abs, max, min
1656 ! ..
1657 ! ..
1658 k = 1
1660 DO l = 1, lspec
1661 dtlsp(l) = huge
1662 END DO
1663 DO l = 1, lpred
1664 IF (qdtc(l)==1) THEN
1665 DO j = jcs, jce
1666 dum(j) = prod(j,l) - loss(j,l)
1667 ! dum(j) = max(abs(dum(j)),epsilc)
1668 dum(j) = max(abs(dum(j)),1.e-30)
1669 dum(j) = .02*vc(j,1,l)/dum(j)
1670 ! DUM(J) = CVMGP(DUM(J),HUGE,VC(J,1,L)-epsilc*100.)
1671 IF (vc(j,1,l)-1.e-10>=0.) THEN
1672 dum(j) = dum(j)
1673 ELSE
1674 dum(j) = huge
1675 END IF
1677 END DO
1678 DO j = jcs, jce
1679 dtlsp(l) = min(dtlsp(l),dum(j))
1680 END DO
1681 END IF
1682 END DO
1683 ! IF (dtc<=dtcmax*.9) THEN
1684 ! dtc = dtc*1.1
1685 ! ELSE
1686 ! dtc = dtcmax
1687 ! END IF
1688 dtc = dtcmax
1689 DO l = 1, lpred
1690 IF (qdtc(l)==1) THEN
1691 IF (dtlsp(l)<dtc) THEN
1692 dtc = dtlsp(l)
1693 END IF
1694 END IF
1695 END DO
1696 IF (dtc<dtcmin) THEN
1697 dtc = dtcmin
1698 END IF
1699 IF ((timenow+dtc)>dt60) dtc = dt60 - timenow
1700 RETURN
1701 END SUBROUTINE setdtc
1702 SUBROUTINE chemin
1703 implicit none
1704 ! .. Scalar Arguments ..
1705 RETURN
1706 END SUBROUTINE chemin
1708 END MODULE module_radm