module_mosaic_cloudchem.F
References to this file elsewhere.
1 !**********************************************************************************
2 ! This computer software was prepared by Battelle Memorial Institute, hereinafter
3 ! the Contractor, under Contract No. DE-AC05-76RL0 1830 with the Department of
4 ! Energy (DOE). NEITHER THE GOVERNMENT NOR THE CONTRACTOR MAKES ANY WARRANTY,
5 ! EXPRESS OR IMPLIED, OR ASSUMES ANY LIABILITY FOR THE USE OF THIS SOFTWARE.
6 !
7 ! MOSAIC module: see module_mosaic_driver.F for information and terms of use
8 !**********************************************************************************
9
10 module module_mosaic_cloudchem
11
12
13
14 integer, parameter :: l_so4_aqyy = 1
15 integer, parameter :: l_no3_aqyy = 2
16 integer, parameter :: l_cl_aqyy = 3
17 integer, parameter :: l_nh4_aqyy = 4
18 integer, parameter :: l_na_aqyy = 5
19 integer, parameter :: l_oin_aqyy = 6
20 integer, parameter :: l_bc_aqyy = 7
21 integer, parameter :: l_oc_aqyy = 8
22
23 integer, parameter :: nyyy = 8
24
25
26
27 contains
28
29
30
31 !-----------------------------------------------------------------------
32 subroutine mosaic_cloudchem_driver( &
33 id, ktau, ktauc, dtstepc, config_flags, &
34 p_phy, t_phy, rho_phy, alt, &
35 cldfra, ph_no2, &
36 moist, chem, &
37 gas_aqfrac, numgas_aqfrac, &
38 ids,ide, jds,jde, kds,kde, &
39 ims,ime, jms,jme, kms,kme, &
40 its,ite, jts,jte, kts,kte )
41
42 use module_state_description, only: &
43 num_moist, num_chem, p_qc
44
45 use module_configure, only: grid_config_rec_type
46
47 use module_data_mosaic_asect, only: cw_phase, nphase_aer
48
49 use module_data_mosaic_other, only: k_pegbegin, name
50
51 use module_mosaic_driver, only: mapaer_tofrom_host
52
53
54 implicit none
55
56 ! subr arguments
57 integer, intent(in) :: &
58 id, ktau, ktauc, &
59 numgas_aqfrac, &
60 ids, ide, jds, jde, kds, kde, &
61 ims, ime, jms, jme, kms, kme, &
62 its, ite, jts, jte, kts, kte
63 ! id - domain index
64 ! ktau - time step number
65 ! ktauc - gas and aerosol chemistry time step number
66 ! numgas_aqfrac - last dimension of gas_aqfrac
67
68 ! [ids:ide, kds:kde, jds:jde] - spatial (x,z,y) indices for 'domain'
69 ! [ims:ime, kms:kme, jms:jme] - spatial (x,z,y) indices for 'memory'
70 ! Most arrays that are arguments to chem_driver
71 ! are dimensioned with these spatial indices.
72 ! [its:ite, kts:kte, jts:jte] - spatial (x,z,y) indices for 'tile'
73 ! chem_driver and routines under it do calculations
74 ! over these spatial indices.
75
76 type(grid_config_rec_type), intent(in) :: config_flags
77 ! config_flags - configuration and control parameters
78
79 real, intent(in) :: &
80 dtstepc
81 ! dtstepc - time step for gas and aerosol chemistry(s)
82
83 real, intent(in), &
84 dimension( ims:ime, kms:kme, jms:jme ) :: &
85 p_phy, t_phy, rho_phy, alt, cldfra, ph_no2
86 ! p_phy - air pressure (Pa)
87 ! t_phy - temperature (K)
88 ! rho_phy - moist air density (kg/m^3)
89 ! alt - dry air specific volume (m^3/kg)
90 ! cldfra - cloud fractional area (0-1)
91 ! ph_no2 - no2 photolysis rate (1/min)
92
93 real, intent(in), &
94 dimension( ims:ime, kms:kme, jms:jme, 1:num_moist ) :: &
95 moist
96 ! moist - mixing ratios of moisture species (water vapor,
97 ! cloud water, ...) (kg/kg for mass species, #/kg for number species)
98
99 real, intent(inout), &
100 dimension( ims:ime, kms:kme, jms:jme, 1:num_chem ) :: &
101 chem
102 ! chem - mixing ratios of trace gas and aerosol species (ppm for gases,
103 ! ug/kg for aerosol mass species, #/kg for aerosol number species)
104
105 real, intent(inout), &
106 dimension( ims:ime, kms:kme, jms:jme, numgas_aqfrac ) :: &
107 gas_aqfrac
108 ! gas_aqfrac - fraction (0-1) of gas that is dissolved in cloud water
109
110
111 ! local variables
112 integer :: it, jt, kt, kpeg, k_pegshift, l, mpeg
113 integer :: icase
114 integer :: igaschem_onoff, iphotol_onoff, iradical_onoff
115
116 real :: gas_aqfrac_box(numgas_aqfrac)
117 real :: ph_aq_box
118 real, parameter :: qcldwtr_cutoff = 1.0e-6
119 real :: qcldwtr
120
121
122 ! check that cw_phase is active
123 if ((cw_phase .le. 0) .or. (cw_phase .gt. nphase_aer)) then
124 print *, '*** mosaic_cloudchem_driver - cw_phase not active'
125 return
126 end if
127
128 print 93010, 'entering mosaic_cloudchem_driver - ktau =', ktau
129
130 icase = 0
131
132 ! iphotol_onoff = 1 if photolysis rate calcs are on; 0 if off
133 iphotol_onoff = 0
134 if (config_flags%phot_opt .gt. 0) iphotol_onoff = 1
135 ! igaschem_onoff = 1 if gas-phase chemistry is on; 0 if off
136 igaschem_onoff = 0
137 if (config_flags%gaschem_onoff .gt. 0) igaschem_onoff = 1
138
139 ! iradical_onoff turns aqueous radical chemistry on/off
140 ! set iradical_onoff=0 if either photolysis or gas-phase chem are off
141 if ((igaschem_onoff .le. 0) .or. (iphotol_onoff .le. 0)) then
142 iradical_onoff = 0
143 else
144 iradical_onoff = 1
145 end if
146 ! following line turns aqueous radical chem off unconditionally
147 iradical_onoff = 0
148
149
150 do 3920 jt = jts, jte
151 do 3910 it = its, ite
152
153 do 3800 kt = kts, kte-1
154
155 qcldwtr = moist(it,kt,jt,p_qc)
156 if (qcldwtr .le. qcldwtr_cutoff) goto 3800
157
158
159 k_pegshift = k_pegbegin - kts
160 kpeg = kt + k_pegshift
161 mpeg = 1
162 icase = icase + 1
163
164 ! detailed dump for debugging
165 if (ktau .eq. -13579) then
166 ! if ((ktau .eq. 30) .and. (it .eq. 23) .and. &
167 ! (jt .eq. 1) .and. (kt .eq. 11)) then
168 call mosaic_cloudchem_dumpaa( &
169 id, ktau, ktauc, dtstepc, config_flags, &
170 p_phy, t_phy, rho_phy, alt, &
171 cldfra, ph_no2, &
172 moist, chem, &
173 gas_aqfrac, numgas_aqfrac, &
174 ids,ide, jds,jde, kds,kde, &
175 ims,ime, jms,jme, kms,kme, &
176 its,ite, jts,jte, kts,kte, &
177 qcldwtr_cutoff, &
178 it, jt, kt )
179 end if
180
181 ! map from wrf-chem 3d arrays to pegasus clm & sub arrays
182 call mapaer_tofrom_host( 0, &
183 ims,ime, jms,jme, kms,kme, &
184 its,ite, jts,jte, kts,kte, &
185 it, jt, kt, kt, &
186 num_moist, num_chem, moist, chem, &
187 t_phy, p_phy, rho_phy )
188
189 ! make call '1box' cloudchem routine
190 ! print 93010, 'calling mosaic_cloudchem_1 at ijk =', it, jt, kt
191 call mosaic_cloudchem_1box( &
192 id, ktau, ktauc, dtstepc, &
193 iphotol_onoff, iradical_onoff, &
194 ph_no2(it,kt,jt), &
195 ph_aq_box, gas_aqfrac_box, &
196 numgas_aqfrac, it, jt, kt, kpeg, mpeg, icase )
197
198 ! map back to wrf-chem 3d arrays
199 call mapaer_tofrom_host( 1, &
200 ims,ime, jms,jme, kms,kme, &
201 its,ite, jts,jte, kts,kte, &
202 it, jt, kt, kt, &
203 num_moist, num_chem, moist, chem, &
204 t_phy, p_phy, rho_phy )
205
206 gas_aqfrac(it,kt,jt,:) = gas_aqfrac_box(:)
207
208
209 3800 continue
210
211 3910 continue
212 3920 continue
213
214 print 93010, 'leaving mosaic_cloudchem_driver - ktau =', ktau, icase
215 93010 format( a, 8(1x,i6) )
216
217 return
218 end subroutine mosaic_cloudchem_driver
219
220
221
222 !-----------------------------------------------------------------------
223 subroutine mosaic_cloudchem_1box( &
224 id, ktau, ktauc, dtstepc, &
225 iphotol_onoff, iradical_onoff, &
226 photol_no2_box, &
227 ph_aq_box, gas_aqfrac_box, &
228 numgas_aqfrac, it, jt, kt, kpeg, mpeg, icase )
229
230 use module_state_description, only: &
231 num_moist, num_chem
232
233 use module_data_mosaic_asect, only: &
234 msectional, &
235 maxd_asize, maxd_atype, &
236 cw_phase, nsize_aer, ntype_aer, &
237 lptr_so4_aer, lptr_no3_aer, lptr_cl_aer, lptr_co3_aer, &
238 lptr_msa_aer, lptr_nh4_aer, lptr_na_aer, lptr_ca_aer, &
239 lptr_oin_aer, lptr_bc_aer, lptr_oc_aer
240
241 use module_data_mosaic_other, only: &
242 l2maxd, ltot2, rsub
243
244 use module_data_cmu_bulkaqchem, only: &
245 meqn1max
246
247
248 implicit none
249
250 ! subr arguments
251 integer, intent(in) :: &
252 id, ktau, ktauc, &
253 numgas_aqfrac, it, jt, kt, kpeg, mpeg, &
254 icase, iphotol_onoff, iradical_onoff
255
256 real, intent(in) :: &
257 dtstepc, photol_no2_box
258
259 real, intent(inout) :: ph_aq_box
260
261 real, intent(inout), dimension( numgas_aqfrac ) :: gas_aqfrac_box
262
263 ! local variables
264 integer :: iphase
265 integer :: icase_in, idecomp_hmsa_hso5, &
266 iradical_in, istat_aqop
267
268 integer :: lptr_yyy_cwaer(maxd_asize,maxd_atype,nyyy)
269
270 real :: co2_mixrat_in
271 real :: ph_cmuaq_cur
272 real :: photol_no2_in
273
274 real :: yaq_beg(meqn1max), yaq_end(meqn1max)
275 real :: rbox(l2maxd), rbox_sv1(l2maxd)
276 real :: rbulk_cwaer(nyyy,2)
277
278 real, dimension( maxd_asize, maxd_atype ) :: fr_partit_cw
279
280
281 !
282 ! set the lptr_yyy_cwaer
283 !
284 iphase = cw_phase
285 lptr_yyy_cwaer(:,:,l_so4_aqyy) = lptr_so4_aer(:,:,iphase)
286 lptr_yyy_cwaer(:,:,l_no3_aqyy) = lptr_no3_aer(:,:,iphase)
287 lptr_yyy_cwaer(:,:,l_cl_aqyy ) = lptr_cl_aer( :,:,iphase)
288 lptr_yyy_cwaer(:,:,l_nh4_aqyy) = lptr_nh4_aer(:,:,iphase)
289 lptr_yyy_cwaer(:,:,l_na_aqyy ) = lptr_na_aer( :,:,iphase)
290 lptr_yyy_cwaer(:,:,l_oin_aqyy) = lptr_oin_aer(:,:,iphase)
291 lptr_yyy_cwaer(:,:,l_bc_aqyy ) = lptr_bc_aer( :,:,iphase)
292 lptr_yyy_cwaer(:,:,l_oc_aqyy ) = lptr_oc_aer( :,:,iphase)
293
294 !
295 ! xfer from rsub to rbox
296 !
297 rbox(1:ltot2) = max( 0.0, rsub(1:ltot2,kpeg,mpeg) )
298 rbox_sv1(1:ltot2) = rbox(1:ltot2)
299
300 !
301 !
302 ! do bulk cloud-water chemistry
303 !
304 !
305 icase_in = icase
306 iradical_in = 1
307 idecomp_hmsa_hso5 = 1
308
309 co2_mixrat_in = 350.0
310
311 photol_no2_in = photol_no2_box
312
313 ! turn off aqueous phase photolytic and radical chemistry
314 ! if either of the iphotol_onoff and iradical_onoff flags are 0
315 if ((iphotol_onoff .le. 0) .or. (iradical_onoff .le. 0)) then
316 photol_no2_in = 0.0
317 iradical_in = 0
318 end if
319
320 #if defined ( ccboxtest_box_testing_active)
321 ! following is for off-line box testing only
322 call ccboxtest_extra_args_aa( 'get', &
323 co2_mixrat_in, iradical_in, &
324 idecomp_hmsa_hso5, icase_in )
325 #endif
326
327 gas_aqfrac_box(:) = 0.0
328
329
330 ! make call to interface_to_aqoperator1
331 call interface_to_aqoperator1( &
332 istat_aqop, &
333 dtstepc, &
334 rbox, gas_aqfrac_box, &
335 rbulk_cwaer, lptr_yyy_cwaer, &
336 co2_mixrat_in, photol_no2_in, &
337 iradical_in, idecomp_hmsa_hso5, &
338 yaq_beg, yaq_end, ph_cmuaq_cur, &
339 numgas_aqfrac, id, it, jt, kt, kpeg, mpeg, ktau, icase_in )
340
341 ph_aq_box = ph_cmuaq_cur
342
343
344 #if defined ( ccboxtest_box_testing_active)
345 ! following is for off-line box testing only
346 call ccboxtest_extra_args_bb( 'put', &
347 yaq_beg, yaq_end, ph_cmuaq_cur )
348 #endif
349
350
351 !
352 !
353 ! calculate fraction of cloud-water associated with each activated aerosol bin
354 !
355 !
356 call partition_cldwtr( &
357 rbox, fr_partit_cw, &
358 it, jt, kt, kpeg, mpeg, icase_in )
359
360 !
361 !
362 ! distribute changes in bulk cloud-water composition among size bins
363 !
364 !
365 call distribute_bulk_changes( &
366 rbox, rbox_sv1, fr_partit_cw, &
367 rbulk_cwaer, lptr_yyy_cwaer, &
368 it, jt, kt, kpeg, mpeg, icase_in )
369
370
371 !
372 ! xfer back to rsub
373 !
374 rsub(1:ltot2,kpeg,mpeg) = max( 0.0, rbox(1:ltot2) )
375
376
377 !
378 ! do move-sections
379 !
380 if (msectional .lt. 1000000000) then
381 call cloudchem_apply_move_sections( &
382 rbox, rbox_sv1, &
383 it, jt, kt, kpeg, mpeg, icase_in )
384 end if
385
386
387
388 return
389 end subroutine mosaic_cloudchem_1box
390
391
392
393 !-----------------------------------------------------------------------
394 subroutine interface_to_aqoperator1( &
395 istat_aqop, &
396 dtstepc, &
397 rbox, gas_aqfrac_box, &
398 rbulk_cwaer, lptr_yyy_cwaer, &
399 co2_mixrat_in, photol_no2_in, iradical_in, idecomp_hmsa_hso5, &
400 yaq_beg, yaq_end, ph_cmuaq_cur, &
401 numgas_aqfrac, id, it, jt, kt, kpeg, mpeg, ktau, icase )
402
403 use module_state_description, only: &
404 num_chem, param_first_scalar, p_qc, &
405 p_nh3, p_hno3, p_hcl, p_sulf, p_hcho, &
406 p_ora1, p_so2, p_h2o2, p_o3, p_ho, &
407 p_ho2, p_no3, p_no, p_no2, p_hono, &
408 p_pan, p_ch3o2, p_ch3oh, p_op1
409
410 use module_data_cmu_bulkaqchem, only: &
411 meqn1max, naers, ngas, &
412 na4, naa, nac, nae, nah, nahmsa, nahso5, &
413 nan, nao, nar, nas, naw, &
414 ng4, nga, ngc, ngch3co3h, ngch3o2, ngch3o2h, ngch3oh, &
415 ngh2o2, nghcho, nghcooh, nghno2, ngho2, &
416 ngn, ngno, ngno2, ngno3, ngo3, ngoh, ngpan, ngso2
417
418 use module_cmu_bulkaqchem, only: aqoperator1
419
420 use module_data_mosaic_asect, only: &
421 maxd_asize, maxd_atype, &
422 cw_phase, nsize_aer, ntype_aer, &
423 lptr_so4_aer, lptr_no3_aer, lptr_cl_aer, lptr_co3_aer, &
424 lptr_msa_aer, lptr_nh4_aer, lptr_na_aer, lptr_ca_aer, &
425 lptr_oin_aer, lptr_bc_aer, lptr_oc_aer, &
426 mw_cl_aer, mw_na_aer, mw_nh4_aer, mw_no3_aer, mw_so4_aer
427
428 use module_data_mosaic_other, only: &
429 aboxtest_units_convert, cairclm, &
430 ktemp, l2maxd, ptotclm, rcldwtr_sub
431
432
433 implicit none
434
435 ! subr arguments
436 integer, intent(in) :: &
437 iradical_in, idecomp_hmsa_hso5, &
438 numgas_aqfrac, id, it, jt, kt, kpeg, mpeg, ktau, icase
439 integer, intent(inout) :: &
440 istat_aqop
441
442 integer, intent(in) :: lptr_yyy_cwaer(maxd_asize,maxd_atype,nyyy)
443
444 real, intent(in) :: &
445 dtstepc, co2_mixrat_in, &
446 photol_no2_in
447
448 real, intent(inout) :: ph_cmuaq_cur
449
450 real, intent(inout), dimension( 1:l2maxd ) :: rbox
451
452 real, intent(inout), dimension( nyyy, 2 ) :: rbulk_cwaer
453
454 real, intent(inout), dimension( 1:numgas_aqfrac ) :: gas_aqfrac_box
455
456 real, intent(inout), dimension( meqn1max ) :: yaq_beg, yaq_end
457
458
459 ! local variables
460 integer :: i, iphase, isize, itype
461 integer :: iaq, istat_fatal, istat_warn
462 integer :: l, lunxx, lyyy
463 integer :: p1st
464
465 real, parameter :: eps=0.622 ! (mw h2o)/(mw air)
466
467
468 real :: cair_moleperm3
469 real :: dum, dumb
470 real :: factgas, factlwc, factpatm, factphoto
471 real :: factaerbc, factaercl, factaerna, factaernh4, &
472 factaerno3, factaeroc, factaeroin, factaerso4
473 real :: lwc
474 real :: p_atm, photo_in
475 real :: rh
476 real :: temp, tstep_beg_sec, tstep_end_sec
477 real :: totsulf_beg, totsulf_end
478 real :: gas(ngas), aerosol(naers)
479 real :: gas_aqfrac_cmu(ngas)
480
481 double precision tstep_beg_sec_dp, tstep_end_sec_dp, &
482 temp_dp, p_atm_dp, lwc_dp, rh_dp, &
483 co2_mixrat_in_dp, photo_in_dp, ph_cmuaq_cur_dp
484 double precision gas_dp(ngas), gas_aqfrac_cmu_dp(ngas), &
485 aerosol_dp(naers), yaq_beg_dp(meqn1max), yaq_end_dp(meqn1max)
486
487
488
489 p1st = param_first_scalar
490
491 !
492 ! units conversion factors
493 ! 'cmuaq-bulk' value = pegasus value X factor
494 !
495 ! [pres in atmospheres] = [pres in dynes/cm2] * factpatm
496 factpatm = 1.0/1.01325e6
497 ! [cldwtr in g-h2o/m3-air] = [cldwtr in mole-h2o/mole-air] * factlwc
498 factlwc = 28.966*eps*1.0e6*cairclm(kpeg)
499 ! [aq photolysis rate scaling factor in --] = [jno2 in 1/min] * factphoto
500 factphoto = 1.6
501
502 ! [gas in ppm] = [gas in mole/mole-air] * factgas
503 factgas = 1.0e6
504
505 ! [aerosol in ug/m3-air] = [aerosol in mole/mole-air] * factaer
506 dum = cairclm(kpeg)*1.0e12
507 factaerso4 = dum*mw_so4_aer
508 factaerno3 = dum*mw_no3_aer
509 factaercl = dum*mw_cl_aer
510 factaernh4 = dum*mw_nh4_aer
511 factaerna = dum*mw_na_aer
512 factaeroin = dum
513 factaeroc = dum
514 factaerbc = dum
515
516 ! rce 2005-jul-11 - use same molecular weights here as in cmu code
517 ! factaerso4 = dum*96.0
518 ! factaerno3 = dum*62.0
519 ! factaercl = dum*35.5
520 ! factaernh4 = dum*18.0
521 ! factaerna = dum*23.0
522
523 ! If aboxtest_units_convert=10, turn off units conversions both here
524 ! and in module_mosaic. This is for testing, to allow exact agreements.
525 if (aboxtest_units_convert .eq. 10) then
526 factpatm = 1.0
527 factlwc = 1.0
528 factphoto = 1.0
529 factgas = 1.0
530 factaerso4 = 1.0
531 factaerno3 = 1.0
532 factaercl = 1.0
533 factaernh4 = 1.0
534 factaerna = 1.0
535 factaeroin = 1.0
536 factaeroc = 1.0
537 factaerbc = 1.0
538 end if
539
540 !
541 ! map from rbox to gas,aerosol
542 !
543 temp = rbox(ktemp)
544
545 lwc = rcldwtr_sub(kpeg,mpeg) * factlwc
546 p_atm = ptotclm(kpeg) * factpatm
547
548 ! rce 2005-jul-11 - set p_atm so that cmu code's cair will match cairclm
549 p_atm = cairclm(kpeg)*1.0e3*0.082058e0*temp
550
551 photo_in = photol_no2_in * factphoto
552
553 rh = 1.0
554 iaq = 1
555
556 tstep_beg_sec = 0.0
557 tstep_end_sec = dtstepc
558
559 ! map gases and convert to ppm
560 gas(:) = 0.0
561
562 gas(nga ) = rbox(p_nh3 ) * factgas
563 gas(ngn ) = rbox(p_hno3 ) * factgas
564 gas(ngc ) = rbox(p_hcl ) * factgas
565 gas(ng4 ) = rbox(p_sulf ) * factgas
566
567 gas(nghcho ) = rbox(p_hcho ) * factgas
568 gas(nghcooh ) = rbox(p_ora1 ) * factgas
569 gas(ngso2 ) = rbox(p_so2 ) * factgas
570 gas(ngh2o2 ) = rbox(p_h2o2 ) * factgas
571 gas(ngo3 ) = rbox(p_o3 ) * factgas
572 gas(ngoh ) = rbox(p_ho ) * factgas
573 gas(ngho2 ) = rbox(p_ho2 ) * factgas
574 gas(ngno3 ) = rbox(p_no3 ) * factgas
575
576 gas(ngno ) = rbox(p_no ) * factgas
577 gas(ngno2 ) = rbox(p_no2 ) * factgas
578 gas(nghno2 ) = rbox(p_hono ) * factgas
579 gas(ngpan ) = rbox(p_pan ) * factgas
580 gas(ngch3o2 ) = rbox(p_ch3o2 ) * factgas
581 gas(ngch3oh ) = rbox(p_ch3oh ) * factgas
582 gas(ngch3o2h) = rbox(p_op1 ) * factgas
583
584 ! compute bulk activated-aerosol mixing ratios
585 aerosol(:) = 0.0
586 rbulk_cwaer(:,:) = 0.0
587
588 iphase = cw_phase
589 do itype = 1, ntype_aer
590 do isize = 1, nsize_aer(itype)
591
592 do lyyy = 1, nyyy
593
594 l = lptr_yyy_cwaer(isize,itype,lyyy)
595 if (l .ge. p1st) rbulk_cwaer(lyyy,1) = rbulk_cwaer(lyyy,1) + rbox(l)
596
597 end do
598
599 end do
600 end do
601
602 ! map them to 'aerosol' array and convert to ug/m3
603 aerosol(na4) = rbulk_cwaer(l_so4_aqyy,1) * factaerso4
604 aerosol(nan) = rbulk_cwaer(l_no3_aqyy,1) * factaerno3
605 aerosol(nac) = rbulk_cwaer(l_cl_aqyy, 1) * factaercl
606 aerosol(naa) = rbulk_cwaer(l_nh4_aqyy,1) * factaernh4
607 aerosol(nas) = rbulk_cwaer(l_na_aqyy, 1) * factaerna
608 aerosol(nar) = rbulk_cwaer(l_oin_aqyy,1) * factaeroin
609 aerosol(nae) = rbulk_cwaer(l_bc_aqyy, 1) * factaerbc
610 aerosol(nao) = rbulk_cwaer(l_oc_aqyy, 1) * factaeroc
611
612
613 !
614 ! make call to aqoperator1
615 !
616 #if defined ( ccboxtest_box_testing_active)
617 lunxx = 87
618 lunxx = -1
619 if (lunxx .gt. 0) then
620 write(lunxx,*)
621 write(lunxx,*)
622 write(lunxx,*) 'interface_to_aqoperator1 - icase, irad, idecomp'
623 write(lunxx,9870) icase, iradical_in, idecomp_hmsa_hso5
624 write(lunxx,*) 'it, jt, kt, kpeg, mpeg, ktau'
625 write(lunxx,9870) it, jt, kt, kpeg, mpeg, ktau
626 write(lunxx,*) 'temp, p_atm, lwc, photo, co2'
627 write(lunxx,9875) temp, p_atm, lwc, photo_in, co2_mixrat_in
628 write(lunxx,*) 'ptot, cair, rcldwtr_clm, dt_sec'
629 write(lunxx,9875) ptotclm(kpeg), cairclm(kpeg), &
630 rcldwtr_sub(kpeg,mpeg), (tstep_end_sec-tstep_beg_sec)
631 write(lunxx,*) 'gas (1=nh3, 2=hno3, 3=hcl, 4=h2so4, 11=so2, 12=h2o2, 18=o3)'
632 write(lunxx,9875) gas
633 write(lunxx,*) 'rbox(nh3, hno3, hcl, h2so4, so2, h2o2, o3)'
634 write(lunxx,9875) rbox(p_nh3), rbox(p_hno3), rbox(p_hcl), &
635 rbox(p_sulf), rbox(p_so2), rbox(p_h2o2), rbox(p_o3)
636 write(lunxx,*) 'aerosol (1=na, 3=nh4, 4=no3, 5=cl, 6=so4, 8=ec, 9=oc, 10=crus)'
637 write(lunxx,9875) aerosol
638 write(lunxx,*) 'rbulk_cwaer (1=so4, 2=no3, 3-cl, 4=nh4, 5=na, 6=oin, 7=bc, 8=oc)'
639 write(lunxx,9875) rbulk_cwaer(:,1)
640 ! if (icase .ge. 3) then
641 ! write(*,*) &
642 ! '*** stopping in interface_to_aqop1 at icase =', icase
643 ! stop
644 ! end if
645 end if
646 9870 format( 8i5 )
647 9875 format( 5(1pe14.6) )
648 #endif
649
650 #if 0
651 ! Print outs for debugging of aqoperator1... wig, 26-Oct-2005
652 !!$ if( (id == 1 .and. ktau >= 207 ) .or. &
653 !!$ (id == 2 .and. ktau >= 610 ) .or. &
654 !!$ (id == 3 .and. ktau >= 1830 ) ) then
655 write(6,'(a)') '---Begin input for aqoperator1---'
656 write(6,'(a,4i)') 'id, it, jt, kt =', id, it, jt, kt
657 write(6,'(a,1p,2e20.12)') 'tstep_beg_sec, tstep_end_sec = ',tstep_beg_sec, tstep_end_sec
658 do l=1,ngas
659 write(6,'("gas(",i2,") = ",1p,1e20.12)') l, gas(l)
660 end do
661 do l=1,naers
662 write(6,'("aerosol(",i2,") = ",1p,1e20.12)') l, aerosol(l)
663 end do
664 write(6,'(a,1p,4e20.12)') "temp, p_atm, lwc, rh = ", temp, p_atm, lwc, rh
665 write(6,'(a,1p,2e20.12)') "co2_mixrat_in, photo_in = ", co2_mixrat_in, photo_in
666 write(6,'(a,3i)') " iradical_in, idecomp_hmsa_hso5, iaq = ", iradical_in, idecomp_hmsa_hso5, iaq
667 write(6,'(a)') "---End input for aqoperator1---"
668 !!$ end if
669 #endif
670
671
672 ! convert arguments to double prec
673 tstep_beg_sec_dp = 0.0d0
674 if (tstep_beg_sec .ne. 0.0) tstep_beg_sec_dp = tstep_beg_sec
675 tstep_end_sec_dp = 0.0d0
676 if (tstep_end_sec .ne. 0.0) tstep_end_sec_dp = tstep_end_sec
677 temp_dp = 0.0d0
678 if (temp .ne. 0.0) temp_dp = temp
679 p_atm_dp = 0.0d0
680 if (p_atm .ne. 0.0) p_atm_dp = p_atm
681 lwc_dp = 0.0d0
682 if (lwc .ne. 0.0) lwc_dp = lwc
683 rh_dp = 0.0d0
684 if (rh .ne. 0.0) rh_dp = rh
685 co2_mixrat_in_dp = 0.0d0
686 if (co2_mixrat_in .ne. 0.0) co2_mixrat_in_dp = co2_mixrat_in
687 photo_in_dp = 0.0d0
688 if (photo_in .ne. 0.0) photo_in_dp = photo_in
689 ph_cmuaq_cur_dp = 0.0d0
690 if (ph_cmuaq_cur .ne. 0.0) ph_cmuaq_cur_dp = ph_cmuaq_cur
691
692 do i = 1, ngas
693 gas_dp(i) = 0.0d0
694 if (gas(i) .ne. 0.0) gas_dp(i) = gas(i)
695 end do
696 do i = 1, naers
697 aerosol_dp(i) = 0.0d0
698 if (aerosol(i) .ne. 0.0) aerosol_dp(i) = aerosol(i)
699 end do
700 do i = 1, ngas
701 gas_aqfrac_cmu_dp(i) = 0.0d0
702 if (gas_aqfrac_cmu(i) .ne. 0.0) gas_aqfrac_cmu_dp(i) = gas_aqfrac_cmu(i)
703 end do
704 do i = 1, meqn1max
705 yaq_beg_dp(i) = 0.0d0
706 if (yaq_beg(i) .ne. 0.0) yaq_beg_dp(i) = yaq_beg(i)
707 end do
708 do i = 1, meqn1max
709 yaq_end_dp(i) = 0.0d0
710 if (yaq_end(i) .ne. 0.0) yaq_end_dp(i) = yaq_end(i)
711 end do
712
713
714 ! total sulfur species conc as sulfate (ug/m3)
715 cair_moleperm3 = 1.0e3*p_atm_dp/(0.082058e0*temp_dp)
716 totsulf_beg = ( aerosol_dp(na4)/96. &
717 + aerosol_dp(nahso5)/113. + aerosol_dp(nahmsa)/111. &
718 + (gas_dp(ngso2) + gas_dp(ng4))*cair_moleperm3 )*96.0
719
720 ! call aqoperator1( &
721 ! istat_fatal, istat_warn, &
722 ! tstep_beg_sec, tstep_end_sec, &
723 ! gas, aerosol, gas_aqfrac_cmu, &
724 ! temp, p_atm, lwc, rh, &
725 ! co2_mixrat_in, photo_in, iradical_in, idecomp_hmsa_hso5, iaq, &
726 ! yaq_beg, yaq_end, ph_cmuaq_cur )
727
728 call aqoperator1( &
729 istat_fatal, istat_warn, &
730 tstep_beg_sec_dp, tstep_end_sec_dp, &
731 gas_dp, aerosol_dp, gas_aqfrac_cmu_dp, &
732 temp_dp, p_atm_dp, lwc_dp, rh_dp, &
733 co2_mixrat_in_dp, photo_in_dp, iradical_in, idecomp_hmsa_hso5, iaq, &
734 yaq_beg_dp, yaq_end_dp, ph_cmuaq_cur_dp )
735
736 totsulf_end = ( aerosol_dp(na4)/96. &
737 + aerosol_dp(nahso5)/113. + aerosol_dp(nahmsa)/111. &
738 + (gas_dp(ngso2) + gas_dp(ng4))*cair_moleperm3 )*96.0
739
740
741 ! convert arguments back to single prec
742 tstep_beg_sec = tstep_beg_sec_dp
743 tstep_end_sec = tstep_end_sec_dp
744 temp = temp_dp
745 p_atm = p_atm_dp
746 lwc = lwc_dp
747 rh = rh_dp
748 ! co2_mixrat_in = co2_mixrat_in_dp ! this has intent(in)
749 ! photo_in = photo_in_dp ! this has intent(in)
750 ph_cmuaq_cur = ph_cmuaq_cur_dp
751
752 do i = 1, ngas
753 gas(i) = gas_dp(i)
754 end do
755 do i = 1, naers
756 aerosol(i) = aerosol_dp(i)
757 end do
758 do i = 1, ngas
759 gas_aqfrac_cmu(i) = gas_aqfrac_cmu_dp(i)
760 end do
761 do i = 1, meqn1max
762 yaq_beg(i) = yaq_beg_dp(i)
763 end do
764 do i = 1, meqn1max
765 yaq_end(i) = yaq_end_dp(i)
766 end do
767
768
769 !
770 ! warning message when status flags are non-zero
771 !
772 istat_aqop = 0
773 if (istat_fatal .ne. 0) then
774 write(6,*) &
775 '*** mosaic_cloudchem_driver, subr interface_to_aqoperator1'
776 write(6,'(a,4i5,2i10)') &
777 ' id,it,jt,kt, istat_fatal, warn =', &
778 id, it, jt, kt, istat_fatal, istat_warn
779 istat_aqop = -10
780 end if
781
782 !
783 ! warning message when sulfur mass balance error exceeds the greater
784 ! of (1.0e-3 ug/m3) OR (1.0e-3 X total sulfur mixing ratio)
785 !
786 dum = totsulf_end - totsulf_beg
787 dumb = max( totsulf_beg, totsulf_end )
788 if (abs(dum) .gt. max(1.0e-3,1.0e-3*dumb)) then
789 write(6,*) &
790 '*** mosaic_cloudchem_driver, sulfur balance warning'
791 write(6,'(a,4i5,1p,3e12.4)') &
792 ' id,it,jt,kt, total_sulfur_beg, _end, _error =', &
793 id, it, jt, kt, totsulf_beg, totsulf_end, dum
794 end if
795
796 !
797 ! map from gas,aerosol to rbox
798 !
799 rbox(p_nh3 ) = gas(nga ) / factgas
800 rbox(p_hno3 ) = gas(ngn ) / factgas
801 rbox(p_hcl ) = gas(ngc ) / factgas
802 rbox(p_sulf ) = gas(ng4 ) / factgas
803
804 rbox(p_hcho ) = gas(nghcho ) / factgas
805 rbox(p_ora1 ) = gas(nghcooh ) / factgas
806 rbox(p_so2 ) = gas(ngso2 ) / factgas
807 rbox(p_h2o2 ) = gas(ngh2o2 ) / factgas
808 rbox(p_o3 ) = gas(ngo3 ) / factgas
809 rbox(p_ho ) = gas(ngoh ) / factgas
810 rbox(p_ho2 ) = gas(ngho2 ) / factgas
811 rbox(p_no3 ) = gas(ngno3 ) / factgas
812
813 rbox(p_no ) = gas(ngno ) / factgas
814 rbox(p_no2 ) = gas(ngno2 ) / factgas
815 rbox(p_hono ) = gas(nghno2 ) / factgas
816 rbox(p_pan ) = gas(ngpan ) / factgas
817 rbox(p_ch3o2 ) = gas(ngch3o2 ) / factgas
818 rbox(p_ch3oh ) = gas(ngch3oh ) / factgas
819 rbox(p_op1 ) = gas(ngch3o2h) / factgas
820
821 gas_aqfrac_box(:) = 0.0
822
823 if (p_nh3 .le. numgas_aqfrac) &
824 gas_aqfrac_box(p_nh3 ) = gas_aqfrac_cmu(nga )
825 if (p_hno3 .le. numgas_aqfrac) &
826 gas_aqfrac_box(p_hno3 ) = gas_aqfrac_cmu(ngn )
827 if (p_hcl .le. numgas_aqfrac) &
828 gas_aqfrac_box(p_hcl ) = gas_aqfrac_cmu(ngc )
829 if (p_sulf .le. numgas_aqfrac) &
830 gas_aqfrac_box(p_sulf ) = gas_aqfrac_cmu(ng4 )
831
832 if (p_hcho .le. numgas_aqfrac) &
833 gas_aqfrac_box(p_hcho ) = gas_aqfrac_cmu(nghcho )
834 if (p_ora1 .le. numgas_aqfrac) &
835 gas_aqfrac_box(p_ora1 ) = gas_aqfrac_cmu(nghcooh )
836 if (p_so2 .le. numgas_aqfrac) &
837 gas_aqfrac_box(p_so2 ) = gas_aqfrac_cmu(ngso2 )
838 if (p_h2o2 .le. numgas_aqfrac) &
839 gas_aqfrac_box(p_h2o2 ) = gas_aqfrac_cmu(ngh2o2 )
840 if (p_o3 .le. numgas_aqfrac) &
841 gas_aqfrac_box(p_o3 ) = gas_aqfrac_cmu(ngo3 )
842 if (p_ho .le. numgas_aqfrac) &
843 gas_aqfrac_box(p_ho ) = gas_aqfrac_cmu(ngoh )
844 if (p_ho2 .le. numgas_aqfrac) &
845 gas_aqfrac_box(p_ho2 ) = gas_aqfrac_cmu(ngho2 )
846 if (p_no3 .le. numgas_aqfrac) &
847 gas_aqfrac_box(p_no3 ) = gas_aqfrac_cmu(ngno3 )
848
849 if (p_no .le. numgas_aqfrac) &
850 gas_aqfrac_box(p_no ) = gas_aqfrac_cmu(ngno )
851 if (p_no2 .le. numgas_aqfrac) &
852 gas_aqfrac_box(p_no2 ) = gas_aqfrac_cmu(ngno2 )
853 if (p_hono .le. numgas_aqfrac) &
854 gas_aqfrac_box(p_hono ) = gas_aqfrac_cmu(nghno2 )
855 if (p_pan .le. numgas_aqfrac) &
856 gas_aqfrac_box(p_pan ) = gas_aqfrac_cmu(ngpan )
857 if (p_ch3o2 .le. numgas_aqfrac) &
858 gas_aqfrac_box(p_ch3o2 ) = gas_aqfrac_cmu(ngch3o2 )
859 if (p_ch3oh .le. numgas_aqfrac) &
860 gas_aqfrac_box(p_ch3oh ) = gas_aqfrac_cmu(ngch3oh )
861 if (p_op1 .le. numgas_aqfrac) &
862 gas_aqfrac_box(p_op1 ) = gas_aqfrac_cmu(ngch3o2h)
863
864 rbulk_cwaer(l_so4_aqyy,2) = aerosol(na4) / factaerso4
865 rbulk_cwaer(l_no3_aqyy,2) = aerosol(nan) / factaerno3
866 rbulk_cwaer(l_cl_aqyy, 2) = aerosol(nac) / factaercl
867 rbulk_cwaer(l_nh4_aqyy,2) = aerosol(naa) / factaernh4
868 rbulk_cwaer(l_na_aqyy, 2) = aerosol(nas) / factaerna
869 rbulk_cwaer(l_oin_aqyy,2) = aerosol(nar) / factaeroin
870 rbulk_cwaer(l_bc_aqyy, 2) = aerosol(nae) / factaerbc
871 rbulk_cwaer(l_oc_aqyy, 2) = aerosol(nao) / factaeroc
872
873
874 #if defined ( ccboxtest_box_testing_active)
875 lunxx = 87
876 lunxx = -1
877 if (lunxx .gt. 0) then
878 write(lunxx,*)
879 write(lunxx,*) 'interface_to_aqoperator1 - after call'
880 write(lunxx,*) 'gas (1=nh3, 2=hno3, 3=hcl, 4=h2so4, 11=so2, 12=h2o2, 18=o3)'
881 write(lunxx,9875) gas
882 write(lunxx,*) 'rbox(nh3, hno3, hcl, h2so4, so2, h2o2, o3)'
883 write(lunxx,9875) rbox(p_nh3), rbox(p_hno3), rbox(p_hcl), &
884 rbox(p_sulf), rbox(p_so2), rbox(p_h2o2), rbox(p_o3)
885 write(lunxx,*) 'aerosol (1=na, 3=nh4, 4=no3, 5=cl, 6=so4, 8=ec, 9=oc, 10=crus)'
886 write(lunxx,9875) aerosol
887 write(lunxx,*) 'rbulk_cwaer (1=so4, 2=no3, 3-cl, 4=nh4, 5=na, 6=oin, 7=bc, 8=oc)'
888 write(lunxx,9875) rbulk_cwaer(:,2)
889 write(lunxx,*) 'ph_cmuaq_cur'
890 write(lunxx,9875) ph_cmuaq_cur
891 if (icase .ge. 10) then
892 write(*,*) &
893 '*** stopping in interface_to_aqop1 at icase =', icase
894 stop
895 end if
896 end if
897 #endif
898
899
900 return
901 end subroutine interface_to_aqoperator1
902
903
904
905 !-----------------------------------------------------------------------
906 subroutine partition_cldwtr( &
907 rbox, fr_partit_cw, &
908 it, jt, kt, kpeg, mpeg, icase )
909
910 use module_state_description, only: &
911 param_first_scalar
912
913 use module_data_mosaic_asect, only: &
914 maxd_asize, maxd_atype, &
915 cw_phase, nsize_aer, ntype_aer, ncomp_aer, &
916 massptr_aer, numptr_aer, &
917 dens_aer, mw_aer, volumlo_sect, volumhi_sect
918
919 use module_data_mosaic_other, only: &
920 aboxtest_units_convert, cairclm, &
921 ktemp, l2maxd, ptotclm, rcldwtr_sub
922
923
924 implicit none
925
926 ! subr arguments
927 integer, intent(in) :: it, jt, kt, kpeg, mpeg, icase
928
929 real, intent(inout), dimension( 1:l2maxd ) :: rbox
930
931 real, intent(inout), dimension( maxd_asize, maxd_atype ) :: &
932 fr_partit_cw
933
934 ! local variables
935 integer :: iphase, isize, itype
936 integer :: jdone_mass, jdone_numb, jpos, jpos_mass, jpos_numb
937 integer :: l, ll, lunxx
938 integer :: p1st
939
940 real, parameter :: partit_wght_mass = 0.5
941
942 real :: dum, duma, dumb, dumc, dummass, dumnumb, dumvolu
943 real :: tmass, tnumb, umass, unumb, wmass, wnumb
944 real, dimension( maxd_asize, maxd_atype ) :: fmass, fnumb, xmass, xnumb
945
946 p1st = PARAM_FIRST_SCALAR
947
948 iphase = cw_phase
949 tmass = 0.0
950 tnumb = 0.0
951 umass = 0.0
952 unumb = 0.0
953
954 ! compute
955 ! xmass, xnumb = mass, number mixing ratio for a bin
956 ! tmass, tnumb = sum over all bins of xmass, xnumb
957 ! umass, unumb = max over all bins of xmass, xnumb
958 ! set xmass, xnumb = 0.0 if bin mass, numb < 1.0e-37
959 ! constrain xnumb so that mean particle volume is
960 ! within bin boundaries
961 do itype = 1, ntype_aer
962 do isize = 1, nsize_aer(itype)
963 dummass = 0.0
964 dumvolu = 0.0
965 do ll = 1, ncomp_aer(itype)
966 l = massptr_aer(ll,isize,itype,iphase)
967 if (l .ge. p1st) then
968 dum = max( 0.0, rbox(l) )*mw_aer(ll,itype)
969 dummass = dummass + dum
970 dumvolu = dumvolu + dum/dens_aer(ll,itype)
971 end if
972 end do
973
974 l = numptr_aer(isize,itype,iphase)
975 dumnumb = max( 0.0, rbox(l) )
976 if (dumnumb .gt. dumvolu/volumlo_sect(isize,itype)) then
977 dumnumb = dumvolu/volumlo_sect(isize,itype)
978 rbox(l) = dumnumb
979 else if (dumnumb .lt. dumvolu/volumhi_sect(isize,itype)) then
980 dumnumb = dumvolu/volumhi_sect(isize,itype)
981 rbox(l) = dumnumb
982 end if
983
984 if (dummass .lt. 1.0e-37) dummass = 0.0
985 xmass(isize,itype) = dummass
986 if (dumnumb .lt. 1.0e-37) dumnumb = 0.0
987 xnumb(isize,itype) = dumnumb
988
989 tmass = tmass + xmass(isize,itype)
990 tnumb = tnumb + xnumb(isize,itype)
991 umass = max( umass, xmass(isize,itype) )
992 unumb = max( unumb, xnumb(isize,itype) )
993 end do
994 end do
995
996 ! compute
997 ! fmass, fnumb = fraction of total mass, number that is in a bin
998 ! if tmass<1e-35 and umass>0, set fmass=1 for bin with largest xmass
999 ! if tmass<1e-35 and umass=0, set fmass=0 for all
1000 jdone_mass = 0
1001 jdone_numb = 0
1002 jpos_mass = 0
1003 jpos_numb = 0
1004 do itype = 1, ntype_aer
1005 do isize = 1, nsize_aer(itype)
1006 fmass(isize,itype) = 0.0
1007 if (tmass .ge. 1.0e-35) then
1008 fmass(isize,itype) = xmass(isize,itype)/tmass
1009 else if (umass .gt. 0.0) then
1010 if ( (jdone_mass .eq. 0) .and. &
1011 (xmass(isize,itype) .eq. umass) ) then
1012 jdone_mass = 1
1013 fmass(isize,itype) = 1.0
1014 end if
1015 end if
1016 if (fmass(isize,itype) .gt. 0) jpos_mass = jpos_mass + 1
1017
1018 fnumb(isize,itype) = 0.0
1019 if (tnumb .ge. 1.0e-35) then
1020 fnumb(isize,itype) = xnumb(isize,itype)/tnumb
1021 else if (unumb .gt. 0.0) then
1022 if ( (jdone_numb .eq. 0) .and. &
1023 (xnumb(isize,itype) .eq. unumb) ) then
1024 jdone_numb = 1
1025 fnumb(isize,itype) = 1.0
1026 end if
1027 end if
1028 if (fnumb(isize,itype) .gt. 0) jpos_numb = jpos_numb + 1
1029 end do
1030 end do
1031
1032 ! if only 1 bin has fmass or fnumb > 0, set value to 1.0 exactly
1033 if ((jpos_mass .eq. 1) .or. (jpos_numb .eq. 1)) then
1034 do itype = 1, ntype_aer
1035 do isize = 1, nsize_aer(itype)
1036 if (jpos_mass .eq. 1) then
1037 if (fmass(isize,itype) .gt. 0) fmass(isize,itype) = 1.0
1038 end if
1039 if (jpos_numb .eq. 1) then
1040 if (fnumb(isize,itype) .gt. 0) fnumb(isize,itype) = 1.0
1041 end if
1042 end do
1043 end do
1044 end if
1045
1046 !
1047 ! compute fr_partit_cw as weighted average of fmass & fnumb, except
1048 ! if tmass<1e-35 and umass=0, use only fnumb
1049 ! if tnumb<1e-35 and unumb=0, use only fmass
1050 ! if tmass,tnumb<1e-35 and umass,unumb=0,
1051 ! set fr_partit_cw=1 for center bin of itype=1
1052 !
1053 fr_partit_cw(:,:) = 0.0
1054 if ((jpos_mass .eq. 0) .and. (jpos_numb .eq. 0)) then
1055 itype = 1
1056 isize = (nsize_aer(itype)+1)/2
1057 fr_partit_cw(isize,itype) = 1.0
1058
1059 else if (jpos_mass .eq. 0) then
1060 fr_partit_cw(:,:) = fnumb(:,:)
1061
1062 else if (jpos_numb .eq. 0) then
1063 fr_partit_cw(:,:) = fmass(:,:)
1064
1065 else
1066 wmass = max( 0.0, min( 1.0, partit_wght_mass ) )
1067 wnumb = 1.0 - wmass
1068 fr_partit_cw(:,:) = wmass*fmass(:,:) + wnumb*fnumb(:,:)
1069
1070 jpos = 0
1071 do itype = 1, ntype_aer
1072 do isize = 1, nsize_aer(itype)
1073 if (fr_partit_cw(isize,itype) .gt. 0.0) jpos = jpos + 1
1074 end do
1075 end do
1076
1077 ! if only 1 bin has fr_partit_cw > 0, set value to 1.0 exactly
1078 if (jpos .eq. 1) then
1079 do itype = 1, ntype_aer
1080 do isize = 1, nsize_aer(itype)
1081 if (fr_partit_cw(isize,itype) .gt. 0.0) &
1082 fr_partit_cw(isize,itype) = 1.0
1083 end do
1084 end do
1085 end if
1086 end if
1087
1088
1089 #if defined ( ccboxtest_box_testing_active)
1090 ! diagnostics when lunxx > 0
1091 lunxx = 86
1092 lunxx = -1
1093 if (lunxx .gt. 0) then
1094 if (icase .le. 9) then
1095 write(lunxx,9800)
1096 write(lunxx,9800)
1097 write(lunxx,9800) &
1098 'partition_cldwtr - icase, jpos, jpos_mass, jpos_numb'
1099 write(lunxx,9810) icase, jpos, jpos_mass, jpos_numb
1100 write(lunxx,9800) 'tmass, umass, wmass'
1101 write(lunxx,9820) tmass, umass, wmass
1102 write(lunxx,9800) 'tnumb, unumb, wnumb'
1103 write(lunxx,9820) tnumb, unumb, wnumb
1104 write(lunxx,9800) 'xmass, fmass, xnumb, fnumb, fr_partit_cw'
1105 duma = 0.0
1106 dumb = 0.0
1107 dumc = 0.0
1108 do itype = 1, ntype_aer
1109 do isize = 1, nsize_aer(itype)
1110 write(lunxx,9820) xmass(isize,itype), fmass(isize,itype), &
1111 xnumb(isize,itype), fnumb(isize,itype), &
1112 fr_partit_cw(isize,itype)
1113 duma = duma + fmass(isize,itype)
1114 dumb = dumb + fnumb(isize,itype)
1115 dumc = dumc + fr_partit_cw(isize,itype)
1116 end do
1117 end do
1118 write(lunxx,9800) &
1119 'sum_fmass-1.0, sum_fnumb-1.0, sum_fr_partit-1.0'
1120 write(lunxx,9820) (duma-1.0), (dumb-1.0), (dumc-1.0)
1121 if (icase .eq. -2) then
1122 write(*,*) '*** stopping in partition_cldwtr at icase =', icase
1123 stop
1124 end if
1125 9800 format( a )
1126 9810 format( 5i10 )
1127 9820 format( 5(1pe10.2) )
1128 end if
1129 end if
1130 #endif
1131
1132
1133 return
1134 end subroutine partition_cldwtr
1135
1136
1137
1138 !-----------------------------------------------------------------------
1139 subroutine distribute_bulk_changes( &
1140 rbox, rbox_sv1, fr_partit_cw, &
1141 rbulk_cwaer, lptr_yyy_cwaer, &
1142 it, jt, kt, kpeg, mpeg, icase )
1143
1144 use module_state_description, only: &
1145 param_first_scalar
1146
1147 use module_data_mosaic_asect, only: &
1148 maxd_asize, maxd_atype, &
1149 cw_phase, nsize_aer, ntype_aer, &
1150 lptr_so4_aer, lptr_no3_aer, lptr_cl_aer, lptr_co3_aer, &
1151 lptr_msa_aer, lptr_nh4_aer, lptr_na_aer, lptr_ca_aer, &
1152 lptr_oin_aer, lptr_bc_aer, lptr_oc_aer
1153
1154 use module_data_mosaic_other, only: l2maxd, lunout, name
1155
1156
1157 implicit none
1158
1159 ! subr arguments
1160 integer, intent(in) :: it, jt, kt, kpeg, mpeg, icase
1161
1162 integer, intent(in) :: lptr_yyy_cwaer(maxd_asize,maxd_atype,nyyy)
1163
1164 real, intent(inout), dimension( 1:l2maxd ) :: rbox, rbox_sv1
1165
1166 real, intent(in), dimension( maxd_asize, maxd_atype ) :: &
1167 fr_partit_cw
1168
1169 real, intent(in), dimension( nyyy, 2 ) :: rbulk_cwaer
1170
1171
1172 ! local variables
1173 integer :: iphase, isize, itype
1174 integer :: idone, icount, ncount
1175 integer :: jpos, jpos_sv
1176 integer :: l, lunxx, lunxxaa, lunxxbb, lyyy
1177 integer :: p1st
1178
1179 real :: duma, dumb, dumc
1180 real :: fr, frsum_cur
1181 real :: fr_cur(maxd_asize,maxd_atype)
1182 real :: del_r_current, del_r_remain
1183 real :: del_rbulk_cwaer(nyyy)
1184
1185
1186 p1st = param_first_scalar
1187
1188 do lyyy = 1, nyyy
1189 del_rbulk_cwaer(lyyy) = rbulk_cwaer(lyyy,2) - rbulk_cwaer(lyyy,1)
1190 end do
1191
1192 iphase = cw_phase
1193
1194
1195 jpos = 0
1196 do itype = 1, ntype_aer
1197 do isize = 1, nsize_aer(itype)
1198 if (fr_partit_cw(isize,itype) .gt. 0) jpos = jpos + 1
1199 end do
1200 end do
1201 jpos_sv = jpos
1202
1203 !
1204 ! distribution is trivial when only 1 bin has fr_partit_cw > 0
1205 !
1206 if (jpos_sv .eq. 1) then
1207 do lyyy = 1, nyyy
1208
1209 do itype = 1, ntype_aer
1210 do isize = 1, nsize_aer(itype)
1211 fr = fr_partit_cw(isize,itype)
1212 if (fr .eq. 1.0) then
1213 l = lptr_yyy_cwaer(isize,itype,lyyy)
1214 if (l .ge. p1st) rbox(l) = rbulk_cwaer(lyyy,2)
1215 end if
1216 end do
1217 end do
1218
1219 end do
1220 goto 7900
1221 end if
1222
1223
1224 do 3900 lyyy = 1, nyyy
1225
1226 !
1227 ! distribution is simple when del_rbulk_cwaer(lyyy) >= 0
1228 !
1229 if (del_rbulk_cwaer(lyyy) .eq. 0.0) then
1230 goto 3900
1231 else if (del_rbulk_cwaer(lyyy) .gt. 0.0) then
1232 do itype = 1, ntype_aer
1233 do isize = 1, nsize_aer(itype)
1234 fr = fr_partit_cw(isize,itype)
1235 if (fr .gt. 0.0) then
1236 l = lptr_yyy_cwaer(isize,itype,lyyy)
1237 if (l .ge. p1st) then
1238 rbox(l) = rbox(l) + fr*del_rbulk_cwaer(lyyy)
1239 end if
1240 end if
1241 end do
1242 end do
1243
1244 goto 3900
1245 end if
1246
1247 !
1248 ! distribution is complicated when del_rbulk_cwaer(lyyy) < 0,
1249 ! because you cannot produce any negative mixrats
1250 !
1251 del_r_remain = del_rbulk_cwaer(lyyy)
1252 fr_cur(:,:) = fr_partit_cw(:,:)
1253
1254 ncount = max( 1, jpos_sv*2 )
1255 icount = 0
1256
1257 ! iteration loop
1258 do while (icount .le. ncount)
1259
1260 icount = icount + 1
1261 del_r_current = del_r_remain
1262 jpos = 0
1263 frsum_cur = 0.0
1264
1265 do itype = 1, ntype_aer
1266 do isize = 1, nsize_aer(itype)
1267
1268 fr = fr_cur(isize,itype)
1269
1270 if (fr .gt. 0.0) then
1271 l = lptr_yyy_cwaer(isize,itype,lyyy)
1272 if (l .ge. p1st) then
1273 duma = fr*del_r_current
1274 dumb = rbox(l) + duma
1275 if (dumb .gt. 0.0) then
1276 jpos = jpos + 1
1277 else if (dumb .eq. 0.0) then
1278 fr_cur(isize,itype) = 0.0
1279 else
1280 duma = -rbox(l)
1281 dumb = 0.0
1282 fr_cur(isize,itype) = 0.0
1283 end if
1284 del_r_remain = del_r_remain - duma
1285 rbox(l) = dumb
1286 frsum_cur = frsum_cur + fr_cur(isize,itype)
1287 else
1288 fr_cur(isize,itype) = 0.0
1289 end if
1290 end if
1291
1292 end do ! isize = 1, nsize_aer
1293 end do ! itype = 1, ntype_aer
1294
1295 ! done if jpos = jpos_sv, because bins reached zero mixrat
1296 if (jpos .eq. jpos_sv) then
1297 idone = 1
1298 ! del_r_remain starts as negative, so done if non-negative
1299 else if (del_r_remain .ge. 0.0) then
1300 idone = 2
1301 ! del_r_remain starts as negative, so done if non-negative
1302 else if (abs(del_r_remain) .le. 1.0e-7*abs(del_rbulk_cwaer(lyyy))) then
1303 idone = 3
1304 ! done if all bins have fr_cur = 0
1305 else if (frsum_cur .le. 0.0) then
1306 idone = 4
1307 ! same thing basically
1308 else if (jpos .le. 0) then
1309 idone = 5
1310 else
1311 idone = 0
1312 end if
1313
1314 ! check for done, and (conditionally) print message
1315 if (idone .gt. 0) then
1316 lunxxaa = 6
1317 #if defined ( ccboxtest_box_testing_active)
1318 lunxxaa = 86
1319 #endif
1320 if ((lunxxaa .gt. 0) .and. (icount .gt. (1+jpos_sv)/2)) then
1321 write(lunxxaa,9800) &
1322 'distribute_bulk_changes - icount>jpos_sv/2 - i,j,k'
1323 write(lunxxaa,9810) it, jt, kt
1324 write(lunxxaa,9800) 'icase, lyyy, idone, icount, jpos, jpos_sv'
1325 write(lunxxaa,9810) icase, lyyy, idone, icount, jpos, jpos_sv
1326 end if
1327 goto 3900
1328 end if
1329
1330 ! rescale fr_cur for next iteration
1331 fr_cur(:,:) = fr_cur(:,:)/frsum_cur
1332
1333 end do ! while (icount .le. ncount)
1334
1335
1336 ! icount > ncount, so print message
1337 lunxxbb = 6
1338 #if defined ( ccboxtest_box_testing_active)
1339 lunxxbb = 86
1340 #endif
1341 if (lunxxbb .gt. 0) then
1342 write(lunxxbb,9800)
1343 write(lunxxbb,9800) &
1344 'distribute_bulk_changes - icount>ncount - i,j,k'
1345 write(lunxxbb,9810) it, jt, kt
1346 write(lunxxbb,9800) 'icase, lyyy, icount, ncount, jpos_sv, jpos'
1347 write(lunxxbb,9810) icase, lyyy, icount, ncount, jpos_sv, jpos
1348 write(lunxxbb,9800) 'rbulk_cwaer(1), del_rbulk_cwaer, del_r_remain, frsum_cur, (frsum_cur-1.0)'
1349 write(lunxxbb,9820) rbulk_cwaer(lyyy,1), del_rbulk_cwaer(lyyy), &
1350 del_r_remain, frsum_cur, (frsum_cur-1.0)
1351 end if
1352 9800 format( a )
1353 9810 format( 7i10 )
1354 9820 format( 7(1pe10.2) )
1355 9840 format( 2i3, 5(1pe14.6) )
1356
1357
1358 3900 continue
1359
1360 7900 continue
1361
1362
1363 #if defined ( ccboxtest_box_testing_active)
1364 ! diagnostics for testing
1365 lunxx = 88
1366 lunxx = -1
1367 if (lunxx .gt. 0) then
1368 icount = 0
1369 do lyyy = 1, nyyy
1370 duma = del_rbulk_cwaer(lyyy)
1371 if ( abs(duma) .gt. &
1372 max( 1.0e-35, 1.0e-5*abs(rbulk_cwaer(lyyy,1)) ) ) then
1373 icount = icount + 1
1374 if (icount .eq. 1) write(lunxx,9800)
1375 if (icount .eq. 1) write(lunxx,9800)
1376 write(lunxx,9800)
1377 write(lunxx,9801) 'distribute_bulk_changes - ',name(lptr_yyy_cwaer(1,1,lyyy)), ' - icase, lyyy'
1378 write(lunxx,9810) icase, lyyy
1379 write(lunxx,9800) ' tp sz rbox_sv1, rbox, del_rbox, del_rbox/del_rbulk_cwaer, (...-fr_partit_cw)'
1380 do itype = 1, ntype_aer
1381 do isize = 1, nsize_aer(itype)
1382 l = lptr_yyy_cwaer(isize,itype,lyyy)
1383 dumb = rbox(l) - rbox_sv1(l)
1384 dumc = dumb/max( abs(duma), 1.0e-35 )
1385 if (duma .lt. 0.0) dumc = -dumc
1386 write(lunxx,9840) itype, isize, rbox_sv1(l), rbox(l), &
1387 dumb, dumc, (dumc-fr_partit_cw(isize,itype))
1388 end do
1389 end do
1390 end if
1391 end do
1392 if (icase .ge. 5) then
1393 write(*,*) &
1394 '*** stop in distribute_bulk_changes diags, icase =', icase
1395 stop
1396 end if
1397 end if
1398 #endif
1399
1400
1401 return
1402 end subroutine distribute_bulk_changes
1403
1404
1405
1406 !-----------------------------------------------------------------------
1407 subroutine cloudchem_apply_move_sections( &
1408 rbox, rbox_sv1, &
1409 it, jt, kt, kpeg, mpeg, icase )
1410
1411 use module_state_description, only: &
1412 param_first_scalar
1413
1414 use module_data_mosaic_asect, only: &
1415 msectional, &
1416 maxd_asize, maxd_atype, &
1417 cw_phase, nsize_aer, ntype_aer, ncomp_aer, &
1418 massptr_aer, numptr_aer, mw_aer, dens_aer, &
1419 lptr_so4_aer, lptr_no3_aer, lptr_cl_aer, lptr_co3_aer, &
1420 lptr_msa_aer, lptr_nh4_aer, lptr_na_aer, lptr_ca_aer, &
1421 lptr_oin_aer, lptr_bc_aer, lptr_oc_aer, &
1422 drymass_aftgrow, drymass_pregrow, &
1423 drydens_aftgrow, drydens_pregrow
1424
1425 use module_data_mosaic_other, only: l2maxd, name, rsub
1426
1427 use module_mosaic_movesect, only: move_sections
1428
1429
1430 implicit none
1431
1432 ! subr arguments
1433 integer, intent(in) :: it, jt, kt, kpeg, mpeg, icase
1434
1435 real, intent(inout), dimension( 1:l2maxd ) :: rbox, rbox_sv1
1436
1437
1438 ! local variables
1439 integer :: idum_msect
1440 integer :: iphase, isize, itype
1441 integer :: l, ll, lunxx
1442 integer :: p1st
1443 integer :: lptr_dum(maxd_asize,maxd_atype)
1444
1445 real :: densdefault
1446 real :: dmaft, dmpre, dvaft, dvpre
1447 real :: duma, dumb, dumc
1448 real :: smallmassbb
1449
1450
1451 p1st = param_first_scalar
1452 iphase = cw_phase
1453
1454 !
1455 ! compute drymass before and after growth
1456 ! set drydens = -1.0, so it will be calculated
1457 !
1458 densdefault = 2.0
1459 smallmassbb = 1.0e-30
1460
1461 do 1800 itype = 1, ntype_aer
1462 do 1800 isize = 1, nsize_aer(itype)
1463 dmaft = 0.0
1464 dmpre = 0.0
1465 dvaft = 0.0
1466 dvpre = 0.0
1467
1468 do ll = 1, ncomp_aer(itype)
1469 l = massptr_aer(ll,isize,itype,iphase)
1470 if (l .ge. p1st) then
1471 duma = mw_aer(ll,itype)
1472 dmaft = dmaft + duma*rbox(l)
1473 dmpre = dmpre + duma*rbox_sv1(l)
1474
1475 duma = duma/dens_aer(ll,itype)
1476 dvaft = dvaft + duma*rbox(l)
1477 dvpre = dvpre + duma*rbox_sv1(l)
1478 end if
1479 end do
1480
1481 drymass_aftgrow(isize,itype) = dmaft
1482 drymass_pregrow(isize,itype) = dmpre
1483
1484 if (min(dmaft,dvaft) .le. smallmassbb) then
1485 drydens_aftgrow(isize,itype) = densdefault
1486 else
1487 drydens_aftgrow(isize,itype) = dmaft/dvaft
1488 end if
1489 if (min(dmpre,dvpre) .le. smallmassbb) then
1490 drydens_pregrow(isize,itype) = densdefault
1491 else
1492 drydens_pregrow(isize,itype) = dmpre/dvpre
1493 end if
1494
1495 1800 continue
1496
1497
1498 ! apply move sections routine
1499 ! (and conditionally turn on movesect diagnostics)
1500 idum_msect = msectional
1501
1502 #if defined ( ccboxtest_box_testing_active )
1503 lunxx = 88
1504 lunxx = -1
1505 if (lunxx .gt. 0) then
1506 if (msectional .lt. 7000) msectional = msectional + 7000
1507 end if
1508 #endif
1509
1510 call move_sections( 2, it, jt, kpeg, mpeg )
1511
1512 msectional = idum_msect
1513
1514
1515 #if defined ( ccboxtest_box_testing_active)
1516 ! diagnostics for testing
1517 if (lunxx .gt. 0) then
1518 do ll = 1, 4
1519 if (ll .eq. 1) then
1520 lptr_dum(:,:) = lptr_so4_aer(:,:,iphase)
1521 else if (ll .eq. 2) then
1522 lptr_dum(:,:) = lptr_nh4_aer(:,:,iphase)
1523 else if (ll .eq. 3) then
1524 lptr_dum(:,:) = lptr_oin_aer(:,:,iphase)
1525 else if (ll .eq. 4) then
1526 lptr_dum(:,:) = numptr_aer(:,:,iphase)
1527 end if
1528
1529 if (ll .eq. 1) write(lunxx,9800)
1530 if (ll .eq. 1) write(lunxx,9800)
1531 write(lunxx,9800)
1532 write(lunxx,9800)
1533 write(lunxx,9801) 'cloudchem_apply_move_sections - ', &
1534 name(lptr_dum(1,1)), ' - icase, ll'
1535 write(lunxx,9810) icase, ll
1536 write(lunxx,9800) ' tp sz rbox_sv1, rbox, rsub'
1537
1538 do itype = 1, ntype_aer
1539 do isize = 1, nsize_aer(itype)
1540 l = lptr_dum(isize,itype)
1541 dumb = rbox(l) - rbox_sv1(l)
1542 dumc = dumb/max( abs(duma), 1.0e-35 )
1543 if (duma .lt. 0.0) dumc = -dumc
1544 write(lunxx,9840) itype, isize, rbox_sv1(l), rbox(l), &
1545 rsub(l,kpeg,mpeg)
1546 end do
1547 end do
1548 end do
1549
1550 if (icase .ge. 5) then
1551 write(*,*) &
1552 '*** stop in cloudchem_apply_move_sections diags, icase =', &
1553 icase
1554 stop
1555 end if
1556 end if
1557 9800 format( a )
1558 9801 format( 3a )
1559 9810 format( 7i10 )
1560 9820 format( 7(1pe10.2) )
1561 9840 format( 2i3, 5(1pe14.6) )
1562 #endif
1563
1564
1565 return
1566 end subroutine cloudchem_apply_move_sections
1567
1568
1569
1570 !-----------------------------------------------------------------------
1571 subroutine mosaic_cloudchem_dumpaa( &
1572 id, ktau, ktauc, dtstepc, config_flags, &
1573 p_phy, t_phy, rho_phy, alt, &
1574 cldfra, ph_no2, &
1575 moist, chem, &
1576 gas_aqfrac, numgas_aqfrac, &
1577 ids,ide, jds,jde, kds,kde, &
1578 ims,ime, jms,jme, kms,kme, &
1579 its,ite, jts,jte, kts,kte, &
1580 qcldwtr_cutoff, &
1581 itcur, jtcur, ktcur )
1582
1583 use module_state_description, only: &
1584 num_moist, num_chem, p_qc
1585
1586 use module_configure, only: grid_config_rec_type
1587
1588 use module_data_mosaic_asect
1589
1590 use module_data_mosaic_other, only: k_pegbegin, name
1591
1592 use module_mosaic_driver, only: mapaer_tofrom_host
1593
1594
1595 implicit none
1596
1597 ! subr arguments
1598 integer, intent(in) :: &
1599 id, ktau, ktauc, &
1600 numgas_aqfrac, &
1601 ids, ide, jds, jde, kds, kde, &
1602 ims, ime, jms, jme, kms, kme, &
1603 its, ite, jts, jte, kts, kte, &
1604 itcur, jtcur, ktcur
1605 ! id - domain index
1606 ! ktau - time step number
1607 ! ktauc - gas and aerosol chemistry time step number
1608 ! numgas_aqfrac - last dimension of gas_aqfrac
1609
1610 ! [ids:ide, kds:kde, jds:jde] - spatial (x,z,y) indices for 'domain'
1611 ! [ims:ime, kms:kme, jms:jme] - spatial (x,z,y) indices for 'memory'
1612 ! Most arrays that are arguments to chem_driver
1613 ! are dimensioned with these spatial indices.
1614 ! [its:ite, kts:kte, jts:jte] - spatial (x,z,y) indices for 'tile'
1615 ! chem_driver and routines under it do calculations
1616 ! over these spatial indices.
1617
1618 type(grid_config_rec_type), intent(in) :: config_flags
1619 ! config_flags - configuration and control parameters
1620
1621 real, intent(in) :: &
1622 dtstepc, qcldwtr_cutoff
1623 ! dtstepc - time step for gas and aerosol chemistry(s)
1624
1625 real, intent(in), &
1626 dimension( ims:ime, kms:kme, jms:jme ) :: &
1627 p_phy, t_phy, rho_phy, alt, cldfra, ph_no2
1628 ! p_phy - air pressure (Pa)
1629 ! t_phy - temperature (K)
1630 ! rho_phy - moist air density (kg/m^3)
1631 ! alt - dry air specific volume (m^3/kg)
1632 ! cldfra - cloud fractional area (0-1)
1633 ! ph_no2 - no2 photolysis rate (1/min)
1634
1635 real, intent(in), &
1636 dimension( ims:ime, kms:kme, jms:jme, 1:num_moist ) :: &
1637 moist
1638 ! moist - mixing ratios of moisture species (water vapor,
1639 ! cloud water, ...) (kg/kg for mass species, #/kg for number species)
1640
1641 real, intent(inout), &
1642 dimension( ims:ime, kms:kme, jms:jme, 1:num_chem ) :: &
1643 chem
1644 ! chem - mixing ratios of trace gas and aerosol species (ppm for gases,
1645 ! ug/kg for aerosol mass species, #/kg for aerosol number species)
1646
1647 real, intent(inout), &
1648 dimension( ims:ime, kms:kme, jms:jme, numgas_aqfrac ) :: &
1649 gas_aqfrac
1650 ! gas_aqfrac - fraction (0-1) of gas that is dissolved in cloud water
1651
1652
1653 ! local variables
1654 integer :: it, jt, kt, l, ll, n
1655 integer :: isize, itype
1656
1657 real :: dumai, dumcw
1658 real :: qcldwtr
1659
1660
1661 it = itcur
1662 jt = jtcur
1663 kt = ktcur
1664
1665 write(*,*)
1666 write(*,*)
1667 write(*,*)
1668 write(*,9100)
1669 write(*,9102) ktau, it, jt, kt
1670 9100 format( 7('----------') )
1671 9102 format( &
1672 'mosaic_cloudchem_dumpaa - ktau, i, j, k =', 4i5 )
1673
1674 itype = 1
1675 do 2900 isize = 1, nsize_aer(itype)
1676
1677 write(*,9110) isize
1678 9110 format( / 'isize =', i3 / &
1679 ' k cldwtr mass-ai numb-ai mass-cw numb-cw' )
1680
1681 do 2800 kt = kte-1, kts, -1
1682
1683 dumai = 0.0
1684 dumcw = 0.0
1685 do ll = 1, ncomp_aer(itype)
1686 l = massptr_aer(ll,isize,itype,1)
1687 dumai = dumai + chem(it,kt,jt,l)
1688 l = massptr_aer(ll,isize,itype,2)
1689 dumcw = dumcw + chem(it,kt,jt,l)
1690 end do
1691 write(*,9120) kt, &
1692 moist(it,kt,jt,p_qc), &
1693 dumai, chem(it,kt,jt,numptr_aer(isize,itype,1)), &
1694 dumcw, chem(it,kt,jt,numptr_aer(isize,itype,2))
1695 9120 format( i3, 1p, e10.2, 2(3x, 2e10.2) )
1696
1697 2800 continue
1698 2900 continue
1699
1700 write(*,*)
1701 write(*,9100)
1702 write(*,*)
1703
1704 ! map from wrf-chem 3d arrays to pegasus clm & sub arrays
1705 kt = ktcur
1706 if ((ktau .eq. 30) .and. (it .eq. 23) .and. &
1707 (jt .eq. 1) .and. (kt .eq. 11)) then
1708 qcldwtr = moist(it,kt,jt,p_qc)
1709 write(*,*)
1710 write(*,*)
1711 write(*,9102) ktau, it, jt, kt
1712 write(*,*)
1713 write( *, '(3(1pe10.2,3x,a))' ) &
1714 (chem(it,kt,jt,l), name(l)(1:10), l=1,num_chem)
1715 write(*,*)
1716 write( *, '(3(1pe10.2,3x,a))' ) &
1717 p_phy(it,kt,jt), 'p_phy ', &
1718 t_phy(it,kt,jt), 't_phy ', &
1719 rho_phy(it,kt,jt), 'rho_phy ', &
1720 alt(it,kt,jt), 'alt ', &
1721 qcldwtr, 'qcldwtr ', &
1722 qcldwtr_cutoff, 'qcldwtrcut'
1723 write(*,*)
1724 write(*,9100)
1725 write(*,*)
1726 end if
1727
1728
1729 return
1730 end subroutine mosaic_cloudchem_dumpaa
1731
1732
1733
1734 !-----------------------------------------------------------------------
1735 end module module_mosaic_cloudchem