| blob | 453ac3e0429c254fc1ce6c70c312a0eedc90112a |
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 references and terms of use
8 !**********************************************************************************
9 !CPP directives to control ic/bc conditions...
10 !(The directive in module_input_chem_data also needs to be set.)
11 ! CASENAME = 0 Uses Texas August 2004 case values and profiles
12 ! 1 Uses same concentrations as TX, but uses different
13 ! profiles depending on the species. (NEAQS2004 case)
14 ! 2 Mexico City
15 ! 3 Mexico City for Testbed, to be consistent with MADE/SORGAM too
16 ! 4 bc,ic values from wrfinput and wrfbdy files
17 #define CASENAME 4
19 module module_mosaic_initmixrats
21 ! rce 2005-feb-18 - several fixes for indices of dlo/hi_sect, volumlo/hi_sect,
22 ! which are now (isize,itype)
24 ! rce 2004-dec-03 - many changes associated with the new aerosol "pointer"
25 ! variables in module_data_mosaic_asect
27 USE module_state_description
29 integer, parameter :: mosaic_init_wrf_mixrats_flagaa = 1
30 ! turns subr mosaic_init_wrf_mixrats on/off
32 contains
35 !-----------------------------------------------------------------------
36 subroutine mosaic_init_wrf_mixrats( &
37 iflagaa, config_flags, &
38 chem, alt, z_at_w, g, &
39 ids,ide, jds,jde, kds,kde, &
40 ims,ime, jms,jme, kms,kme, &
41 its,ite, jts,jte, kts,kte )
43 !
44 ! initializes the species and number mixing ratios for each section
45 !
46 ! this top level routine simply calls other routines depending on value
47 ! of config_flags%aer_ic_opt
48 !
49 use module_configure, only: grid_config_rec_type
50 use module_state_description, only: num_chem, param_first_scalar, &
51 aer_ic_pnnl
52 use module_data_mosaic_asect
53 use module_data_mosaic_other
54 use module_peg_util, only: peg_message, peg_error_fatal
56 implicit none
59 ! subr arguments
60 type(grid_config_rec_type), intent(in) :: config_flags
62 integer, intent(in) :: &
63 iflagaa, &
64 ids, ide, jds, jde, kds, kde, &
65 ims, ime, jms, jme, kms, kme, &
66 its, ite, jts, jte, kts, kte
68 real, intent(inout), &
69 dimension( ims:ime, kms:kme, jms:jme, 1:num_chem ) :: &
70 chem
72 real, intent(in), &
73 dimension( ims:ime, kms:kme, jms:jme ) :: &
74 alt, z_at_w
75 real, intent(in) :: g
77 ! local variables
78 integer :: i, ic, j, k
80 #if (CASENAME == 0) || (CASENAME == 1) || (CASENAME == 2) || (CASENAME == 3)
81 if (config_flags%aer_ic_opt == aer_ic_pnnl) then
82 call mosaic_init_wrf_mixrats_opt2( &
83 iflagaa, config_flags, &
84 chem, z_at_w, g, &
85 ids,ide, jds,jde, kds,kde, &
86 ims,ime, jms,jme, kms,kme, &
87 its,ite, jts,jte, kts,kte )
89 else
90 call mosaic_init_wrf_mixrats_opt1( &
91 iflagaa, config_flags, &
92 chem, &
93 ids,ide, jds,jde, kds,kde, &
94 ims,ime, jms,jme, kms,kme, &
95 its,ite, jts,jte, kts,kte )
97 end if
99 ! Aerosol species are returned from above in concentration units. Convert
100 ! them to mixing ratio for use in advection, etc.
101 do ic = p_so4_a01,num_chem
102 do j = jts,jte
103 do k = kts,kte-1
104 do i = its,ite
105 chem(i,k,j,ic) = chem(i,k,j,ic)*alt(i,k,j)
106 end do
107 end do
108 end do
109 end do
111 ! Fill the top z-staggered location to prevent trouble during advection.
112 do ic = p_so4_a01,num_chem
113 do j = jts,jte
114 do i = its,ite
115 chem(i,kte,j,ic) = chem(i,kte-1,j,ic)
116 end do
117 end do
118 end do
119 #endif
121 return
122 end subroutine mosaic_init_wrf_mixrats
125 !-----------------------------------------------------------------------
126 subroutine mosaic_init_wrf_mixrats_opt1( &
127 iflagaa, config_flags, &
128 chem, &
129 ids,ide, jds,jde, kds,kde, &
130 ims,ime, jms,jme, kms,kme, &
131 its,ite, jts,jte, kts,kte )
133 !
134 ! initializes the species and number mixing ratios for each section
135 ! based on user-specified lognormal modes that span the size distribution
136 !
137 ! rce 11-sep-2004 - eliminated use of the _wrfch pointers (lptr_xxx_a_wrfch,
138 ! lwaterptr_wrfch, numptr_wrfch); use only the _amode pointers now;
139 ! added l1dum
140 !
141 use module_configure, only: grid_config_rec_type
142 use module_state_description, only: num_chem, param_first_scalar
143 use module_data_mosaic_asect
144 use module_data_mosaic_other
145 use module_peg_util, only: peg_message, peg_error_fatal
147 implicit none
149 ! subr arguments
150 type(grid_config_rec_type), intent(in) :: config_flags
152 integer, intent(in) :: &
153 iflagaa, &
154 ids, ide, jds, jde, kds, kde, &
155 ims, ime, jms, jme, kms, kme, &
156 its, ite, jts, jte, kts, kte
158 real, intent(inout), &
159 dimension( ims:ime, kms:kme, jms:jme, 1:num_chem ) :: &
160 chem
162 ! local variables
163 integer i, j, k, l, ll, l1, l3, l1dum, m, mdum, nsm
164 integer it, jt, kt
166 real dum, dumdp, dumrsfc, dumvol, &
167 xlo, xhi, &
168 dumvol1p, &
169 pdummb, zdumkm, zscalekm, zfactor
171 real vtot_nsm_ofmode(maxd_asize)
172 real dumarr(maxd_acomp+5)
174 real erfc
176 ! double precision fracnum, fracvol, tlo, thi
177 real fracvol, tlo, thi
179 integer nmaxd_nsm
180 parameter (nmaxd_nsm = 4)
182 integer iphase, itype, ntot_nsm
183 integer iiprof_nsm(nmaxd_nsm)
184 integer lldum_so4, lldum_nh4, lldum_oc, lldum_bc, &
185 lldum_oin, lldum_na, lldum_cl, lldum_hysw
187 real sx_nsm(nmaxd_nsm), sxr2_nsm(nmaxd_nsm), &
188 x0_nsm(nmaxd_nsm), x3_nsm(nmaxd_nsm), &
189 rtot_nsm(maxd_acomp,nmaxd_nsm), &
190 vtot_nsm(nmaxd_nsm), xntot_nsm(nmaxd_nsm)
192 real dgnum_nsm(nmaxd_nsm), sigmag_nsm(nmaxd_nsm)
193 real aaprof_nsm(maxd_acomp+1,nmaxd_nsm)
194 real bbprof_nsm(nmaxd_nsm)
196 character*80 msg
197 character*10 dumname
200 ! check for on/off
201 if (mosaic_init_wrf_mixrats_flagaa .le. 0) return
204 ! *** currently only works for ntype_aer = 1
205 itype = 1
206 iphase = ai_phase
207 m = 1
209 ! set values for initialization modes
210 iiprof_nsm(:) = 1
211 aaprof_nsm(:,:) = 0.0
212 bbprof_nsm(:) = 0.0
214 ntot_nsm = 4
215 ntot_nsm = min( ntot_nsm, nsize_aer(itype) )
217 lldum_so4 = 0
218 lldum_nh4 = 0
219 lldum_oc = 0
220 lldum_bc = 0
221 lldum_oin = 0
222 lldum_na = 0
223 lldum_cl = 0
224 lldum_hysw = 0
226 do ll = 1, ncomp_plustracer_aer(itype)
227 if (massptr_aer(ll,m,itype,iphase) .eq. &
228 lptr_so4_aer(m,itype,iphase)) lldum_so4 = ll
229 if (massptr_aer(ll,m,itype,iphase) .eq. &
230 lptr_nh4_aer(m,itype,iphase)) lldum_nh4 = ll
231 if (massptr_aer(ll,m,itype,iphase) .eq. &
232 lptr_oc_aer(m,itype,iphase)) lldum_oc = ll
233 if (massptr_aer(ll,m,itype,iphase) .eq. &
234 lptr_bc_aer(m,itype,iphase)) lldum_bc = ll
235 if (massptr_aer(ll,m,itype,iphase) .eq. &
236 lptr_oin_aer(m,itype,iphase)) lldum_oin = ll
237 if (massptr_aer(ll,m,itype,iphase) .eq. &
238 lptr_na_aer(m,itype,iphase)) lldum_na = ll
239 if (massptr_aer(ll,m,itype,iphase) .eq. &
240 lptr_cl_aer(m,itype,iphase)) lldum_cl = ll
241 end do
242 if (hyswptr_aer(m,itype) .gt. 0) &
243 lldum_hysw = ncomp_plustracer_aer(itype) + 1
245 msg = ' '
246 if (lldum_so4 .le. 0) &
247 msg = '*** subr mosaic_init_wrf_mixrats - lldum_so4 = 0'
248 if (lldum_nh4 .le. 0) &
249 msg = '*** subr mosaic_init_wrf_mixrats - lldum_nh4 = 0'
250 if (lldum_oc .le. 0) &
251 msg = '*** subr mosaic_init_wrf_mixrats - lldum_oc = 0'
252 if (lldum_bc .le. 0) &
253 msg = '*** subr mosaic_init_wrf_mixrats - lldum_bc = 0'
254 if (lldum_oin .le. 0) &
255 msg = '*** subr mosaic_init_wrf_mixrats - lldum_oin = 0'
256 if (lldum_na .le. 0) &
257 msg = '*** subr mosaic_init_wrf_mixrats - lldum_na = 0'
258 if (lldum_cl .le. 0) &
259 msg = '*** subr mosaic_init_wrf_mixrats - lldum_cl = 0'
260 if (lldum_hysw .le. 0) &
261 msg = '*** subr mosaic_init_wrf_mixrats - lldum_hysw = 0'
262 if (msg .ne. ' ') call peg_error_fatal( lunerr, msg )
265 do nsm = 1, ntot_nsm
267 if (nsm .eq. 1) then
268 ! accum mode with so4, nh4, oc, bc
269 dgnum_nsm( nsm) = 0.15e-4
270 sigmag_nsm(nsm) = 2.0
271 aaprof_nsm(lldum_so4,nsm) = 0.50
272 aaprof_nsm(lldum_nh4,nsm) = aaprof_nsm(lldum_so4,nsm) &
273 * (mw_nh4_aer/mw_so4_aer)
274 aaprof_nsm(lldum_oc,nsm) = 0.25
275 aaprof_nsm(lldum_bc,nsm) = 0.05
276 aaprof_nsm(lldum_hysw,nsm) = aaprof_nsm(lldum_so4,nsm) * 0.2
278 else if (nsm .eq. 2) then
279 ! aitken mode with so4, nh4, oc, bc
280 dgnum_nsm( nsm) = 0.03e-4
281 sigmag_nsm(nsm) = 2.0
282 aaprof_nsm(lldum_so4,nsm) = 0.50 * 0.020
283 aaprof_nsm(lldum_nh4,nsm) = aaprof_nsm(lldum_so4,nsm) &
284 * (mw_nh4_aer/mw_so4_aer)
285 aaprof_nsm(lldum_oc,nsm) = 0.25 * 0.020
286 aaprof_nsm(lldum_bc,nsm) = 0.05 * 0.020
287 aaprof_nsm(lldum_hysw,nsm) = aaprof_nsm(lldum_so4,nsm) * 0.2
289 else if (nsm .eq. 3) then
290 ! coarse dust mode with oin
291 dgnum_nsm( nsm) = 1.0e-4
292 sigmag_nsm(nsm) = 2.0
293 aaprof_nsm(lldum_oin,nsm) = 0.5
294 aaprof_nsm(lldum_hysw,nsm) = aaprof_nsm( 9,nsm) * 1.0e-3
296 else if (nsm .eq. 4) then
297 ! coarse seasalt mode with na, cl
298 dgnum_nsm( nsm) = 2.0e-4
299 sigmag_nsm(nsm) = 2.0
300 aaprof_nsm(lldum_cl,nsm) = 0.1
301 aaprof_nsm(lldum_na,nsm) = aaprof_nsm(lldum_cl,nsm) &
302 * (mw_na_aer/mw_cl_aer)
303 aaprof_nsm(lldum_hysw,nsm) = aaprof_nsm(lldum_cl,nsm) * 0.2
305 end if
307 end do
309 ! when iflagaa = 1/2/3/4, use only the nsm-mode = iflagaa
310 if (iflagaa .gt. 0) then
311 do nsm = 1, ntot_nsm
312 if (nsm .ne. iflagaa) aaprof_nsm(:,nsm) = 0.0
313 end do
314 end if
318 !
319 ! do the initialization now
320 !
322 ! calculate mode parameters
323 do nsm = 1, ntot_nsm
324 sx_nsm(nsm) = alog( sigmag_nsm(nsm) )
325 sxr2_nsm(nsm) = sx_nsm(nsm) * sqrt(2.0)
326 x0_nsm(nsm) = alog( dgnum_nsm(nsm) )
327 x3_nsm(nsm) = x0_nsm(nsm) + 3.0*sx_nsm(nsm)*sx_nsm(nsm)
328 end do
330 ! initialize rclm array to zero
331 rclm(:,:) = 0.
332 ! rclmsvaa(:,:) = 0.
334 !
335 ! loop over all vertical levels
336 !
337 ! do 12900 k = 1, ktot
338 do 12900 k = 1, 1
340 ! pdummb = 1013.*scord(k)
341 ! zdumkm = ptoz( pdummb ) * 1.e-3
342 zdumkm = 0.0
345 ! for each species and nsm mode, define total mixing ratio
346 ! (for all sizes) at level k
347 !
348 ! iiprof_nsm = +1 gives constant mixing ratio
349 ! aaprof_nsm(l,nsm) = constant mixing ratio (ppb)
350 ! iiprof_nsm = +2 gives exponential profile
351 ! aaprof_nsm(l,nsm) = surface mixing ratio (ppb)
352 ! bbprof_nsm(l) = scale height (km)
353 ! iiprof_nsm = +3 gives linear profile (then zero at z > zscalekm)
354 ! aaprof_nsm(l,nsm) = surface mixing ratio (ppb)
355 ! bbprof_nsm(l) = height (km) at which mixing ratio = 0
357 do nsm = 1, ntot_nsm
359 if (iiprof_nsm(nsm) .eq. 2) then
360 zscalekm = bbprof_nsm(nsm)
361 zfactor = exp( -zdumkm/zscalekm )
362 else if (iiprof_nsm(nsm) .eq. 3) then
363 zscalekm = bbprof_nsm(nsm)
364 zfactor = max( 0., (1. - zdumkm/zscalekm) )
365 else
366 zfactor = 1.0
367 end if
369 do ll = 1, ncomp_plustracer_aer(itype) + 1
370 rtot_nsm(ll,nsm) = aaprof_nsm(ll,nsm) * zfactor
371 end do
373 end do
375 ! compute total volume and number mixing ratios for each nsm mode
376 ! rtot_nsm is ug/m3, 1.0e-6*rtot is g/m3, vtot_nsm is cm3/m3
377 do nsm = 1, ntot_nsm
378 dumvol = 0.
379 do ll = 1, ncomp_aer(itype)
380 dum = 1.0e-6*rtot_nsm(ll,nsm)/dens_aer(ll,itype)
381 dumvol = dumvol + max( 0., dum )
382 end do
383 vtot_nsm(nsm) = dumvol
384 end do
386 ! now compute species and number mixing ratios for each bin
387 do 12700 m = 1, nsize_aer(itype)
389 vtot_nsm_ofmode(m) = 0.0
391 do 12500 nsm = 1, ntot_nsm
393 ! for nsm_mode = n, compute fraction of number and volume
394 ! that is in bin m
395 xlo = alog( dlo_sect(m,itype) )
396 xhi = alog( dhi_sect(m,itype) )
398 tlo = (xlo - x3_nsm(nsm))/sxr2_nsm(nsm)
399 thi = (xhi - x3_nsm(nsm))/sxr2_nsm(nsm)
400 if (tlo .ge. 0.) then
401 ! fracvol = 0.5*( erfc(tlo) - erfc(thi) )
402 fracvol = 0.5*( erfc_num_recipes(tlo) - erfc_num_recipes(thi) )
403 else
404 ! fracvol = 0.5*( erfc(-thi) - erfc(-tlo) )
405 fracvol = 0.5*( erfc_num_recipes(-thi) - erfc_num_recipes(-tlo) )
406 end if
407 fracvol = max( fracvol, 0.0 )
409 vtot_nsm_ofmode(m) = vtot_nsm_ofmode(m) + vtot_nsm(nsm)*fracvol
411 ! now load that fraction of species-mixing-ratio
412 ! into the appropriate rclm location
413 do ll = 1, ncomp_plustracer_aer(itype)
414 rclm( k, massptr_aer(ll,m,itype,iphase) ) = &
415 rclm( k, massptr_aer(ll,m,itype,iphase) ) + &
416 fracvol*rtot_nsm(ll,nsm)
417 end do
419 if ((iphase .eq. ai_phase) .and. &
420 (lldum_hysw .gt. 0) .and. &
421 (hyswptr_aer(m,itype) .gt. 0)) then
423 rclm( k, hyswptr_aer(m,itype) ) = &
424 rclm( k, hyswptr_aer(m,itype) ) + &
425 fracvol*rtot_nsm(lldum_hysw,nsm)
426 end if
428 12500 continue
430 ! now compute number from volume
431 dum = sqrt( dlo_sect(m,itype)*dhi_sect(m,itype) )
432 dumvol1p = (pi/6.0)*(dum**3)
433 rclm( k, numptr_aer(m,itype,iphase) ) = vtot_nsm_ofmode(m)/dumvol1p
435 ! set water = hyswatr
436 if ((iphase .eq. ai_phase) .and. &
437 (lldum_hysw .gt. 0) .and. &
438 (hyswptr_aer(m,itype) .gt. 0) .and. &
439 (waterptr_aer(m,itype) .gt. 0)) then
441 rclm( k, waterptr_aer(m,itype) ) = &
442 rclm( k, hyswptr_aer(m,itype) )
443 end if
445 12700 continue
447 12900 continue
450 !
451 ! do diagnostic output
452 !
454 ! temporary
455 ! temporary
457 dumarr(:) = 0.0
458 msg = ' '
459 call peg_message( lunout, msg )
460 msg = '*** subr mosaic_init_wrf_mixrats_opt1 results'
461 call peg_message( lunout, msg )
462 msg = ' mass in ug/m3 number in #/m3 volume in cm3/m3'
463 call peg_message( lunout, msg )
464 msg = ' '
465 call peg_message( lunout, msg )
466 msg = ' mode l l1 species conc'
467 call peg_message( lunout, msg )
469 do 14390 mdum = 1, nsize_aer(itype)+1
470 m = min( mdum, nsize_aer(itype) )
471 msg = ' '
472 call peg_message( lunout, msg )
473 do 14350 l = 1, ncomp_plustracer_aer(itype)+4
475 if (l .le. ncomp_plustracer_aer(itype)) then
476 l1 = massptr_aer(l,m,itype,iphase)
477 dumname = name_aer(l,itype)
478 dum = rclm(1,l1)
479 else if (l .eq. ncomp_plustracer_aer(itype)+1) then
480 l1 = hyswptr_aer(m,itype)
481 dumname = 'hystwatr'
482 dum = rclm(1,l1)
483 else if (l .eq. ncomp_plustracer_aer(itype)+2) then
484 l1 = waterptr_aer(m,itype)
485 dumname = 'water'
486 dum = rclm(1,l1)
487 else if (l .eq. ncomp_plustracer_aer(itype)+3) then
488 l1 = numptr_aer(m,itype,iphase)
489 dumname = 'number'
490 dum = rclm(1,l1)
491 else if (l .eq. ncomp_plustracer_aer(itype)+4) then
492 l1 = 0
493 dumname = 'volume'
494 dum = vtot_nsm_ofmode(m)
495 else
496 dumname = '=BADBAD='
497 l1 = -1
498 dum = -1.0
499 end if
501 l1dum = l1
502 if (aboxtest_testmode .gt. 0) l1dum = max( l1-1, 0 )
504 if (mdum .le. nsize_aer(itype)) then
505 dumarr(l) = dumarr(l) + dum
506 write(msg,9620) m, l, l1dum, dumname, dum
507 else
508 write(msg,9625) l, dumname, dumarr(l)
509 end if
510 call peg_message( lunout, msg )
512 14350 continue
513 14390 continue
515 9620 format( 3i4, 2x, a, 3(1pe12.3) )
516 9625 format( ' sum', i4, ' -', 2x, a, 3(1pe12.3) )
519 !
520 ! load the chem array
521 !
522 do 16390 m = 1, nsize_aer(itype)
523 do 16350 l = 1, 15
525 if (l .eq. 1) then
526 l1 = lptr_so4_aer(m,itype,iphase)
527 else if (l .eq. 2) then
528 l1 = lptr_no3_aer(m,itype,iphase)
529 else if (l .eq. 3) then
530 l1 = lptr_cl_aer(m,itype,iphase)
531 else if (l .eq. 4) then
532 l1 = lptr_msa_aer(m,itype,iphase)
533 else if (l .eq. 5) then
534 l1 = lptr_co3_aer(m,itype,iphase)
535 else if (l .eq. 6) then
536 l1 = lptr_nh4_aer(m,itype,iphase)
537 else if (l .eq. 7) then
538 l1 = lptr_na_aer(m,itype,iphase)
539 else if (l .eq. 8) then
540 l1 = lptr_ca_aer(m,itype,iphase)
541 else if (l .eq. 9) then
542 l1 = lptr_oin_aer(m,itype,iphase)
543 else if (l .eq. 10) then
544 l1 = lptr_oc_aer(m,itype,iphase)
545 else if (l .eq. 11) then
546 l1 = lptr_bc_aer(m,itype,iphase)
547 else if (l .eq. 12) then
548 l1 = hyswptr_aer(m,itype)
549 else if (l .eq. 13) then
550 l1 = waterptr_aer(m,itype)
551 else if (l .eq. 14) then
552 l1 = numptr_aer(m,itype,iphase)
553 else
554 goto 16350
555 end if
556 l3 = l1
558 if ((l1 .gt. 0) .and. (l1 .le. lmaxd) .and. &
559 (l3 .ge. param_first_scalar)) then
560 do it = its, ite
561 ! *** note: not sure what the kt limits should be here
562 do kt = kts, kte-1
563 do jt = jts, jte
564 chem(it,kt,jt,l3) = rclm(1,l1)
565 end do
566 end do
567 end do
568 end if
570 16350 continue
571 16390 continue
574 ! all done
575 return
576 end subroutine mosaic_init_wrf_mixrats_opt1
579 !-----------------------------------------------------------------------
580 !wig 10-May-2004, added phb, ph, g
581 subroutine mosaic_init_wrf_mixrats_opt2( &
582 iflagaa, config_flags, &
583 chem, z_at_w, g, &
584 ids,ide, jds,jde, kds,kde, &
585 ims,ime, jms,jme, kms,kme, &
586 its,ite, jts,jte, kts,kte )
588 !
589 ! provides initial values for mosaic aerosol species (mass and number
590 ! mixing ratio) for "Texas August 2000" simulations
591 ! modified to use different profiles for different aerosols for NEAQS case, egc 7/2005
592 ! currently all the initial values are uniform in x, y, AND z
593 !
594 ! rce 11-sep-2004 - eliminated use of the _wrfch pointers (lptr_xxx_a_wrfch,
595 ! lwaterptr_wrfch, numptr_wrfch); use only the _amode pointers now
596 !
597 use module_configure, only: grid_config_rec_type
598 use module_state_description, only: num_chem, param_first_scalar
599 use module_data_mosaic_asect
600 use module_data_mosaic_other
601 use module_peg_util, only: peg_message, peg_error_fatal
603 implicit none
605 ! subr arguments
606 type(grid_config_rec_type), intent(in) :: config_flags
608 integer, intent(in) :: &
609 iflagaa, &
610 ids, ide, jds, jde, kds, kde, &
611 ims, ime, jms, jme, kms, kme, &
612 its, ite, jts, jte, kts, kte
614 real, intent(inout), &
615 dimension( ims:ime, kms:kme, jms:jme, 1:num_chem ) :: &
616 chem
618 real, intent(in), &
619 dimension( ims:ime, kms:kme, jms:jme ) :: &
620 z_at_w
621 real :: g
623 ! local variables
624 integer l, l1, l3, m, mdum
625 integer iphase, itype
626 integer it, jt, kt
628 !wig 10-May-2004, added mult
629 real dum, dumvol1p, mult
630 real qcoar, qfine, qval
632 real :: vtot_ofmode(maxd_asize)
633 real :: dumarr(maxd_acomp+5)
634 real :: fr_coar(8), fr_fine(8)
636 !wig 01-Apr-2005, Updated fractional breakdown between bins. (See also
637 ! bdy_chem_value_mosaic in this file and mosaic_addemiss in
638 ! module_mosaic_addemiss.F) Note that the values here no
639 ! longer match those in mosaic_addemiss.
640 !rce 10-May-2005, changed fr8b_coar to values determined by jdf
641 real, save :: fr8b_coar(8) = &
642 (/ 0.0, 0.0, 0.0, 0.019, 0.0212, 0.1101, 0.2751, 0.7018 /) ! jdf Testbed
643 !jdf (/ 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.300, 0.700 /) ! 10-May-2005
644 ! (/ 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.933, 0.067 /) ! 01-Apr-2005
645 ! (/ 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.16, 0.84 /) ! "old"
647 real, save :: fr8b_fine(8) = &
648 (/ 0.0245, 0.1472, 0.3501, 0.3317, 0.1251, 0.0186, 0.0011, 0.0/) !jdf Testbed
649 !jdf (/ 0.006, 0.024, 0.231, 0.341, 0.140, 0.258, 0., 0./) !5-Apr-2005 values
650 ! (/ 0.0275, 0.0426, 0.2303, 0.3885, 0.1100, 0.2011, 0., 0./) !15-Nov-2004 values
651 ! (/ 0.01, 0.05, 0.145, 0.60, 0.145, 0.05, 0.00, 0.00 /)
652 ! (/ 0.04, 0.10, 0.35, 0.29, 0.15, 0.07, 0.00, 0.00 /)
654 real :: qfine_so4, qfine_no3, qfine_cl, qfine_msa, &
655 qfine_co3, qfine_nh4, qfine_na, qfine_ca, qfine_oin, &
656 qfine_oc, qfine_bc, qfine_hysw, qfine_watr, qfine_vol
657 real :: qcoar_so4, qcoar_no3, qcoar_cl, qcoar_msa, &
658 qcoar_co3, qcoar_nh4, qcoar_na, qcoar_ca, qcoar_oin, &
659 qcoar_oc, qcoar_bc, qcoar_hysw, qcoar_watr, qcoar_vol
661 !wig 10-May-2004, added z
662 real, dimension( ims:ime, kms:kme, jms:jme ) :: z
664 character*80 msg
665 character*10 dumname
668 ! *** currently only works for ntype_aer = 1
669 itype = 1
670 iphase = ai_phase
671 m = 1
673 !wig 10-May-2004, block begin...
674 ! calculate the mid-level heights
675 do jt = jts, min(jte,jde-1)
676 do kt = kts, kte-1
677 do it = its, min(ite,ide-1)
678 z(it,kt,jt) = (z_at_w(it,kt,jt)+z_at_w(it,kt+1,jt))*0.5
679 end do
680 end do
681 end do
682 !wig 10-May-2004, ...end block
684 ! set the partitioning fractions for either 8 or 4 bins
685 if (nsize_aer(itype) == 8) then
686 fr_coar(:) = fr8b_coar(:)
687 fr_fine(:) = fr8b_fine(:)
688 else if (nsize_aer(itype) == 4) then
689 do m = 1, nsize_aer(itype)
690 fr_coar(m) = fr8b_coar(2*m) + fr8b_coar(2*m-1)
691 fr_fine(m) = fr8b_fine(2*m) + fr8b_fine(2*m-1)
692 end do
693 else
694 write(msg,'(a,i5)') &
695 'subr mosaic_init_wrf_mixrats_opt2' // &
696 ' - nsize_aer(itype) must be 4 or 8 but = ', nsize_aer(itype)
697 call peg_error_fatal( lunout, msg )
698 end if
700 !
701 ! compute initial values (currently uniform in x, y, AND z)
702 ! load them into the chem array
703 ! also load them into the rclm array for diagnostic output
704 !
705 rclm(:,:) = 0.0
706 vtot_ofmode(:) = 0.0
708 ! Set values for fine and coarse mass concentrations, and
709 ! compute fine/coarse volume concentrations. These are also set in
710 ! bdy_chem_value_mosaic, below.
711 ! (these are latest values from 1-Apr-2005 discussions)
712 !
713 ! rce 10-may-2005 - changed qfine_cl, _na, _oin to values determined by jdf
714 qfine_so4 = 2.14
715 qcoar_so4 = 0.242
716 qfine_no3 = 0.11
717 qcoar_no3 = 0.03
718 ! qfine_cl = 0.3
719 qfine_cl = 0.14 ! 10-May-2005
720 qcoar_cl = 0.139
721 qfine_msa = 0.0
722 qcoar_msa = 0.0
723 qfine_co3 = 0.0
724 qcoar_co3 = 0.0
725 qfine_nh4 = 0.83
726 qcoar_nh4 = 0.10
727 ! qfine_na = 0.2
728 qfine_na = 0.1 ! 10-May-2005
729 qcoar_na = 0.09
730 qfine_ca = 0.0
731 qcoar_ca = 0.0
732 ! qfine_oin = 6.2
733 qfine_oin = 3.48 ! 10-May-2005
734 qcoar_oin = 0.35
735 qfine_oc = 1.00
736 qcoar_oc = 1.50
737 qfine_bc = 0.2
738 qcoar_bc = 0.075
739 qfine_hysw = 0.0
740 qcoar_hysw = 0.0
741 qfine_watr = 0.0
742 qcoar_watr = 0.0
743 #if (CASENAME == 2)
744 ! qfine_so4 = 0.50
745 ! qcoar_so4 = 0.50/2.0
746 ! qfine_no3 = 0.75
747 ! qcoar_no3 = 0.75/2.0
748 ! qfine_cl = 0.10
749 ! qcoar_cl = 0.10/2.0
750 ! qfine_nh4 = 0.30
751 ! qcoar_nh4 = 0.30/2.0
752 ! qfine_na = 0.20
753 ! qcoar_na = 0.20/2.0
754 ! qfine_oin = 4.00
755 ! qcoar_oin = 4.00/2.0
756 ! qfine_oc = 3.00
757 ! qcoar_oc = 3.00/2.0
758 ! qfine_bc = 0.50
759 ! qcoar_bc = 0.50/2.0
760 qfine_so4 = 0.300
761 qcoar_so4 = 0.300/2.0
762 qfine_no3 = 0.001
763 qcoar_no3 = 0.001/2.0
764 qfine_cl = 1.050
765 qcoar_cl = 1.05/02.0
766 qfine_nh4 = 0.094
767 qcoar_nh4 = 0.094/2.0
768 qfine_na = 0.700
769 qcoar_na = 0.700/2.0
770 qfine_oin = 4.500
771 qcoar_oin = 4.500/2.0
772 qfine_oc = 0.088
773 qcoar_oc = 0.088/2.0
774 qfine_bc = 0.013
775 qcoar_bc = 0.013/2.0
776 #endif
777 #if (CASENAME == 3)
778 qfine_so4 = 0.300
779 qcoar_so4 = 0.0
780 qfine_nh4 = 0.094
781 qcoar_nh4 = 0.0
782 qfine_no3 = 0.001
783 qcoar_no3 = 0.0
784 qfine_cl = 0.0
785 qcoar_cl = 1.050
786 qfine_na = 0.000
787 qcoar_na = 0.700
788 qfine_oin = 4.500
789 qcoar_oin = 4.500/2.0
790 qfine_oc = 0.088
791 qcoar_oc = 0.0
792 qfine_bc = 0.013
793 qcoar_bc = 0.0
794 #endif
796 !!$! old values from 23-Apr-2004:
797 !!$ qfine_so4 = 2.554
798 !!$ qcoar_so4 = 0.242
799 !!$ qfine_no3 = 0.07
800 !!$ qcoar_no3 = 0.03
801 !!$ qfine_cl = 0.324
802 !!$ qcoar_cl = 0.139
803 !!$ qfine_msa = 0.0
804 !!$ qcoar_msa = 0.0
805 !!$ qfine_co3 = 0.0
806 !!$ qcoar_co3 = 0.0
807 !!$ qfine_nh4 = 0.98
808 !!$ qcoar_nh4 = 0.10
809 !!$ qfine_na = 0.21
810 !!$ qcoar_na = 0.09
811 !!$ qfine_ca = 0.0
812 !!$ qcoar_ca = 0.0
813 !!$ qfine_oin = 0.15
814 !!$ qcoar_oin = 0.35
815 !!$ qfine_oc = 1.00
816 !!$ qcoar_oc = 1.50
817 !!$ qfine_bc = 0.175
818 !!$ qcoar_bc = 0.075
819 !!$ qfine_hysw = 0.0
820 !!$ qcoar_hysw = 0.0
821 !!$ qfine_watr = 0.0
822 !!$ qcoar_watr = 0.0
825 ! qfine_so4 ... are ug/m3 so 1.0e-6 factor gives g/m3
826 ! dens_so4 ... are g/cm3; qfine_vol ... are cm3/m3
827 qfine_vol = 1.0e-6 * ( &
828 (qfine_so4/dens_so4_aer) + (qfine_no3/dens_no3_aer) + &
829 (qfine_cl /dens_cl_aer ) + (qfine_msa/dens_msa_aer) + &
830 (qfine_co3/dens_co3_aer) + (qfine_nh4/dens_nh4_aer) + &
831 (qfine_na /dens_na_aer ) + (qfine_ca /dens_ca_aer ) + &
832 (qfine_oin/dens_oin_aer) + (qfine_oc /dens_oc_aer ) + &
833 (qfine_bc /dens_bc_aer ) )
834 qcoar_vol = 1.0e-6 * ( &
835 (qcoar_so4/dens_so4_aer) + (qcoar_no3/dens_no3_aer) + &
836 (qcoar_cl /dens_cl_aer ) + (qcoar_msa/dens_msa_aer) + &
837 (qcoar_co3/dens_co3_aer) + (qcoar_nh4/dens_nh4_aer) + &
838 (qcoar_na /dens_na_aer ) + (qcoar_ca /dens_ca_aer ) + &
839 (qcoar_oin/dens_oin_aer) + (qcoar_oc /dens_oc_aer ) + &
840 (qcoar_bc /dens_bc_aer ) )
842 do 2900 m = 1, nsize_aer(itype)
843 do 2800 l = 1, 15
845 if (l .eq. 1) then
846 l1 = lptr_so4_aer(m,itype,iphase)
847 qfine = qfine_so4
848 qcoar = qcoar_so4
849 else if (l .eq. 2) then
850 l1 = lptr_no3_aer(m,itype,iphase)
851 qfine = qfine_no3
852 qcoar = qcoar_no3
853 else if (l .eq. 3) then
854 l1 = lptr_cl_aer(m,itype,iphase)
855 qfine = qfine_cl
856 qcoar = qcoar_cl
857 else if (l .eq. 4) then
858 l1 = lptr_msa_aer(m,itype,iphase)
859 qfine = qfine_msa
860 qcoar = qcoar_msa
861 else if (l .eq. 5) then
862 l1 = lptr_co3_aer(m,itype,iphase)
863 qfine = qfine_co3
864 qcoar = qcoar_co3
865 else if (l .eq. 6) then
866 l1 = lptr_nh4_aer(m,itype,iphase)
867 qfine = qfine_nh4
868 qcoar = qcoar_nh4
869 else if (l .eq. 7) then
870 l1 = lptr_na_aer(m,itype,iphase)
871 qfine = qfine_na
872 qcoar = qcoar_na
873 else if (l .eq. 8) then
874 l1 = lptr_ca_aer(m,itype,iphase)
875 qfine = qfine_ca
876 qcoar = qcoar_ca
877 else if (l .eq. 9) then
878 l1 = lptr_oin_aer(m,itype,iphase)
879 qfine = qfine_oin
880 qcoar = qcoar_oin
881 else if (l .eq. 10) then
882 l1 = lptr_oc_aer(m,itype,iphase)
883 qfine = qfine_oc
884 qcoar = qcoar_oc
885 else if (l .eq. 11) then
886 l1 = lptr_bc_aer(m,itype,iphase)
887 qfine = qfine_bc
888 qcoar = qcoar_bc
889 else if (l .eq. 12) then
890 l1 = hyswptr_aer(m,itype)
891 qfine = qfine_hysw
892 qcoar = qcoar_hysw
893 else if (l .eq. 13) then
894 l1 = waterptr_aer(m,itype)
895 qfine = qfine_watr
896 qcoar = qcoar_watr
897 else if (l .eq. 14) then
898 l1 = numptr_aer(m,itype,iphase)
899 dumvol1p = sqrt(volumlo_sect(m,itype)*volumhi_sect(m,itype))
900 qfine = qfine_vol/dumvol1p
901 qcoar = qcoar_vol/dumvol1p
902 vtot_ofmode(m) = &
903 qfine_vol*fr_fine(m) + qcoar_vol*fr_coar(m)
904 else
905 goto 2800
906 end if
907 l3 = l1
909 if ((l1 .gt. 0) .and. (l1 .le. lmaxd) .and. &
910 (l3 .ge. param_first_scalar)) then
911 qval = qfine*fr_fine(m) + qcoar*fr_coar(m)
912 rclm(1,l1) = qval
914 do jt = jts, min(jte,jde-1)
915 do kt = kts, kte-1
916 do it = its, min(ite,ide-1)
918 !wig 28-May-2004, begin block...
919 ! Determine height multiplier...
920 ! This should mimic the calculation in sorgam_set_aer_bc_pnnl,
921 ! sorgam_init_aer_ic_pnnl, bdy_chem_value_mosaic
922 !!$! Height(m) Multiplier
923 !!$! --------- ----------
924 !!$! <=2000 1.0
925 !!$! 2000<z<3000 linear transition zone to 0.5
926 !!$! 3000<z<5000 linear transition zone to 0.25
927 !!$! >=5000 0.25
928 !!$!
929 !!$! which translates to:
930 !!$! 2000<z<3000 mult = 1.0 + (z-2000.)*(0.5-1.0)/(3000.-2000.)
931 !!$! 3000<z<5000 mult = 0.5 + (z-3000.)*(0.25-0.5)/(5000.-3000.)
932 !!$! or in reduced form:
933 !!$ if( z(it,kt,jt) <= 2000. ) then
934 !!$ mult = 1.0
935 !!$ elseif( z(it,kt,jt) > 2000. &
936 !!$ .and. z(it,kt,jt) <= 3000. ) then
937 !!$ mult = 1.0 - 0.0005*(z(it,kt,jt)-2000.)
938 !!$ elseif( z(it,kt,jt) > 3000. &
939 !!$ .and. z(it,kt,jt) <= 5000. ) then
940 !!$ mult = 0.5 - 1.25e-4*(z(it,kt,jt)-3000.)
941 !!$ else
942 !!$ mult = 0.25
943 !!$ end if
944 ! Updated 1-Apr-2005:
945 #if (CASENAME == 1)
946 ! further updated 7-20-05 egc: species specific profiles based on MIRAG2 output
947 mult = 1.0
948 if ( (l3 == 1) .or. (l3 == 2) .or. (l3 == 6) ) then
949 ! so4, no3 and nh4 aerosol profiles
950 if ( z(it,kt,jt) <= 1000. ) then
951 mult = 1.0
952 elseif( z(it,kt,jt) > 1000. &
953 .and. z(it,kt,jt) <= 4000. ) then
954 mult = 1.0 - 2.367e-4*(z(it,kt,jt)-1000.)
955 elseif( z(it,kt,jt) > 4000. &
956 .and. z(it,kt,jt) <= 5500. ) then
957 mult = 0.29 - 4.0e-5*(z(it,kt,jt)-4000.)
958 else
959 mult = 0.23
960 end if
961 else if ( ( l3 == 3) .or. (l3 ==7) ) then
962 ! na and cl aerosol profiles
963 if ( z(it,kt,jt) <= 100. ) then
964 mult = 1.0
965 elseif( z(it,kt,jt) > 100. &
966 .and. z(it,kt,jt) <= 265. ) then
967 mult = 1.0 - 2.9e-3*(z(it,kt,jt)-100.)
968 elseif( z(it,kt,jt) > 265. &
969 .and. z(it,kt,jt) <= 2000. ) then
970 mult = 0.52 - 2.94e-4*(z(it,kt,jt)-265.)
971 elseif( z(it,kt,jt) > 2000. &
972 .and. z(it,kt,jt) <= 7000. ) then
973 mult = 0.01
974 else
975 mult = 1.e-10
976 end if
977 else if ( l3 == 10) then
978 ! oc aerosol profile
979 if ( z(it,kt,jt) <= 600. ) then
980 mult = 1.0
981 elseif( z(it,kt,jt) > 600. &
982 .and. z(it,kt,jt) <= 3389. ) then
983 mult = 1.0 - 2.04e-4*(z(it,kt,jt)-600.)
984 elseif( z(it,kt,jt) > 3389. &
985 .and. z(it,kt,jt) <= 8840. ) then
986 mult = 0.43 - 7.045e-5*(z(it,kt,jt)-3389.)
987 else
988 mult = 0.046
989 end if
990 else if ( l3 == 11) then
991 ! bc aerosol profile
992 if ( z(it,kt,jt) <= 100. ) then
993 mult = 1.0
994 elseif( z(it,kt,jt) > 100. &
995 .and. z(it,kt,jt) <= 3400. ) then
996 mult = 1.0 - 2.51e-4*(z(it,kt,jt)-100.)
997 elseif( z(it,kt,jt) > 3400. &
998 .and. z(it,kt,jt) <= 8400. ) then
999 mult = 0.172 -2.64e-5*(z(it,kt,jt)-3400.)
1000 else
1001 mult = 0.04
1002 end if
1003 else if ( l3 == 9) then
1004 ! oin aerosol profile
1005 if ( z(it,kt,jt) <= 1580. ) then
1006 mult = 1.0 + 1.77e-4 *z(it,kt,jt)
1007 elseif( z(it,kt,jt) > 1580. &
1008 .and. z(it,kt,jt) <= 5280. ) then
1009 mult = 1.28 - 2.65e-4*(z(it,kt,jt)-1580.)
1010 else
1011 mult = 0.30
1012 end if
1013 else
1014 ! generic profile for other other species (which should have groundlevel values=0)
1015 #endif
1016 ! Height(m) Multiplier
1017 ! --------- ----------
1018 ! <=2000 1.0
1019 ! 2000<z<3000 linear transition zone to 0.25
1020 ! 3000<z<5000 linear transision zone to 0.125
1021 ! >=5000 0.125
1022 !
1023 ! which translates to:
1024 ! 2000<z<3000 mult = 1.00 + (z-2000.)*(0.25-1.0)/(3000.-2000.)
1025 ! 3000<z<5000 mult = 0.25 + (z-3000.)*(0.125-0.25)/(5000.-3000.)
1026 ! or in reduced form:
1027 if( z(it,kt,jt) <= 2000. ) then
1028 mult = 1.0
1029 elseif( z(it,kt,jt) > 2000. &
1030 .and. z(it,kt,jt) <= 3000. ) then
1031 mult = 1.0 - 0.00075*(z(it,kt,jt)-2000.)
1032 elseif( z(it,kt,jt) > 3000. &
1033 .and. z(it,kt,jt) <= 5000. ) then
1034 mult = 0.25 - 4.166666667e-5*(z(it,kt,jt)-3000.)
1035 else
1036 mult = 0.125
1037 end if
1038 #if (CASENAME == 1)
1039 end if !close l3 comparison block
1040 #endif
1041 #if (CASENAME == 2)
1042 ! if( z(it,kt,jt) <= 3000. ) then
1043 ! mult = 1.0
1044 ! elseif( z(it,kt,jt) > 3000. &
1045 ! .and. z(it,kt,jt) <= 5000. ) then
1046 ! mult = 1.0 - 0.0004375*(z(it,kt,jt)-3000.)
1047 ! else
1048 ! mult = 0.125
1049 ! end if
1050 !
1051 ! if( z(it,kt,jt) <= 5000. ) then
1052 ! mult = 1.0
1053 ! elseif( z(it,kt,jt) > 5000. &
1054 ! .and. z(it,kt,jt) <= 7000. ) then
1055 ! mult = 1.0 - 0.0004375*(z(it,kt,jt)-5000.)
1056 ! else
1057 ! mult = 0.125
1058 ! end if
1059 mult=1.0
1060 if ( (l3 == lptr_so4_aer(m,itype,iphase)) .or. &
1061 (l3 == lptr_no3_aer(m,itype,iphase)) .or. &
1062 (l3 == lptr_nh4_aer(m,itype,iphase)) ) then
1063 ! so4, no3 and nh4 aerosol profiles
1064 if( z(it,kt,jt) <= 500. ) then
1065 mult = 1.000
1066 elseif( z(it,kt,jt) > 500. &
1067 .and. z(it,kt,jt) <= 1000. ) then
1068 mult = 1.000 - 0.0010740*(z(it,kt,jt)-500.)
1069 elseif( z(it,kt,jt) > 1000. &
1070 .and. z(it,kt,jt) <= 5000. ) then
1071 mult = 0.463 - 0.0001110*(z(it,kt,jt)-1000.)
1072 elseif( z(it,kt,jt) > 5000. &
1073 .and. z(it,kt,jt) <= 20000. ) then
1074 mult = 0.019
1075 else
1076 mult = 0.019
1077 end if
1078 else if ( (l3 == lptr_na_aer(m,itype,iphase)) .or. &
1079 (l3 == lptr_cl_aer(m,itype,iphase)) ) then
1080 ! na and cl aerosol profiles
1081 if( z(it,kt,jt) <= 500. ) then
1082 mult = 1.000
1083 elseif( z(it,kt,jt) > 500. &
1084 .and. z(it,kt,jt) <= 1000. ) then
1085 mult = 1.000 - 0.0014260*(z(it,kt,jt)-500.)
1086 elseif( z(it,kt,jt) > 1000. &
1087 .and. z(it,kt,jt) <= 5000. ) then
1088 mult = 0.287 - 0.0000643*(z(it,kt,jt)-1000.)
1089 elseif( z(it,kt,jt) > 5000. &
1090 .and. z(it,kt,jt) <= 10000. ) then
1091 mult = 0.030
1092 elseif( z(it,kt,jt) > 10000. &
1093 .and. z(it,kt,jt) <= 20000. ) then
1094 mult = 0.030 - 0.0000028*(z(it,kt,jt)-10000.)
1095 else
1096 mult = 0.002
1097 end if
1098 else if ( (l3 == lptr_oc_aer(m,itype,iphase)) .or. &
1099 (l3 == lptr_bc_aer(m,itype,iphase)) ) then
1100 ! oc and bc aerosol profile
1101 if( z(it,kt,jt) <= 500. ) then
1102 mult = 1.000
1103 elseif( z(it,kt,jt) > 500. &
1104 .and. z(it,kt,jt) <= 1000. ) then
1105 mult = 1.000 - 0.0002500*(z(it,kt,jt)-500.)
1106 elseif( z(it,kt,jt) > 1000. &
1107 .and. z(it,kt,jt) <= 5000. ) then
1108 mult = 0.875 - 0.0001843*(z(it,kt,jt)-1000.)
1109 elseif( z(it,kt,jt) > 5000. &
1110 .and. z(it,kt,jt) <= 20000. ) then
1111 mult = 0.138
1112 else
1113 mult = 0.138
1114 end if
1115 else if ( l3 == lptr_oin_aer(m,itype,iphase) ) then
1116 ! oin aerosol profile
1117 if( z(it,kt,jt) <= 500. ) then
1118 mult = 1.000
1119 elseif( z(it,kt,jt) > 500. &
1120 .and. z(it,kt,jt) <= 1000. ) then
1121 mult = 1.000 - 0.0007340*(z(it,kt,jt)-500.)
1122 elseif( z(it,kt,jt) > 1000. &
1123 .and. z(it,kt,jt) <= 5000. ) then
1124 mult = 0.633 - 0.0001103*(z(it,kt,jt)-1000.)
1125 elseif( z(it,kt,jt) > 5000. &
1126 .and. z(it,kt,jt) <= 10000. ) then
1127 mult = 0.192 - 0.0000080*(z(it,kt,jt)-5000.)
1128 elseif( z(it,kt,jt) > 10000. &
1129 .and. z(it,kt,jt) <= 20000. ) then
1130 mult = 0.152 - 0.0000137*(z(it,kt,jt)-10000.)
1131 else
1132 mult = 0.015
1133 end if
1134 else
1135 if( z(it,kt,jt) <= 5000. ) then
1136 mult = 1.0
1137 elseif( z(it,kt,jt) > 5000. &
1138 .and. z(it,kt,jt) <= 7000. ) then
1139 mult = 1.0 - 0.0004375*(z(it,kt,jt)-5000.)
1140 else
1141 mult = 0.125
1142 end if
1143 end if
1144 #endif
1145 #if (CASENAME == 3)
1146 mult = 1.000
1147 if( z(it,kt,jt) <= 500. ) then
1148 mult = 1.000
1149 elseif( z(it,kt,jt) > 500. &
1150 .and. z(it,kt,jt) <= 1000. ) then
1151 mult = 1.000 - 0.0010740*(z(it,kt,jt)-500.)
1152 elseif( z(it,kt,jt) > 1000. &
1153 .and. z(it,kt,jt) <= 5000. ) then
1154 mult = 0.463 - 0.0001110*(z(it,kt,jt)-1000.)
1155 else
1156 mult = 0.019
1157 end if
1158 #endif
1160 chem(it,kt,jt,l3) = mult*rclm(1,l1)
1161 !wig 28-May-2004, ...end block
1162 ! chem(it,kt,jt,l3) = rclm(1,l1)
1163 end do
1164 end do
1165 end do
1166 end if
1169 2800 continue
1170 2900 continue
1172 !
1173 ! do diagnostic output
1174 !
1175 dumarr(:) = 0.0
1176 msg = ' '
1177 call peg_message( lunout, msg )
1178 msg = '*** subr mosaic_init_wrf_mixrats_opt2 results'
1179 call peg_message( lunout, msg )
1180 msg = ' mass in ug/m3 number in #/m3 volume in cm3/m3'
1181 call peg_message( lunout, msg )
1182 msg = ' '
1183 call peg_message( lunout, msg )
1184 msg = ' mode l l1 species conc'
1185 call peg_message( lunout, msg )
1187 do 3190 mdum = 1, nsize_aer(itype)+1
1188 m = min( mdum, nsize_aer(itype) )
1189 msg = ' '
1190 call peg_message( lunout, msg )
1191 do 3150 l = 1, ncomp_plustracer_aer(itype)+4
1193 if (l .le. ncomp_plustracer_aer(itype)) then
1194 l1 = massptr_aer(l,m,itype,iphase)
1195 dumname = name_aer(l,itype)
1196 dum = rclm(1,l1)
1197 else if (l .eq. ncomp_plustracer_aer(itype)+1) then
1198 l1 = hyswptr_aer(m,itype)
1199 dumname = 'hystwatr'
1200 dum = rclm(1,l1)
1201 else if (l .eq. ncomp_plustracer_aer(itype)+2) then
1202 l1 = waterptr_aer(m,itype)
1203 dumname = 'water'
1204 dum = rclm(1,l1)
1205 else if (l .eq. ncomp_plustracer_aer(itype)+3) then
1206 l1 = numptr_aer(m,itype,iphase)
1207 dumname = 'number'
1208 dum = rclm(1,l1)
1209 else if (l .eq. ncomp_plustracer_aer(itype)+4) then
1210 l1 = 0
1211 dumname = 'volume'
1212 dum = vtot_ofmode(m)
1213 else
1214 dumname = '=BADBAD='
1215 l1 = -1
1216 dum = -1.0
1217 end if
1219 if (mdum .le. nsize_aer(itype)) then
1220 dumarr(l) = dumarr(l) + dum
1221 write(msg,9620) m, l, l1, dumname, dum
1222 else
1223 write(msg,9625) l, dumname, dumarr(l)
1224 end if
1225 call peg_message( lunout, msg )
1227 3150 continue
1228 3190 continue
1230 9620 format( 3i4, 2x, a, 3(1pe12.3) )
1231 9625 format( ' sum', i4, ' -', 2x, a, 3(1pe12.3) )
1234 ! all done
1235 return
1236 end subroutine mosaic_init_wrf_mixrats_opt2
1239 !-----------------------------------------------------------------------
1240 real function erfc_num_recipes( x )
1241 !
1242 ! from press et al, numerical recipes, 1990, page 164
1243 !
1244 implicit none
1245 real x
1246 double precision erfc_dbl, dum, t, zz
1248 zz = abs(x)
1249 t = 1.0/(1.0 + 0.5*zz)
1251 ! erfc_num_recipes =
1252 ! & t*exp( -zz*zz - 1.26551223 + t*(1.00002368 + t*(0.37409196 +
1253 ! & t*(0.09678418 + t*(-0.18628806 + t*(0.27886807 +
1254 ! & t*(-1.13520398 +
1255 ! & t*(1.48851587 + t*(-0.82215223 + t*0.17087277 )))))))))
1257 dum = ( -zz*zz - 1.26551223 + t*(1.00002368 + t*(0.37409196 + &
1258 t*(0.09678418 + t*(-0.18628806 + t*(0.27886807 + &
1259 t*(-1.13520398 + &
1260 t*(1.48851587 + t*(-0.82215223 + t*0.17087277 )))))))))
1262 erfc_dbl = t * exp(dum)
1263 if (x .lt. 0.0) erfc_dbl = 2.0d0 - erfc_dbl
1265 erfc_num_recipes = erfc_dbl
1267 return
1268 end function erfc_num_recipes
1271 !-----------------------------------------------------------------------
1272 end module module_mosaic_initmixrats
1277 !-----------------------------------------------------------------------
1278 subroutine bdy_chem_value_mosaic ( chem_bv, alt, zz, nch, config_flags )
1279 !
1280 ! provides boundary values for the mosaic aerosol species
1281 !
1282 ! it is outside of the module declaration because of potential
1283 ! module1 --> module2 --> module1 use conflicts
1284 !
1285 ! rce 11-sep-2004 - eliminated use of the _wrfch pointers (lptr_xxx_a_wrfch,
1286 ! lwaterptr_wrfch, numptr_wrfch); use only the _amode pointers now
1287 !
1288 use module_configure, only: grid_config_rec_type
1289 use module_state_description, only: param_first_scalar, &
1290 aer_bc_pnnl
1291 use module_data_mosaic_asect
1292 use module_data_mosaic_other
1293 implicit none
1295 ! arguments
1296 REAL, intent(OUT) :: chem_bv ! boundary value for chem(-,-,-,nch)
1297 REAL, intent(IN) :: alt ! inverse density
1298 REAL, intent(IN) :: zz ! height
1299 INTEGER, intent(IN) :: nch ! index number of chemical species
1300 TYPE (grid_config_rec_type), intent(in) :: config_flags
1302 ! local variables
1303 integer :: iphase, itype, m
1304 logical :: foundit
1306 real, parameter :: chem_bv_def = 1.0e-20
1307 !wig 10-May-2004, added mult
1308 real :: dumvol1p, mult
1309 integer :: iflg
1310 real :: qcoar, qfine, qval
1312 real :: fr_coar(8), fr_fine(8)
1314 !wig 1-Apr-2005, Updated fractional breakdown between bins. (See also
1315 ! mosaic_init_wrf_mixrats_opt2, above, and mosaic_addemiss
1316 ! in module_mosaic_addemiss.F). Note that these values no
1317 ! longer match those in mosaic_addemiss.
1318 real, save :: fr8b_coar(8) = &
1319 (/ 0.0, 0.0, 0.0, 0.019, 0.0212, 0.1101, 0.2751, 0.7018 /) ! jdf Testbed
1320 !jdf (/ 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.300, 0.700 /) ! 10-May-2005
1321 ! (/ 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.933, 0.067 /) ! 01-Apr-2005
1322 ! (/ 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.16, 0.84 /)
1323 real, save :: fr8b_fine(8) = &
1324 (/ 0.0245, 0.1472, 0.3501, 0.3317, 0.1251, 0.0186, 0.0011, 0.0/) !jdf Testbed
1325 !jdf (/ 0.006, 0.024, 0.231, 0.341, 0.140, 0.258, 0., 0./) !5-Apr-2005 values
1326 ! (/ 0.0275, 0.0426, 0.2303, 0.3885, 0.1100, 0.2011, 0., 0./) !15-Nov-2004 values
1327 ! (/ 0.01, 0.05, 0.145, 0.60, 0.145, 0.05, 0.00, 0.00 /)
1328 ! (/ 0.04, 0.10, 0.35, 0.29, 0.15, 0.07, 0.00, 0.00 /)
1330 real :: qfine_so4, qfine_no3, qfine_cl, qfine_msa, &
1331 qfine_co3, qfine_nh4, qfine_na, qfine_ca, qfine_oin, &
1332 qfine_oc, qfine_bc, qfine_hysw, qfine_watr, qfine_vol
1333 real :: qcoar_so4, qcoar_no3, qcoar_cl, qcoar_msa, &
1334 qcoar_co3, qcoar_nh4, qcoar_na, qcoar_ca, qcoar_oin, &
1335 qcoar_oc, qcoar_bc, qcoar_hysw, qcoar_watr, qcoar_vol
1337 character*80 msg
1341 ! when aer_bc_opt /= aer_bc_pnnl,
1342 ! just chem_bv=near-zero for all species
1343 chem_bv = chem_bv_def
1344 if (config_flags%aer_bc_opt /= aer_bc_pnnl) return
1345 if (nch < param_first_scalar) return
1348 ! *** currently only works for ntype_aer = 1
1349 itype = 1
1350 iphase = ai_phase
1351 m = 1 !This is here for size, but is overridden by loop below.
1354 !
1355 ! following is for aer_bc_opt == aer_bc_pnnl
1356 ! when boundary values are set for Texas August 2000 simulations
1357 !
1358 ! set the partitioning fractions for either 8 or 4 bins
1359 if (nsize_aer(itype) == 8) then
1360 fr_coar(:) = fr8b_coar(:)
1361 fr_fine(:) = fr8b_fine(:)
1362 else if (nsize_aer(itype) == 4) then
1363 do m = 1, nsize_aer(itype)
1364 fr_coar(m) = fr8b_coar(2*m) + fr8b_coar(2*m-1)
1365 fr_fine(m) = fr8b_fine(2*m) + fr8b_fine(2*m-1)
1366 end do
1367 else
1368 write(msg,'(a,i5)') &
1369 'subr bdy_chem_value_mosaic' // &
1370 ' - nsize_aer(itype) must be 4 or 8 but = ', nsize_aer(itype)
1371 call wrf_error_fatal( msg )
1372 end if
1374 ! Determine height multiplier...
1375 ! This should mimic the calculation in sorgam_set_aer_bc_pnnl,
1376 ! sorgam_init_aer_ic_pnnl, and mosaic_init_wrf_mixrats_opt2
1377 !!$! Height(m) Multiplier
1378 !!$! --------- ----------
1379 !!$! <=2000 1.0
1380 !!$! 2000<z<3000 linear transition zone to 0.5
1381 !!$! 3000<z<5000 linear transision zone to 0.25
1382 !!$! >=5000 0.25
1383 !!$!
1384 !!$! which translates to:
1385 !!$! 2000<z<3000 mult = 1.0 + (z-2000.)*(0.5-1.0)/(3000.-2000.)
1386 !!$! 3000<z<5000 mult = 0.5 + (z-3000.)*(0.25-0.5)/(5000.-3000.)
1387 !!$! or in reduced form:
1388 !!$ if( zz <= 2000. ) then
1389 !!$ mult = 1.0
1390 !!$ elseif( zz > 2000. &
1391 !!$ .and. zz <= 3000. ) then
1392 !!$ mult = 1.0 - 0.0005*(zz-2000.)
1393 !!$ elseif( zz > 3000. &
1394 !!$ .and. zz <= 5000. ) then
1395 !!$ mult = 0.5 - 1.25e-4*(zz-3000.)
1396 !!$ else
1397 !!$ mult = 0.25
1398 !!$ end if
1399 #if (CASENAME == 1)
1400 mult = 1.0
1401 SIZE_LOOP : do m=1,nsize_aer(itype)
1402 if( ( nch .eq. lptr_so4_aer(m,itype,iphase) ) .or. &
1403 (nch .eq. lptr_no3_aer(m,itype,iphase) ) .or. &
1404 (nch .eq. lptr_nh4_aer(m,itype,iphase) ) )then
1405 ! so4, no3 and nh4 aerosol profiles
1406 if ( zz <= 1000. ) then
1407 mult = 1.0
1408 elseif( zz > 1000. &
1409 .and. zz <= 4000. ) then
1410 mult = 1.0 - 2.367e-4*(zz-1000.)
1411 elseif( zz > 4000. &
1412 .and. zz <= 5500. ) then
1413 mult = 0.29 - 4.0e-5*(zz-4000.)
1414 else
1415 mult = 0.23
1416 end if
1417 exit SIZE_LOOP
1418 else if ( (nch .eq. lptr_na_aer(m,itype,iphase) ) .or. &
1419 (nch .eq. lptr_cl_aer(m,itype,iphase) ) ) then
1420 ! na and cl aerosol profiles
1421 if ( zz <= 100. ) then
1422 mult = 1.0
1423 elseif( zz > 100. &
1424 .and. zz <= 265. ) then
1425 mult = 1.0 - 2.9e-3*(zz-100.)
1426 elseif( zz > 265. &
1427 .and. zz <= 2000. ) then
1428 mult = 0.52 - 2.94e-4*(zz-265.)
1429 elseif( zz > 2000. &
1430 .and. zz <= 7000. ) then
1431 mult = 0.01
1432 else
1433 mult = 1.e-10
1434 end if
1435 exit SIZE_LOOP
1436 else if (nch .eq. lptr_oc_aer(m,itype,iphase) ) then
1437 ! oc aerosol profile
1438 if ( zz <= 600. ) then
1439 mult = 1.0
1440 elseif( zz > 600. &
1441 .and. zz <= 3389. ) then
1442 mult = 1.0 - 2.04e-4*(zz-600.)
1443 elseif( zz > 3389. &
1444 .and. zz <= 8840. ) then
1445 mult = 0.43 - 7.045e-5*(zz-3389.)
1446 else
1447 mult = 0.046
1448 end if
1449 exit SIZE_LOOP
1450 else if (nch .eq. lptr_bc_aer(m,itype,iphase) ) then
1451 ! bc aerosol profile
1452 if ( zz <= 100. ) then
1453 mult = 1.0
1454 elseif( zz > 100. &
1455 .and. zz <= 3400. ) then
1456 mult = 1.0 - 2.51e-4*(zz-100.)
1457 elseif( zz > 3400. &
1458 .and. zz <= 8400. ) then
1459 mult = 0.172 -2.64e-5*(zz-3400.)
1460 else
1461 mult = 0.04
1462 end if
1463 exit SIZE_LOOP
1464 else if (nch .eq. lptr_oin_aer(m,itype,iphase)) then
1465 ! oin aerosol profile
1466 if ( zz <= 1580. ) then
1467 mult = 1.0 + 1.77e-4 *zz
1468 elseif( zz > 1580. &
1469 .and. zz <= 5280. ) then
1470 mult = 1.28 - 2.65e-4*(zz-1580.)
1471 else
1472 mult = 0.30
1473 end if
1474 exit SIZE_LOOP
1475 else
1476 ! generic profile
1477 #endif
1478 ! Updated aerosol profile multiplier 1-Apr-2005:
1479 ! Height(m) Multiplier
1480 ! --------- ----------
1481 ! <=2000 1.0
1482 ! 2000<z<3000 linear transition zone to 0.25
1483 ! 3000<z<5000 linear transision zone to 0.125
1484 ! >=5000 0.125
1485 !
1486 ! which translates to:
1487 ! 2000<z<3000 mult = 1.00 + (z-2000.)*(0.25-1.0)/(3000.-2000.)
1488 ! 3000<z<5000 mult = 0.25 + (z-3000.)*(0.125-0.25)/(5000.-3000.)
1489 ! or in reduced form:
1490 if( zz <= 2000. ) then
1491 mult = 1.0
1492 elseif( zz > 2000. &
1493 .and. zz <= 3000. ) then
1494 mult = 1.0 - 0.00075*(zz-2000.)
1495 elseif( zz > 3000. &
1496 .and. zz <= 5000. ) then
1497 mult = 0.25 - 4.166666667e-5*(zz-3000.)
1498 else
1499 mult = 0.125
1500 end if
1501 #if (CASENAME == 1)
1502 end if ! end nc block comparison
1503 end do SIZE_LOOP
1504 #endif
1505 #if (CASENAME == 2)
1506 ! if( zz <= 3000. ) then
1507 ! mult = 1.0
1508 ! elseif( zz > 3000. &
1509 ! .and. zz <= 5000. ) then
1510 ! mult = 1.0 - 0.0004375*(zz-3000.)
1511 ! else
1512 ! mult = 0.125
1513 ! end if
1514 ! if( zz <= 5000. ) then
1515 ! mult = 1.0
1516 ! elseif( zz > 5000. &
1517 ! .and. zz <= 7000. ) then
1518 ! mult = 1.0 - 0.0004375*(zz-5000.)
1519 ! else
1520 ! mult = 0.125
1521 ! end if
1522 mult=1.0
1523 iflg=0
1524 do 2901 m = 1, nsize_aer(itype)
1525 if(iflg.eq.1) goto 2902
1526 if( ( nch .eq. lptr_so4_aer(m,itype,iphase) ) .or. &
1527 (nch .eq. lptr_no3_aer(m,itype,iphase) ) .or. &
1528 (nch .eq. lptr_nh4_aer(m,itype,iphase) ) )then
1529 ! so4, no3 and nh4 aerosol profiles
1530 if( zz <= 500. ) then
1531 mult = 1.000
1532 elseif( zz > 500. &
1533 .and. zz <= 1000. ) then
1534 mult = 1.000 - 0.0010740*(zz-500.)
1535 elseif( zz > 1000. &
1536 .and. zz <= 5000. ) then
1537 mult = 0.463 - 0.0001110*(zz-1000.)
1538 elseif( zz > 5000. &
1539 .and. zz <= 20000. ) then
1540 mult = 0.019
1541 else
1542 mult = 0.019
1543 end if
1544 iflg=1
1545 else if ( (nch .eq. lptr_na_aer(m,itype,iphase) ) .or. &
1546 (nch .eq. lptr_cl_aer(m,itype,iphase) ) ) then
1547 ! na and cl aerosol profiles
1548 if( zz <= 500. ) then
1549 mult = 1.000
1550 elseif( zz > 500. &
1551 .and. zz <= 1000. ) then
1552 mult = 1.000 - 0.0014260*(zz-500.)
1553 elseif( zz > 1000. &
1554 .and. zz <= 5000. ) then
1555 mult = 0.287 - 0.0000643*(zz-1000.)
1556 elseif( zz > 5000. &
1557 .and. zz <= 10000. ) then
1558 mult = 0.030
1559 elseif( zz > 10000. &
1560 .and. zz <= 20000. ) then
1561 mult = 0.030 - 0.0000028*(zz-10000.)
1562 else
1563 mult = 0.002
1564 end if
1565 iflg=1
1566 else if ( (nch .eq. lptr_oc_aer(m,itype,iphase) ) .or. &
1567 (nch .eq. lptr_bc_aer(m,itype,iphase) ) ) then
1568 ! oc and bc aerosol profile
1569 if( zz <= 500. ) then
1570 mult = 1.000
1571 elseif( zz > 500. &
1572 .and. zz <= 1000. ) then
1573 mult = 1.000 - 0.0002500*(zz-500.)
1574 elseif( zz > 1000. &
1575 .and. zz <= 5000. ) then
1576 mult = 0.875 - 0.0001843*(zz-1000.)
1577 elseif( zz > 5000. &
1578 .and. zz <= 20000. ) then
1579 mult = 0.138
1580 else
1581 mult = 0.138
1582 end if
1583 iflg=1
1584 else if (nch .eq. lptr_oin_aer(m,itype,iphase)) then
1585 ! oin aerosol profile
1586 if( zz <= 500. ) then
1587 mult = 1.000
1588 elseif( zz > 500. &
1589 .and. zz <= 1000. ) then
1590 mult = 1.000 - 0.0007340*(zz-500.)
1591 elseif( zz > 1000. &
1592 .and. zz <= 5000. ) then
1593 mult = 0.633 - 0.0001103*(zz-1000.)
1594 elseif( zz > 5000. &
1595 .and. zz <= 10000. ) then
1596 mult = 0.192 - 0.0000080*(zz-5000.)
1597 elseif( zz > 10000. &
1598 .and. zz <= 20000. ) then
1599 mult = 0.152 - 0.0000137*(zz-10000.)
1600 else
1601 mult = 0.015
1602 end if
1603 iflg=1
1604 else
1605 if( zz <= 5000. ) then
1606 mult = 1.0
1607 elseif( zz > 5000. &
1608 .and. zz <= 7000. ) then
1609 mult = 1.0 - 0.0004375*(zz-5000.)
1610 else
1611 mult = 0.125
1612 end if
1613 end if
1614 2901 continue
1615 2902 continue
1616 #endif
1617 #if (CASENAME == 3)
1618 if( zz <= 500. ) then
1619 mult = 1.000
1620 elseif( zz > 500. &
1621 .and. zz <= 1000. ) then
1622 mult = 1.000 - 0.0010740*(zz-500.)
1623 elseif( zz > 1000. &
1624 .and. zz <= 5000. ) then
1625 mult = 0.463 - 0.0001110*(zz-1000.)
1626 else
1627 mult = 0.019
1628 end if
1629 #endif
1630 ! Set values for fine and coarse mass concentrations, and
1631 ! compute fine/coarse volume concentrations. These are also set in
1632 ! mosaic_init_wrf_mixrats_opt2.
1633 ! (these are latest values from 1-Apr-2005 discussions)
1634 !wig 10-May-2004, added height multiplier (mult*) to each species...
1635 qfine_so4 = mult*2.14
1636 qcoar_so4 = mult*0.242
1637 qfine_no3 = mult*0.11
1638 qcoar_no3 = mult*0.03
1639 ! qfine_cl = mult*0.3
1640 qfine_cl = mult*0.14 ! 10-May-2005
1641 qcoar_cl = mult*0.139
1642 qfine_msa = mult*0.0
1643 qcoar_msa = mult*0.0
1644 qfine_co3 = mult*0.0
1645 qcoar_co3 = mult*0.0
1646 qfine_nh4 = mult*0.83
1647 qcoar_nh4 = mult*0.10
1648 ! qfine_na = mult*0.2
1649 qfine_na = mult*0.1 ! 10-May-2005
1650 qcoar_na = mult*0.09
1651 qfine_ca = mult*0.0
1652 qcoar_ca = mult*0.0
1653 ! qfine_oin = mult*6.2
1654 ! qfine_oin = 3.48 ! 10-May-2005
1655 qfine_oin = mult*3.48 ! not sure why mult was dropped from qfine_oin
1656 qcoar_oin = mult*0.35
1657 qfine_oc = mult*1.00
1658 qcoar_oc = mult*1.50
1659 qfine_bc = mult*0.2
1660 qcoar_bc = mult*0.075
1661 qfine_hysw = mult*0.0
1662 qcoar_hysw = mult*0.0
1663 qfine_watr = mult*0.0
1664 qcoar_watr = mult*0.0
1665 #if (CASENAME == 2)
1666 ! qfine_so4 = mult*0.50
1667 ! qcoar_so4 = mult*0.50/2.0
1668 ! qfine_no3 = mult*0.75
1669 ! qcoar_no3 = mult*0.75/2.0
1670 ! qfine_cl = mult*0.10
1671 ! qcoar_cl = mult*0.10/2.0
1672 ! qfine_nh4 = mult*0.30
1673 ! qcoar_nh4 = mult*0.30/2.0
1674 ! qfine_na = mult*0.20
1675 ! qcoar_na = mult*0.20/2.0
1676 ! qfine_oin = mult*4.00
1677 ! qcoar_oin = mult*4.00/2.0
1678 ! qfine_oc = mult*3.00
1679 ! qcoar_oc = mult*3.00/2.0
1680 ! qfine_bc = mult*0.50
1681 ! qcoar_bc = mult*0.50/2.0
1682 qfine_so4 = mult*0.300
1683 qcoar_so4 = mult*0.300/2.0
1684 qfine_no3 = mult*0.001
1685 qcoar_no3 = mult*0.001/2.0
1686 qfine_cl = mult*1.050
1687 qcoar_cl = mult*1.05/02.0
1688 qfine_nh4 = mult*0.094
1689 qcoar_nh4 = mult*0.094/2.0
1690 qfine_na = mult*0.700
1691 qcoar_na = mult*0.700/2.0
1692 qfine_oin = mult*4.500
1693 qcoar_oin = mult*4.500/2.0
1694 qfine_oc = mult*0.088
1695 qcoar_oc = mult*0.088/2.0
1696 qfine_bc = mult*0.013
1697 qcoar_bc = mult*0.013/2.0
1698 #endif
1699 #if (CASENAME == 3)
1700 qfine_so4 = mult*0.300
1701 qcoar_so4 = mult*0.0
1702 qfine_nh4 = mult*0.094
1703 qcoar_nh4 = mult*0.0
1704 qfine_no3 = mult*0.001
1705 qcoar_no3 = mult*0.0
1706 qfine_cl = mult*0.0
1707 qcoar_cl = mult*1.050
1708 qfine_na = mult*0.000
1709 qcoar_na = mult*0.700
1710 qfine_oin = mult*4.500
1711 qcoar_oin = mult*4.500/2.0
1712 qfine_oc = mult*0.088
1713 qcoar_oc = mult*0.0
1714 qfine_bc = mult*0.013
1715 qcoar_bc = mult*0.0
1716 #endif
1717 !!$! old values from 23-Apr-2004:
1718 !!$ qfine_so4 = mult*2.554
1719 !!$ qcoar_so4 = mult*0.242
1720 !!$ qfine_no3 = mult*0.07
1721 !!$ qcoar_no3 = mult*0.03
1722 !!$ qfine_cl = mult*0.324
1723 !!$ qcoar_cl = mult*0.139
1724 !!$ qfine_msa = mult*0.0
1725 !!$ qcoar_msa = mult*0.0
1726 !!$ qfine_co3 = mult*0.0
1727 !!$ qcoar_co3 = mult*0.0
1728 !!$ qfine_nh4 = mult*0.98
1729 !!$ qcoar_nh4 = mult*0.10
1730 !!$ qfine_na = mult*0.21
1731 !!$ qcoar_na = mult*0.09
1732 !!$ qfine_ca = mult*0.0
1733 !!$ qcoar_ca = mult*0.0
1734 !!$ qfine_oin = mult*0.15
1735 !!$ qcoar_oin = mult*0.35
1736 !!$ qfine_oc = mult*1.00
1737 !!$ qcoar_oc = mult*1.50
1738 !!$ qfine_bc = mult*0.175
1739 !!$ qcoar_bc = mult*0.075
1740 !!$ qfine_hysw = mult*0.0
1741 !!$ qcoar_hysw = mult*0.0
1742 !!$ qfine_watr = mult*0.0
1743 !!$ qcoar_watr = mult*0.0
1746 ! qfine_so4 ... are ug/m3 so 1.0e-6 factor gives g/m3
1747 ! dens_so4 ... are g/cm3; qfine_vol ... are cm3/m3
1748 qfine_vol = 1.0e-6 * ( &
1749 (qfine_so4/dens_so4_aer) + (qfine_no3/dens_no3_aer) + &
1750 (qfine_cl /dens_cl_aer ) + (qfine_msa/dens_msa_aer) + &
1751 (qfine_co3/dens_co3_aer) + (qfine_nh4/dens_nh4_aer) + &
1752 (qfine_na /dens_na_aer ) + (qfine_ca /dens_ca_aer ) + &
1753 (qfine_oin/dens_oin_aer) + (qfine_oc /dens_oc_aer ) + &
1754 (qfine_bc /dens_bc_aer ) )
1755 qcoar_vol = 1.0e-6 * ( &
1756 (qcoar_so4/dens_so4_aer) + (qcoar_no3/dens_no3_aer) + &
1757 (qcoar_cl /dens_cl_aer ) + (qcoar_msa/dens_msa_aer) + &
1758 (qcoar_co3/dens_co3_aer) + (qcoar_nh4/dens_nh4_aer) + &
1759 (qcoar_na /dens_na_aer ) + (qcoar_ca /dens_ca_aer ) + &
1760 (qcoar_oin/dens_oin_aer) + (qcoar_oc /dens_oc_aer ) + &
1761 (qcoar_bc /dens_bc_aer ) )
1763 qfine = -1.0e30
1764 qcoar = -1.0e30
1766 ! identify the species by checking "nch" against the "lptr_xxx_a_amode(m)"
1767 do 2900 m = 1, nsize_aer(itype)
1769 if (nch .eq. lptr_so4_aer(m,itype,iphase)) then
1770 qfine = qfine_so4
1771 qcoar = qcoar_so4
1772 else if (nch .eq. lptr_no3_aer(m,itype,iphase)) then
1773 qfine = qfine_no3
1774 qcoar = qcoar_no3
1775 else if (nch .eq. lptr_cl_aer(m,itype,iphase)) then
1776 qfine = qfine_cl
1777 qcoar = qcoar_cl
1778 else if (nch .eq. lptr_msa_aer(m,itype,iphase)) then
1779 qfine = qfine_msa
1780 qcoar = qcoar_msa
1781 else if (nch .eq. lptr_co3_aer(m,itype,iphase)) then
1782 qfine = qfine_co3
1783 qcoar = qcoar_co3
1784 else if (nch .eq. lptr_nh4_aer(m,itype,iphase)) then
1785 qfine = qfine_nh4
1786 qcoar = qcoar_nh4
1787 else if (nch .eq. lptr_na_aer(m,itype,iphase)) then
1788 qfine = qfine_na
1789 qcoar = qcoar_na
1790 else if (nch .eq. lptr_ca_aer(m,itype,iphase)) then
1791 qfine = qfine_ca
1792 qcoar = qcoar_ca
1793 else if (nch .eq. lptr_oin_aer(m,itype,iphase)) then
1794 qfine = qfine_oin
1795 qcoar = qcoar_oin
1796 else if (nch .eq. lptr_oc_aer(m,itype,iphase)) then
1797 qfine = qfine_oc
1798 qcoar = qcoar_oc
1799 else if (nch .eq. lptr_bc_aer(m,itype,iphase)) then
1800 qfine = qfine_bc
1801 qcoar = qcoar_bc
1802 else if (nch .eq. hyswptr_aer(m,itype)) then
1803 qfine = qfine_hysw
1804 qcoar = qcoar_hysw
1805 else if (nch .eq. waterptr_aer(m,itype)) then
1806 qfine = qfine_watr
1807 qcoar = qcoar_watr
1808 else if (nch .eq. numptr_aer(m,itype,iphase)) then
1809 dumvol1p = sqrt(volumlo_sect(m,itype)*volumhi_sect(m,itype))
1810 qfine = qfine_vol/dumvol1p
1811 qcoar = qcoar_vol/dumvol1p
1812 end if
1814 if ((qfine >= 0.0) .and. (qcoar >= 0.0)) then
1815 qval = qfine*fr_fine(m) + qcoar*fr_coar(m)
1816 chem_bv = qval*alt
1817 goto 2910
1818 end if
1820 2900 continue
1821 2910 continue
1823 return
1824 end subroutine bdy_chem_value_mosaic