module_mosaic_initmixrats.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 module module_mosaic_initmixrats
10
11 ! rce 2005-feb-18 - several fixes for indices of dlo/hi_sect, volumlo/hi_sect,
12 ! which are now (isize,itype)
13
14 ! rce 2004-dec-03 - many changes associated with the new aerosol "pointer"
15 ! variables in module_data_mosaic_asect
16
17 USE module_state_description
18
19 integer, parameter :: mosaic_init_wrf_mixrats_flagaa = 1
20 ! turns subr mosaic_init_wrf_mixrats on/off
21
22 contains
23
24
25 !-----------------------------------------------------------------------
26 !wig 10-May-2004, added phb, ph, g
27 subroutine mosaic_init_wrf_mixrats( &
28 iflagaa, config_flags, &
29 chem, alt, z_at_w, g, &
30 ids,ide, jds,jde, kds,kde, &
31 ims,ime, jms,jme, kms,kme, &
32 its,ite, jts,jte, kts,kte )
33
34 !
35 ! initializes the species and number mixing ratios for each section
36 !
37 ! this top level routine simply calls other routines depending on value
38 ! of config_flags%aer_ic_opt
39 !
40 use module_configure, only: grid_config_rec_type
41 use module_state_description, only: num_chem, param_first_scalar, &
42 aer_ic_pnnl
43 use module_data_mosaic_asect
44 use module_data_mosaic_other
45 use module_peg_util, only: peg_message, peg_error_fatal
46
47 implicit none
48
49
50 ! subr arguments
51 type(grid_config_rec_type), intent(in) :: config_flags
52
53 integer, intent(in) :: &
54 iflagaa, &
55 ids, ide, jds, jde, kds, kde, &
56 ims, ime, jms, jme, kms, kme, &
57 its, ite, jts, jte, kts, kte
58
59 real, intent(inout), &
60 dimension( ims:ime, kms:kme, jms:jme, 1:num_chem ) :: &
61 chem
62
63 !wig 10-May-2004, added phb, ph, g
64 real, intent(in), &
65 dimension( ims:ime, kms:kme, jms:jme ) :: &
66 alt, z_at_w
67 real, intent(in) :: g
68
69 ! local variables
70 integer :: i, ic, j, k
71
72 if (config_flags%aer_ic_opt == aer_ic_pnnl) then
73 !wig 10-May-2004, added phb, ph, g
74 call mosaic_init_wrf_mixrats_opt2( &
75 iflagaa, config_flags, &
76 chem, z_at_w, g, &
77 ids,ide, jds,jde, kds,kde, &
78 ims,ime, jms,jme, kms,kme, &
79 its,ite, jts,jte, kts,kte )
80
81 else
82 call mosaic_init_wrf_mixrats_opt1( &
83 iflagaa, config_flags, &
84 chem, &
85 ids,ide, jds,jde, kds,kde, &
86 ims,ime, jms,jme, kms,kme, &
87 its,ite, jts,jte, kts,kte )
88
89 end if
90
91 ! Aerosol species are returned from above in concentration units. Convert
92 ! them to mixing ratio for use in advection, etc.
93 do ic = p_so4_a01,num_chem
94 do j = jts,jte
95 do k = kts,kte-1
96 do i = its,ite
97 chem(i,k,j,ic) = chem(i,k,j,ic)*alt(i,k,j)
98 end do
99 end do
100 end do
101 end do
102
103 ! Fill the top z-staggered location to prevent trouble during advection.
104 do ic = p_so4_a01,num_chem
105 do j = jts,jte
106 do i = its,ite
107 chem(i,kte,j,ic) = chem(i,kte-1,j,ic)
108 end do
109 end do
110 end do
111
112 return
113 end subroutine mosaic_init_wrf_mixrats
114
115
116 !-----------------------------------------------------------------------
117 subroutine mosaic_init_wrf_mixrats_opt1( &
118 iflagaa, config_flags, &
119 chem, &
120 ids,ide, jds,jde, kds,kde, &
121 ims,ime, jms,jme, kms,kme, &
122 its,ite, jts,jte, kts,kte )
123
124 !
125 ! initializes the species and number mixing ratios for each section
126 ! based on user-specified lognormal modes that span the size distribution
127 !
128 ! rce 11-sep-2004 - eliminated use of the _wrfch pointers (lptr_xxx_a_wrfch,
129 ! lwaterptr_wrfch, numptr_wrfch); use only the _amode pointers now;
130 ! added l1dum
131 !
132 use module_configure, only: grid_config_rec_type
133 use module_state_description, only: num_chem, param_first_scalar
134 use module_data_mosaic_asect
135 use module_data_mosaic_other
136 use module_peg_util, only: peg_message, peg_error_fatal
137
138 implicit none
139
140 ! subr arguments
141 type(grid_config_rec_type), intent(in) :: config_flags
142
143 integer, intent(in) :: &
144 iflagaa, &
145 ids, ide, jds, jde, kds, kde, &
146 ims, ime, jms, jme, kms, kme, &
147 its, ite, jts, jte, kts, kte
148
149 real, intent(inout), &
150 dimension( ims:ime, kms:kme, jms:jme, 1:num_chem ) :: &
151 chem
152
153 ! local variables
154 integer i, j, k, l, ll, l1, l3, l1dum, m, mdum, nsm
155 integer it, jt, kt
156
157 real dum, dumdp, dumrsfc, dumvol, &
158 xlo, xhi, &
159 dumvol1p, &
160 pdummb, zdumkm, zscalekm, zfactor
161
162 real vtot_nsm_ofmode(maxd_asize)
163 real dumarr(maxd_acomp+5)
164
165 real erfc
166
167 ! double precision fracnum, fracvol, tlo, thi
168 real fracvol, tlo, thi
169
170 integer nmaxd_nsm
171 parameter (nmaxd_nsm = 4)
172
173 integer iphase, itype, ntot_nsm
174 integer iiprof_nsm(nmaxd_nsm)
175 integer lldum_so4, lldum_nh4, lldum_oc, lldum_bc, &
176 lldum_oin, lldum_na, lldum_cl, lldum_hysw
177
178 real sx_nsm(nmaxd_nsm), sxr2_nsm(nmaxd_nsm), &
179 x0_nsm(nmaxd_nsm), x3_nsm(nmaxd_nsm), &
180 rtot_nsm(maxd_acomp,nmaxd_nsm), &
181 vtot_nsm(nmaxd_nsm), xntot_nsm(nmaxd_nsm)
182
183 real dgnum_nsm(nmaxd_nsm), sigmag_nsm(nmaxd_nsm)
184 real aaprof_nsm(maxd_acomp+1,nmaxd_nsm)
185 real bbprof_nsm(nmaxd_nsm)
186
187 character*80 msg
188 character*10 dumname
189
190
191 ! check for on/off
192 if (mosaic_init_wrf_mixrats_flagaa .le. 0) return
193
194
195 ! *** currently only works for ntype_aer = 1
196 itype = 1
197 iphase = ai_phase
198 m = 1
199
200 ! set values for initialization modes
201 iiprof_nsm(:) = 1
202 aaprof_nsm(:,:) = 0.0
203 bbprof_nsm(:) = 0.0
204
205 ntot_nsm = 4
206 ntot_nsm = min( ntot_nsm, nsize_aer(itype) )
207
208 lldum_so4 = 0
209 lldum_nh4 = 0
210 lldum_oc = 0
211 lldum_bc = 0
212 lldum_oin = 0
213 lldum_na = 0
214 lldum_cl = 0
215 lldum_hysw = 0
216
217 do ll = 1, ncomp_plustracer_aer(itype)
218 if (massptr_aer(ll,m,itype,iphase) .eq. &
219 lptr_so4_aer(m,itype,iphase)) lldum_so4 = ll
220 if (massptr_aer(ll,m,itype,iphase) .eq. &
221 lptr_nh4_aer(m,itype,iphase)) lldum_nh4 = ll
222 if (massptr_aer(ll,m,itype,iphase) .eq. &
223 lptr_oc_aer(m,itype,iphase)) lldum_oc = ll
224 if (massptr_aer(ll,m,itype,iphase) .eq. &
225 lptr_bc_aer(m,itype,iphase)) lldum_bc = ll
226 if (massptr_aer(ll,m,itype,iphase) .eq. &
227 lptr_oin_aer(m,itype,iphase)) lldum_oin = ll
228 if (massptr_aer(ll,m,itype,iphase) .eq. &
229 lptr_na_aer(m,itype,iphase)) lldum_na = ll
230 if (massptr_aer(ll,m,itype,iphase) .eq. &
231 lptr_cl_aer(m,itype,iphase)) lldum_cl = ll
232 end do
233 if (hyswptr_aer(m,itype) .gt. 0) &
234 lldum_hysw = ncomp_plustracer_aer(itype) + 1
235
236 msg = ' '
237 if (lldum_so4 .le. 0) &
238 msg = '*** subr mosaic_init_wrf_mixrats - lldum_so4 = 0'
239 if (lldum_nh4 .le. 0) &
240 msg = '*** subr mosaic_init_wrf_mixrats - lldum_nh4 = 0'
241 if (lldum_oc .le. 0) &
242 msg = '*** subr mosaic_init_wrf_mixrats - lldum_oc = 0'
243 if (lldum_bc .le. 0) &
244 msg = '*** subr mosaic_init_wrf_mixrats - lldum_bc = 0'
245 if (lldum_oin .le. 0) &
246 msg = '*** subr mosaic_init_wrf_mixrats - lldum_oin = 0'
247 if (lldum_na .le. 0) &
248 msg = '*** subr mosaic_init_wrf_mixrats - lldum_na = 0'
249 if (lldum_cl .le. 0) &
250 msg = '*** subr mosaic_init_wrf_mixrats - lldum_cl = 0'
251 if (lldum_hysw .le. 0) &
252 msg = '*** subr mosaic_init_wrf_mixrats - lldum_hysw = 0'
253 if (msg .ne. ' ') call peg_error_fatal( lunerr, msg )
254
255
256 do nsm = 1, ntot_nsm
257
258 if (nsm .eq. 1) then
259 ! accum mode with so4, nh4, oc, bc
260 dgnum_nsm( nsm) = 0.15e-4
261 sigmag_nsm(nsm) = 2.0
262 aaprof_nsm(lldum_so4,nsm) = 0.50
263 aaprof_nsm(lldum_nh4,nsm) = aaprof_nsm(lldum_so4,nsm) &
264 * (mw_nh4_aer/mw_so4_aer)
265 aaprof_nsm(lldum_oc,nsm) = 0.25
266 aaprof_nsm(lldum_bc,nsm) = 0.05
267 aaprof_nsm(lldum_hysw,nsm) = aaprof_nsm(lldum_so4,nsm) * 0.2
268
269 else if (nsm .eq. 2) then
270 ! aitken mode with so4, nh4, oc, bc
271 dgnum_nsm( nsm) = 0.03e-4
272 sigmag_nsm(nsm) = 2.0
273 aaprof_nsm(lldum_so4,nsm) = 0.50 * 0.020
274 aaprof_nsm(lldum_nh4,nsm) = aaprof_nsm(lldum_so4,nsm) &
275 * (mw_nh4_aer/mw_so4_aer)
276 aaprof_nsm(lldum_oc,nsm) = 0.25 * 0.020
277 aaprof_nsm(lldum_bc,nsm) = 0.05 * 0.020
278 aaprof_nsm(lldum_hysw,nsm) = aaprof_nsm(lldum_so4,nsm) * 0.2
279
280 else if (nsm .eq. 3) then
281 ! coarse dust mode with oin
282 dgnum_nsm( nsm) = 1.0e-4
283 sigmag_nsm(nsm) = 2.0
284 aaprof_nsm(lldum_oin,nsm) = 0.5
285 aaprof_nsm(lldum_hysw,nsm) = aaprof_nsm( 9,nsm) * 1.0e-3
286
287 else if (nsm .eq. 4) then
288 ! coarse seasalt mode with na, cl
289 dgnum_nsm( nsm) = 2.0e-4
290 sigmag_nsm(nsm) = 2.0
291 aaprof_nsm(lldum_cl,nsm) = 0.1
292 aaprof_nsm(lldum_na,nsm) = aaprof_nsm(lldum_cl,nsm) &
293 * (mw_na_aer/mw_cl_aer)
294 aaprof_nsm(lldum_hysw,nsm) = aaprof_nsm(lldum_cl,nsm) * 0.2
295
296 end if
297
298 end do
299
300 ! when iflagaa = 1/2/3/4, use only the nsm-mode = iflagaa
301 if (iflagaa .gt. 0) then
302 do nsm = 1, ntot_nsm
303 if (nsm .ne. iflagaa) aaprof_nsm(:,nsm) = 0.0
304 end do
305 end if
306
307
308
309 !
310 ! do the initialization now
311 !
312
313 ! calculate mode parameters
314 do nsm = 1, ntot_nsm
315 sx_nsm(nsm) = alog( sigmag_nsm(nsm) )
316 sxr2_nsm(nsm) = sx_nsm(nsm) * sqrt(2.0)
317 x0_nsm(nsm) = alog( dgnum_nsm(nsm) )
318 x3_nsm(nsm) = x0_nsm(nsm) + 3.0*sx_nsm(nsm)*sx_nsm(nsm)
319 end do
320
321 ! initialize rclm array to zero
322 rclm(:,:) = 0.
323 ! rclmsvaa(:,:) = 0.
324
325 !
326 ! loop over all vertical levels
327 !
328 ! do 12900 k = 1, ktot
329 do 12900 k = 1, 1
330
331 ! pdummb = 1013.*scord(k)
332 ! zdumkm = ptoz( pdummb ) * 1.e-3
333 zdumkm = 0.0
334
335
336 ! for each species and nsm mode, define total mixing ratio
337 ! (for all sizes) at level k
338 !
339 ! iiprof_nsm = +1 gives constant mixing ratio
340 ! aaprof_nsm(l,nsm) = constant mixing ratio (ppb)
341 ! iiprof_nsm = +2 gives exponential profile
342 ! aaprof_nsm(l,nsm) = surface mixing ratio (ppb)
343 ! bbprof_nsm(l) = scale height (km)
344 ! iiprof_nsm = +3 gives linear profile (then zero at z > zscalekm)
345 ! aaprof_nsm(l,nsm) = surface mixing ratio (ppb)
346 ! bbprof_nsm(l) = height (km) at which mixing ratio = 0
347
348 do nsm = 1, ntot_nsm
349
350 if (iiprof_nsm(nsm) .eq. 2) then
351 zscalekm = bbprof_nsm(nsm)
352 zfactor = exp( -zdumkm/zscalekm )
353 else if (iiprof_nsm(nsm) .eq. 3) then
354 zscalekm = bbprof_nsm(nsm)
355 zfactor = max( 0., (1. - zdumkm/zscalekm) )
356 else
357 zfactor = 1.0
358 end if
359
360 do ll = 1, ncomp_plustracer_aer(itype) + 1
361 rtot_nsm(ll,nsm) = aaprof_nsm(ll,nsm) * zfactor
362 end do
363
364 end do
365
366 ! compute total volume and number mixing ratios for each nsm mode
367 ! rtot_nsm is ug/m3, 1.0e-6*rtot is g/m3, vtot_nsm is cm3/m3
368 do nsm = 1, ntot_nsm
369 dumvol = 0.
370 do ll = 1, ncomp_aer(itype)
371 dum = 1.0e-6*rtot_nsm(ll,nsm)/dens_aer(ll,itype)
372 dumvol = dumvol + max( 0., dum )
373 end do
374 vtot_nsm(nsm) = dumvol
375 end do
376
377 ! now compute species and number mixing ratios for each bin
378 do 12700 m = 1, nsize_aer(itype)
379
380 vtot_nsm_ofmode(m) = 0.0
381
382 do 12500 nsm = 1, ntot_nsm
383
384 ! for nsm_mode = n, compute fraction of number and volume
385 ! that is in bin m
386 xlo = alog( dlo_sect(m,itype) )
387 xhi = alog( dhi_sect(m,itype) )
388
389 tlo = (xlo - x3_nsm(nsm))/sxr2_nsm(nsm)
390 thi = (xhi - x3_nsm(nsm))/sxr2_nsm(nsm)
391 if (tlo .ge. 0.) then
392 ! fracvol = 0.5*( erfc(tlo) - erfc(thi) )
393 fracvol = 0.5*( erfc_num_recipes(tlo) - erfc_num_recipes(thi) )
394 else
395 ! fracvol = 0.5*( erfc(-thi) - erfc(-tlo) )
396 fracvol = 0.5*( erfc_num_recipes(-thi) - erfc_num_recipes(-tlo) )
397 end if
398 fracvol = max( fracvol, 0.0 )
399
400 vtot_nsm_ofmode(m) = vtot_nsm_ofmode(m) + vtot_nsm(nsm)*fracvol
401
402 ! now load that fraction of species-mixing-ratio
403 ! into the appropriate rclm location
404 do ll = 1, ncomp_plustracer_aer(itype)
405 rclm( k, massptr_aer(ll,m,itype,iphase) ) = &
406 rclm( k, massptr_aer(ll,m,itype,iphase) ) + &
407 fracvol*rtot_nsm(ll,nsm)
408 end do
409
410 if ((iphase .eq. ai_phase) .and. &
411 (lldum_hysw .gt. 0) .and. &
412 (hyswptr_aer(m,itype) .gt. 0)) then
413
414 rclm( k, hyswptr_aer(m,itype) ) = &
415 rclm( k, hyswptr_aer(m,itype) ) + &
416 fracvol*rtot_nsm(lldum_hysw,nsm)
417 end if
418
419 12500 continue
420
421 ! now compute number from volume
422 dum = sqrt( dlo_sect(m,itype)*dhi_sect(m,itype) )
423 dumvol1p = (pi/6.0)*(dum**3)
424 rclm( k, numptr_aer(m,itype,iphase) ) = vtot_nsm_ofmode(m)/dumvol1p
425
426 ! set water = hyswatr
427 if ((iphase .eq. ai_phase) .and. &
428 (lldum_hysw .gt. 0) .and. &
429 (hyswptr_aer(m,itype) .gt. 0) .and. &
430 (waterptr_aer(m,itype) .gt. 0)) then
431
432 rclm( k, waterptr_aer(m,itype) ) = &
433 rclm( k, hyswptr_aer(m,itype) )
434 end if
435
436 12700 continue
437
438 12900 continue
439
440
441 !
442 ! do diagnostic output
443 !
444
445 ! temporary
446 ! temporary
447
448 dumarr(:) = 0.0
449 msg = ' '
450 call peg_message( lunout, msg )
451 msg = '*** subr mosaic_init_wrf_mixrats_opt1 results'
452 call peg_message( lunout, msg )
453 msg = ' mass in ug/m3 number in #/m3 volume in cm3/m3'
454 call peg_message( lunout, msg )
455 msg = ' '
456 call peg_message( lunout, msg )
457 msg = ' mode l l1 species conc'
458 call peg_message( lunout, msg )
459
460 do 14390 mdum = 1, nsize_aer(itype)+1
461 m = min( mdum, nsize_aer(itype) )
462 msg = ' '
463 call peg_message( lunout, msg )
464 do 14350 l = 1, ncomp_plustracer_aer(itype)+4
465
466 if (l .le. ncomp_plustracer_aer(itype)) then
467 l1 = massptr_aer(l,m,itype,iphase)
468 dumname = name_aer(l,itype)
469 dum = rclm(1,l1)
470 else if (l .eq. ncomp_plustracer_aer(itype)+1) then
471 l1 = hyswptr_aer(m,itype)
472 dumname = 'hystwatr'
473 dum = rclm(1,l1)
474 else if (l .eq. ncomp_plustracer_aer(itype)+2) then
475 l1 = waterptr_aer(m,itype)
476 dumname = 'water'
477 dum = rclm(1,l1)
478 else if (l .eq. ncomp_plustracer_aer(itype)+3) then
479 l1 = numptr_aer(m,itype,iphase)
480 dumname = 'number'
481 dum = rclm(1,l1)
482 else if (l .eq. ncomp_plustracer_aer(itype)+4) then
483 l1 = 0
484 dumname = 'volume'
485 dum = vtot_nsm_ofmode(m)
486 else
487 dumname = '????????'
488 l1 = -1
489 dum = -1.0
490 end if
491
492 l1dum = l1
493 if (aboxtest_testmode .gt. 0) l1dum = max( l1-1, 0 )
494
495 if (mdum .le. nsize_aer(itype)) then
496 dumarr(l) = dumarr(l) + dum
497 write(msg,9620) m, l, l1dum, dumname, dum
498 else
499 write(msg,9625) l, dumname, dumarr(l)
500 end if
501 call peg_message( lunout, msg )
502
503 14350 continue
504 14390 continue
505
506 9620 format( 3i4, 2x, a, 3(1pe12.3) )
507 9625 format( ' sum', i4, ' -', 2x, a, 3(1pe12.3) )
508
509
510 !
511 ! load the chem array
512 !
513 do 16390 m = 1, nsize_aer(itype)
514 do 16350 l = 1, 15
515
516 if (l .eq. 1) then
517 l1 = lptr_so4_aer(m,itype,iphase)
518 else if (l .eq. 2) then
519 l1 = lptr_no3_aer(m,itype,iphase)
520 else if (l .eq. 3) then
521 l1 = lptr_cl_aer(m,itype,iphase)
522 else if (l .eq. 4) then
523 l1 = lptr_msa_aer(m,itype,iphase)
524 else if (l .eq. 5) then
525 l1 = lptr_co3_aer(m,itype,iphase)
526 else if (l .eq. 6) then
527 l1 = lptr_nh4_aer(m,itype,iphase)
528 else if (l .eq. 7) then
529 l1 = lptr_na_aer(m,itype,iphase)
530 else if (l .eq. 8) then
531 l1 = lptr_ca_aer(m,itype,iphase)
532 else if (l .eq. 9) then
533 l1 = lptr_oin_aer(m,itype,iphase)
534 else if (l .eq. 10) then
535 l1 = lptr_oc_aer(m,itype,iphase)
536 else if (l .eq. 11) then
537 l1 = lptr_bc_aer(m,itype,iphase)
538 else if (l .eq. 12) then
539 l1 = hyswptr_aer(m,itype)
540 else if (l .eq. 13) then
541 l1 = waterptr_aer(m,itype)
542 else if (l .eq. 14) then
543 l1 = numptr_aer(m,itype,iphase)
544 else
545 goto 16350
546 end if
547 l3 = l1
548
549 if ((l1 .gt. 0) .and. (l1 .le. lmaxd) .and. &
550 (l3 .ge. param_first_scalar)) then
551 do it = its, ite
552 ! *** note: not sure what the kt limits should be here
553 do kt = kts, kte-1
554 do jt = jts, jte
555 chem(it,kt,jt,l3) = rclm(1,l1)
556 end do
557 end do
558 end do
559 end if
560
561 16350 continue
562 16390 continue
563
564
565 ! all done
566 return
567 end subroutine mosaic_init_wrf_mixrats_opt1
568
569
570 !-----------------------------------------------------------------------
571 !wig 10-May-2004, added phb, ph, g
572 subroutine mosaic_init_wrf_mixrats_opt2( &
573 iflagaa, config_flags, &
574 chem, z_at_w, g, &
575 ids,ide, jds,jde, kds,kde, &
576 ims,ime, jms,jme, kms,kme, &
577 its,ite, jts,jte, kts,kte )
578
579 !
580 ! provides initial values for mosaic aerosol species (mass and number
581 ! mixing ratio) for "Texas August 2000" simulations
582 ! currently all the initial values are uniform in x, y, AND z
583 !
584 ! rce 11-sep-2004 - eliminated use of the _wrfch pointers (lptr_xxx_a_wrfch,
585 ! lwaterptr_wrfch, numptr_wrfch); use only the _amode pointers now
586 !
587 use module_configure, only: grid_config_rec_type
588 use module_state_description, only: num_chem, param_first_scalar
589 use module_data_mosaic_asect
590 use module_data_mosaic_other
591 use module_peg_util, only: peg_message, peg_error_fatal
592
593 implicit none
594
595 ! subr arguments
596 type(grid_config_rec_type), intent(in) :: config_flags
597
598 integer, intent(in) :: &
599 iflagaa, &
600 ids, ide, jds, jde, kds, kde, &
601 ims, ime, jms, jme, kms, kme, &
602 its, ite, jts, jte, kts, kte
603
604 real, intent(inout), &
605 dimension( ims:ime, kms:kme, jms:jme, 1:num_chem ) :: &
606 chem
607
608 !wig 10-May-2004, added phb, ph, g
609 real, intent(in), &
610 dimension( ims:ime, kms:kme, jms:jme ) :: &
611 z_at_w
612 real :: g
613
614 ! local variables
615 integer l, l1, l3, m, mdum
616 integer iphase, itype
617 integer it, jt, kt
618
619 !wig 10-May-2004, added mult
620 real dum, dumvol1p, mult
621 real qcoar, qfine, qval
622
623 real :: vtot_ofmode(maxd_asize)
624 real :: dumarr(maxd_acomp+5)
625 real :: fr_coar(8), fr_fine(8)
626
627 !wig 01-Apr-2005, Updated fractional breakdown between bins. (See also
628 ! bdy_chem_value_mosaic in this file and mosaic_addemiss in
629 ! module_mosaic_addemiss.F) Note that the values here no
630 ! longer match those in mosaic_addemiss.
631 !rce 10-May-2005, changed fr8b_coar to values determined by jdf
632 real, save :: fr8b_coar(8) = &
633 (/ 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.300, 0.700 /) ! 10-May-2005
634 ! (/ 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.933, 0.067 /) ! 01-Apr-2005
635 ! (/ 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.16, 0.84 /) ! "old"
636
637 real, save :: fr8b_fine(8) = &
638 (/ 0.006, 0.024, 0.231, 0.341, 0.140, 0.258, 0., 0./) !5-Apr-2005 values
639 ! (/ 0.0275, 0.0426, 0.2303, 0.3885, 0.1100, 0.2011, 0., 0./) !15-Nov-2004 values
640 ! (/ 0.01, 0.05, 0.145, 0.60, 0.145, 0.05, 0.00, 0.00 /)
641 ! (/ 0.04, 0.10, 0.35, 0.29, 0.15, 0.07, 0.00, 0.00 /)
642
643 real :: qfine_so4, qfine_no3, qfine_cl, qfine_msa, &
644 qfine_co3, qfine_nh4, qfine_na, qfine_ca, qfine_oin, &
645 qfine_oc, qfine_bc, qfine_hysw, qfine_watr, qfine_vol
646 real :: qcoar_so4, qcoar_no3, qcoar_cl, qcoar_msa, &
647 qcoar_co3, qcoar_nh4, qcoar_na, qcoar_ca, qcoar_oin, &
648 qcoar_oc, qcoar_bc, qcoar_hysw, qcoar_watr, qcoar_vol
649
650 !wig 10-May-2004, added z
651 real, dimension( ims:ime, kms:kme, jms:jme ) :: z
652
653 character*80 msg
654 character*10 dumname
655
656
657 ! *** currently only works for ntype_aer = 1
658 itype = 1
659 iphase = ai_phase
660 m = 1
661
662 !wig 10-May-2004, block begin...
663 ! calculate the mid-level heights
664 do it = its, ite
665 do kt = kts, kte-1
666 do jt = jts, jte
667 z(it,kt,jt) = (z_at_w(it,kt,jt)+z_at_w(it,kt+1,jt))*0.5
668 end do
669 end do
670 end do
671 !wig 10-May-2004, ...end block
672
673 ! set the partitioning fractions for either 8 or 4 bins
674 if (nsize_aer(itype) == 8) then
675 fr_coar(:) = fr8b_coar(:)
676 fr_fine(:) = fr8b_fine(:)
677 else if (nsize_aer(itype) == 4) then
678 do m = 1, nsize_aer(itype)
679 fr_coar(m) = fr8b_coar(2*m) + fr8b_coar(2*m-1)
680 fr_fine(m) = fr8b_fine(2*m) + fr8b_fine(2*m-1)
681 end do
682 else
683 write(msg,'(a,i5)') &
684 'subr mosaic_init_wrf_mixrats_opt2' // &
685 ' - nsize_aer(itype) must be 4 or 8 but = ', nsize_aer(itype)
686 call peg_error_fatal( lunout, msg )
687 end if
688
689 !
690 ! compute initial values (currently uniform in x, y, AND z)
691 ! load them into the chem array
692 ! also load them into the rclm array for diagnostic output
693 !
694 rclm(:,:) = 0.0
695 vtot_ofmode(:) = 0.0
696
697 ! Set values for fine and coarse mass concentrations, and
698 ! compute fine/coarse volume concentrations. These are also set in
699 ! bdy_chem_value_mosaic, below.
700 ! (these are latest values from 1-Apr-2005 discussions)
701 !
702 ! rce 10-may-2005 - changed qfine_cl, _na, _oin to values determined by jdf
703 qfine_so4 = 2.14
704 qcoar_so4 = 0.242
705 qfine_no3 = 0.11
706 qcoar_no3 = 0.03
707 ! qfine_cl = 0.3
708 qfine_cl = 0.14 ! 10-May-2005
709 qcoar_cl = 0.139
710 qfine_msa = 0.0
711 qcoar_msa = 0.0
712 qfine_co3 = 0.0
713 qcoar_co3 = 0.0
714 qfine_nh4 = 0.83
715 qcoar_nh4 = 0.10
716 ! qfine_na = 0.2
717 qfine_na = 0.1 ! 10-May-2005
718 qcoar_na = 0.09
719 qfine_ca = 0.0
720 qcoar_ca = 0.0
721 ! qfine_oin = 6.2
722 qfine_oin = 3.48 ! 10-May-2005
723 qcoar_oin = 0.35
724 qfine_oc = 1.00
725 qcoar_oc = 1.50
726 qfine_bc = 0.2
727 qcoar_bc = 0.075
728 qfine_hysw = 0.0
729 qcoar_hysw = 0.0
730 qfine_watr = 0.0
731 qcoar_watr = 0.0
732
733 !!$! old values from 23-Apr-2004:
734 !!$ qfine_so4 = 2.554
735 !!$ qcoar_so4 = 0.242
736 !!$ qfine_no3 = 0.07
737 !!$ qcoar_no3 = 0.03
738 !!$ qfine_cl = 0.324
739 !!$ qcoar_cl = 0.139
740 !!$ qfine_msa = 0.0
741 !!$ qcoar_msa = 0.0
742 !!$ qfine_co3 = 0.0
743 !!$ qcoar_co3 = 0.0
744 !!$ qfine_nh4 = 0.98
745 !!$ qcoar_nh4 = 0.10
746 !!$ qfine_na = 0.21
747 !!$ qcoar_na = 0.09
748 !!$ qfine_ca = 0.0
749 !!$ qcoar_ca = 0.0
750 !!$ qfine_oin = 0.15
751 !!$ qcoar_oin = 0.35
752 !!$ qfine_oc = 1.00
753 !!$ qcoar_oc = 1.50
754 !!$ qfine_bc = 0.175
755 !!$ qcoar_bc = 0.075
756 !!$ qfine_hysw = 0.0
757 !!$ qcoar_hysw = 0.0
758 !!$ qfine_watr = 0.0
759 !!$ qcoar_watr = 0.0
760
761
762 ! qfine_so4 ... are ug/m3 so 1.0e-6 factor gives g/m3
763 ! dens_so4 ... are g/cm3; qfine_vol ... are cm3/m3
764 qfine_vol = 1.0e-6 * ( &
765 (qfine_so4/dens_so4_aer) + (qfine_no3/dens_no3_aer) + &
766 (qfine_cl /dens_cl_aer ) + (qfine_msa/dens_msa_aer) + &
767 (qfine_co3/dens_co3_aer) + (qfine_nh4/dens_nh4_aer) + &
768 (qfine_na /dens_na_aer ) + (qfine_ca /dens_ca_aer ) + &
769 (qfine_oin/dens_oin_aer) + (qfine_oc /dens_oc_aer ) + &
770 (qfine_bc /dens_bc_aer ) )
771 qcoar_vol = 1.0e-6 * ( &
772 (qcoar_so4/dens_so4_aer) + (qcoar_no3/dens_no3_aer) + &
773 (qcoar_cl /dens_cl_aer ) + (qcoar_msa/dens_msa_aer) + &
774 (qcoar_co3/dens_co3_aer) + (qcoar_nh4/dens_nh4_aer) + &
775 (qcoar_na /dens_na_aer ) + (qcoar_ca /dens_ca_aer ) + &
776 (qcoar_oin/dens_oin_aer) + (qcoar_oc /dens_oc_aer ) + &
777 (qcoar_bc /dens_bc_aer ) )
778
779 do 2900 m = 1, nsize_aer(itype)
780 do 2800 l = 1, 15
781
782 if (l .eq. 1) then
783 l1 = lptr_so4_aer(m,itype,iphase)
784 qfine = qfine_so4
785 qcoar = qcoar_so4
786 else if (l .eq. 2) then
787 l1 = lptr_no3_aer(m,itype,iphase)
788 qfine = qfine_no3
789 qcoar = qcoar_no3
790 else if (l .eq. 3) then
791 l1 = lptr_cl_aer(m,itype,iphase)
792 qfine = qfine_cl
793 qcoar = qcoar_cl
794 else if (l .eq. 4) then
795 l1 = lptr_msa_aer(m,itype,iphase)
796 qfine = qfine_msa
797 qcoar = qcoar_msa
798 else if (l .eq. 5) then
799 l1 = lptr_co3_aer(m,itype,iphase)
800 qfine = qfine_co3
801 qcoar = qcoar_co3
802 else if (l .eq. 6) then
803 l1 = lptr_nh4_aer(m,itype,iphase)
804 qfine = qfine_nh4
805 qcoar = qcoar_nh4
806 else if (l .eq. 7) then
807 l1 = lptr_na_aer(m,itype,iphase)
808 qfine = qfine_na
809 qcoar = qcoar_na
810 else if (l .eq. 8) then
811 l1 = lptr_ca_aer(m,itype,iphase)
812 qfine = qfine_ca
813 qcoar = qcoar_ca
814 else if (l .eq. 9) then
815 l1 = lptr_oin_aer(m,itype,iphase)
816 qfine = qfine_oin
817 qcoar = qcoar_oin
818 else if (l .eq. 10) then
819 l1 = lptr_oc_aer(m,itype,iphase)
820 qfine = qfine_oc
821 qcoar = qcoar_oc
822 else if (l .eq. 11) then
823 l1 = lptr_bc_aer(m,itype,iphase)
824 qfine = qfine_bc
825 qcoar = qcoar_bc
826 else if (l .eq. 12) then
827 l1 = hyswptr_aer(m,itype)
828 qfine = qfine_hysw
829 qcoar = qcoar_hysw
830 else if (l .eq. 13) then
831 l1 = waterptr_aer(m,itype)
832 qfine = qfine_watr
833 qcoar = qcoar_watr
834 else if (l .eq. 14) then
835 l1 = numptr_aer(m,itype,iphase)
836 dumvol1p = sqrt(volumlo_sect(m,itype)*volumhi_sect(m,itype))
837 qfine = qfine_vol/dumvol1p
838 qcoar = qcoar_vol/dumvol1p
839 vtot_ofmode(m) = &
840 qfine_vol*fr_fine(m) + qcoar_vol*fr_coar(m)
841 else
842 goto 2800
843 end if
844 l3 = l1
845
846 if ((l1 .gt. 0) .and. (l1 .le. lmaxd) .and. &
847 (l3 .ge. param_first_scalar)) then
848 qval = qfine*fr_fine(m) + qcoar*fr_coar(m)
849 rclm(1,l1) = qval
850
851 do it = its, ite
852 do kt = kts, kte-1
853 do jt = jts, jte
854
855 !wig 28-May-2004, begin block...
856 ! Determine height multiplier...
857 ! This should mimic the calculation in sorgam_set_aer_bc_pnnl,
858 ! sorgam_init_aer_ic_pnnl, bdy_chem_value_mosaic
859 !!$! Height(m) Multiplier
860 !!$! --------- ----------
861 !!$! <=2000 1.0
862 !!$! 2000<z<3000 linear transition zone to 0.5
863 !!$! 3000<z<5000 linear transision zone to 0.25
864 !!$! >=5000 0.25
865 !!$!
866 !!$! which translates to:
867 !!$! 2000<z<3000 mult = 1.0 + (z-2000.)*(0.5-1.0)/(3000.-2000.)
868 !!$! 3000<z<5000 mult = 0.5 + (z-3000.)*(0.25-0.5)/(5000.-3000.)
869 !!$! or in reduced form:
870 !!$ if( z(it,kt,jt) <= 2000. ) then
871 !!$ mult = 1.0
872 !!$ elseif( z(it,kt,jt) > 2000. &
873 !!$ .and. z(it,kt,jt) <= 3000. ) then
874 !!$ mult = 1.0 - 0.0005*(z(it,kt,jt)-2000.)
875 !!$ elseif( z(it,kt,jt) > 3000. &
876 !!$ .and. z(it,kt,jt) <= 5000. ) then
877 !!$ mult = 0.5 - 1.25e-4*(z(it,kt,jt)-3000.)
878 !!$ else
879 !!$ mult = 0.25
880 !!$ end if
881 ! Updated 1-Apr-2005:
882 ! Height(m) Multiplier
883 ! --------- ----------
884 ! <=2000 1.0
885 ! 2000<z<3000 linear transition zone to 0.25
886 ! 3000<z<5000 linear transision zone to 0.125
887 ! >=5000 0.125
888 !
889 ! which translates to:
890 ! 2000<z<3000 mult = 1.00 + (z-2000.)*(0.25-1.0)/(3000.-2000.)
891 ! 3000<z<5000 mult = 0.25 + (z-3000.)*(0.125-0.25)/(5000.-3000.)
892 ! or in reduced form:
893 if( z(it,kt,jt) <= 2000. ) then
894 mult = 1.0
895 elseif( z(it,kt,jt) > 2000. &
896 .and. z(it,kt,jt) <= 3000. ) then
897 mult = 1.0 - 0.00075*(z(it,kt,jt)-2000.)
898 elseif( z(it,kt,jt) > 3000. &
899 .and. z(it,kt,jt) <= 5000. ) then
900 mult = 0.25 - 4.166666667e-5*(z(it,kt,jt)-3000.)
901 else
902 mult = 0.125
903 end if
904
905 chem(it,kt,jt,l3) = mult*rclm(1,l1)
906 !wig 28-May-2004, ...end block
907 ! chem(it,kt,jt,l3) = rclm(1,l1)
908 end do
909 end do
910 end do
911 end if
912
913
914 2800 continue
915 2900 continue
916
917 !
918 ! do diagnostic output
919 !
920 dumarr(:) = 0.0
921 msg = ' '
922 call peg_message( lunout, msg )
923 msg = '*** subr mosaic_init_wrf_mixrats_opt2 results'
924 call peg_message( lunout, msg )
925 msg = ' mass in ug/m3 number in #/m3 volume in cm3/m3'
926 call peg_message( lunout, msg )
927 msg = ' '
928 call peg_message( lunout, msg )
929 msg = ' mode l l1 species conc'
930 call peg_message( lunout, msg )
931
932 do 3190 mdum = 1, nsize_aer(itype)+1
933 m = min( mdum, nsize_aer(itype) )
934 msg = ' '
935 call peg_message( lunout, msg )
936 do 3150 l = 1, ncomp_plustracer_aer(itype)+4
937
938 if (l .le. ncomp_plustracer_aer(itype)) then
939 l1 = massptr_aer(l,m,itype,iphase)
940 dumname = name_aer(l,itype)
941 dum = rclm(1,l1)
942 else if (l .eq. ncomp_plustracer_aer(itype)+1) then
943 l1 = hyswptr_aer(m,itype)
944 dumname = 'hystwatr'
945 dum = rclm(1,l1)
946 else if (l .eq. ncomp_plustracer_aer(itype)+2) then
947 l1 = waterptr_aer(m,itype)
948 dumname = 'water'
949 dum = rclm(1,l1)
950 else if (l .eq. ncomp_plustracer_aer(itype)+3) then
951 l1 = numptr_aer(m,itype,iphase)
952 dumname = 'number'
953 dum = rclm(1,l1)
954 else if (l .eq. ncomp_plustracer_aer(itype)+4) then
955 l1 = 0
956 dumname = 'volume'
957 dum = vtot_ofmode(m)
958 else
959 dumname = '????????'
960 l1 = -1
961 dum = -1.0
962 end if
963
964 if (mdum .le. nsize_aer(itype)) then
965 dumarr(l) = dumarr(l) + dum
966 write(msg,9620) m, l, l1, dumname, dum
967 else
968 write(msg,9625) l, dumname, dumarr(l)
969 end if
970 call peg_message( lunout, msg )
971
972 3150 continue
973 3190 continue
974
975 9620 format( 3i4, 2x, a, 3(1pe12.3) )
976 9625 format( ' sum', i4, ' -', 2x, a, 3(1pe12.3) )
977
978
979 ! all done
980 return
981 end subroutine mosaic_init_wrf_mixrats_opt2
982
983
984 !-----------------------------------------------------------------------
985 real function erfc_num_recipes( x )
986 !
987 ! from press et al, numerical recipes, 1990, page 164
988 !
989 implicit none
990 real x
991 double precision erfc_dbl, dum, t, zz
992
993 zz = abs(x)
994 t = 1.0/(1.0 + 0.5*zz)
995
996 ! erfc_num_recipes =
997 ! & t*exp( -zz*zz - 1.26551223 + t*(1.00002368 + t*(0.37409196 +
998 ! & t*(0.09678418 + t*(-0.18628806 + t*(0.27886807 +
999 ! & t*(-1.13520398 +
1000 ! & t*(1.48851587 + t*(-0.82215223 + t*0.17087277 )))))))))
1001
1002 dum = ( -zz*zz - 1.26551223 + t*(1.00002368 + t*(0.37409196 + &
1003 t*(0.09678418 + t*(-0.18628806 + t*(0.27886807 + &
1004 t*(-1.13520398 + &
1005 t*(1.48851587 + t*(-0.82215223 + t*0.17087277 )))))))))
1006
1007 erfc_dbl = t * exp(dum)
1008 if (x .lt. 0.0) erfc_dbl = 2.0d0 - erfc_dbl
1009
1010 erfc_num_recipes = erfc_dbl
1011
1012 return
1013 end function erfc_num_recipes
1014
1015
1016 !-----------------------------------------------------------------------
1017 end module module_mosaic_initmixrats
1018
1019
1020
1021
1022 !-----------------------------------------------------------------------
1023 subroutine bdy_chem_value_mosaic ( chem_bv, alt, zz, nch, config_flags )
1024 !
1025 ! provides boundary values for the mosaic aerosol species
1026 !
1027 ! it is outside of the module declaration because of potential
1028 ! module1 --> module2 --> module1 use conflicts
1029 !
1030 ! rce 11-sep-2004 - eliminated use of the _wrfch pointers (lptr_xxx_a_wrfch,
1031 ! lwaterptr_wrfch, numptr_wrfch); use only the _amode pointers now
1032 !
1033 use module_configure, only: grid_config_rec_type
1034 use module_state_description, only: param_first_scalar, &
1035 aer_bc_pnnl
1036 use module_data_mosaic_asect
1037 use module_data_mosaic_other
1038 implicit none
1039
1040 ! arguments
1041 REAL, intent(OUT) :: chem_bv ! boundary value for chem(-,-,-,nch)
1042 REAL, intent(IN) :: alt ! inverse density
1043 REAL, intent(IN) :: zz ! height
1044 INTEGER, intent(IN) :: nch ! index number of chemical species
1045 TYPE (grid_config_rec_type), intent(in) :: config_flags
1046
1047 ! local variables
1048 integer iphase, itype, m
1049
1050 real, parameter :: chem_bv_def = 1.0e-20
1051 !wig 10-May-2004, added mult
1052 real dumvol1p, mult
1053 real qcoar, qfine, qval
1054
1055 real :: fr_coar(8), fr_fine(8)
1056
1057 !wig 1-Apr-2005, Updated fractional breakdown between bins. (See also
1058 ! mosaic_init_wrf_mixrats_opt2, above, and mosaic_addemiss
1059 ! in module_mosaic_addemiss.F). Note that these values no
1060 ! longer match those in mosaic_addemiss.
1061 real, save :: fr8b_coar(8) = &
1062 (/ 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.933, 0.067 /)
1063 ! (/ 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.16, 0.84 /)
1064 real, save :: fr8b_fine(8) = &
1065 (/ 0.006, 0.024, 0.231, 0.341, 0.140, 0.258, 0., 0./) !5-Apr-2005 values
1066 ! (/ 0.0275, 0.0426, 0.2303, 0.3885, 0.1100, 0.2011, 0., 0./) !15-Nov-2004 values
1067 ! (/ 0.01, 0.05, 0.145, 0.60, 0.145, 0.05, 0.00, 0.00 /)
1068 ! (/ 0.04, 0.10, 0.35, 0.29, 0.15, 0.07, 0.00, 0.00 /)
1069
1070 real :: qfine_so4, qfine_no3, qfine_cl, qfine_msa, &
1071 qfine_co3, qfine_nh4, qfine_na, qfine_ca, qfine_oin, &
1072 qfine_oc, qfine_bc, qfine_hysw, qfine_watr, qfine_vol
1073 real :: qcoar_so4, qcoar_no3, qcoar_cl, qcoar_msa, &
1074 qcoar_co3, qcoar_nh4, qcoar_na, qcoar_ca, qcoar_oin, &
1075 qcoar_oc, qcoar_bc, qcoar_hysw, qcoar_watr, qcoar_vol
1076
1077 character*80 msg
1078
1079
1080
1081 ! when aer_bc_opt /= aer_bc_pnnl,
1082 ! just chem_bv=near-zero for all species
1083 chem_bv = chem_bv_def
1084 if (config_flags%aer_bc_opt /= aer_bc_pnnl) return
1085 if (nch < param_first_scalar) return
1086
1087
1088 ! *** currently only works for ntype_aer = 1
1089 itype = 1
1090 iphase = ai_phase
1091 m = 1
1092
1093
1094 !
1095 ! following is for aer_bc_opt == aer_bc_pnnl
1096 ! when boundary values are set for Texas August 2000 simulations
1097 !
1098 ! set the partitioning fractions for either 8 or 4 bins
1099 if (nsize_aer(itype) == 8) then
1100 fr_coar(:) = fr8b_coar(:)
1101 fr_fine(:) = fr8b_fine(:)
1102 else if (nsize_aer(itype) == 4) then
1103 do m = 1, nsize_aer(itype)
1104 fr_coar(m) = fr8b_coar(2*m) + fr8b_coar(2*m-1)
1105 fr_fine(m) = fr8b_fine(2*m) + fr8b_fine(2*m-1)
1106 end do
1107 else
1108 write(msg,'(a,i5)') &
1109 'subr bdy_chem_value_mosaic' // &
1110 ' - nsize_aer(itype) must be 4 or 8 but = ', nsize_aer(itype)
1111 call wrf_error_fatal( msg )
1112 end if
1113
1114 ! Determine height multiplier...
1115 ! This should mimic the calculation in sorgam_set_aer_bc_pnnl,
1116 ! sorgam_init_aer_ic_pnnl, and mosaic_init_wrf_mixrats_opt2
1117 !!$! Height(m) Multiplier
1118 !!$! --------- ----------
1119 !!$! <=2000 1.0
1120 !!$! 2000<z<3000 linear transition zone to 0.5
1121 !!$! 3000<z<5000 linear transision zone to 0.25
1122 !!$! >=5000 0.25
1123 !!$!
1124 !!$! which translates to:
1125 !!$! 2000<z<3000 mult = 1.0 + (z-2000.)*(0.5-1.0)/(3000.-2000.)
1126 !!$! 3000<z<5000 mult = 0.5 + (z-3000.)*(0.25-0.5)/(5000.-3000.)
1127 !!$! or in reduced form:
1128 !!$ if( zz <= 2000. ) then
1129 !!$ mult = 1.0
1130 !!$ elseif( zz > 2000. &
1131 !!$ .and. zz <= 3000. ) then
1132 !!$ mult = 1.0 - 0.0005*(zz-2000.)
1133 !!$ elseif( zz > 3000. &
1134 !!$ .and. zz <= 5000. ) then
1135 !!$ mult = 0.5 - 1.25e-4*(zz-3000.)
1136 !!$ else
1137 !!$ mult = 0.25
1138 !!$ end if
1139 ! Updated aerosol profile multiplier 1-Apr-2005:
1140 ! Height(m) Multiplier
1141 ! --------- ----------
1142 ! <=2000 1.0
1143 ! 2000<z<3000 linear transition zone to 0.25
1144 ! 3000<z<5000 linear transision zone to 0.125
1145 ! >=5000 0.125
1146 !
1147 ! which translates to:
1148 ! 2000<z<3000 mult = 1.00 + (z-2000.)*(0.25-1.0)/(3000.-2000.)
1149 ! 3000<z<5000 mult = 0.25 + (z-3000.)*(0.125-0.25)/(5000.-3000.)
1150 ! or in reduced form:
1151 if( zz <= 2000. ) then
1152 mult = 1.0
1153 elseif( zz > 2000. &
1154 .and. zz <= 3000. ) then
1155 mult = 1.0 - 0.00075*(zz-2000.)
1156 elseif( zz > 3000. &
1157 .and. zz <= 5000. ) then
1158 mult = 0.25 - 4.166666667e-5*(zz-3000.)
1159 else
1160 mult = 0.125
1161 end if
1162
1163 ! Set values for fine and coarse mass concentrations, and
1164 ! compute fine/coarse volume concentrations. These are also set in
1165 ! mosaic_init_wrf_mixrats_opt2.
1166 ! (these are latest values from 1-Apr-2005 discussions)
1167 !wig 10-May-2004, added height multiplier (mult*) to each species...
1168 qfine_so4 = mult*2.14
1169 qcoar_so4 = mult*0.242
1170 qfine_no3 = mult*0.11
1171 qcoar_no3 = mult*0.03
1172 qfine_cl = mult*0.3
1173 qcoar_cl = mult*0.139
1174 qfine_msa = mult*0.0
1175 qcoar_msa = mult*0.0
1176 qfine_co3 = mult*0.0
1177 qcoar_co3 = mult*0.0
1178 qfine_nh4 = mult*0.83
1179 qcoar_nh4 = mult*0.10
1180 qfine_na = mult*0.2
1181 qcoar_na = mult*0.09
1182 qfine_ca = mult*0.0
1183 qcoar_ca = mult*0.0
1184 qfine_oin = mult*6.2
1185 qcoar_oin = mult*0.35
1186 qfine_oc = mult*1.00
1187 qcoar_oc = mult*1.50
1188 qfine_bc = mult*0.2
1189 qcoar_bc = mult*0.075
1190 qfine_hysw = mult*0.0
1191 qcoar_hysw = mult*0.0
1192 qfine_watr = mult*0.0
1193 qcoar_watr = mult*0.0
1194 !!$! old values from 23-Apr-2004:
1195 !!$ qfine_so4 = mult*2.554
1196 !!$ qcoar_so4 = mult*0.242
1197 !!$ qfine_no3 = mult*0.07
1198 !!$ qcoar_no3 = mult*0.03
1199 !!$ qfine_cl = mult*0.324
1200 !!$ qcoar_cl = mult*0.139
1201 !!$ qfine_msa = mult*0.0
1202 !!$ qcoar_msa = mult*0.0
1203 !!$ qfine_co3 = mult*0.0
1204 !!$ qcoar_co3 = mult*0.0
1205 !!$ qfine_nh4 = mult*0.98
1206 !!$ qcoar_nh4 = mult*0.10
1207 !!$ qfine_na = mult*0.21
1208 !!$ qcoar_na = mult*0.09
1209 !!$ qfine_ca = mult*0.0
1210 !!$ qcoar_ca = mult*0.0
1211 !!$ qfine_oin = mult*0.15
1212 !!$ qcoar_oin = mult*0.35
1213 !!$ qfine_oc = mult*1.00
1214 !!$ qcoar_oc = mult*1.50
1215 !!$ qfine_bc = mult*0.175
1216 !!$ qcoar_bc = mult*0.075
1217 !!$ qfine_hysw = mult*0.0
1218 !!$ qcoar_hysw = mult*0.0
1219 !!$ qfine_watr = mult*0.0
1220 !!$ qcoar_watr = mult*0.0
1221
1222
1223 ! qfine_so4 ... are ug/m3 so 1.0e-6 factor gives g/m3
1224 ! dens_so4 ... are g/cm3; qfine_vol ... are cm3/m3
1225 qfine_vol = 1.0e-6 * ( &
1226 (qfine_so4/dens_so4_aer) + (qfine_no3/dens_no3_aer) + &
1227 (qfine_cl /dens_cl_aer ) + (qfine_msa/dens_msa_aer) + &
1228 (qfine_co3/dens_co3_aer) + (qfine_nh4/dens_nh4_aer) + &
1229 (qfine_na /dens_na_aer ) + (qfine_ca /dens_ca_aer ) + &
1230 (qfine_oin/dens_oin_aer) + (qfine_oc /dens_oc_aer ) + &
1231 (qfine_bc /dens_bc_aer ) )
1232 qcoar_vol = 1.0e-6 * ( &
1233 (qcoar_so4/dens_so4_aer) + (qcoar_no3/dens_no3_aer) + &
1234 (qcoar_cl /dens_cl_aer ) + (qcoar_msa/dens_msa_aer) + &
1235 (qcoar_co3/dens_co3_aer) + (qcoar_nh4/dens_nh4_aer) + &
1236 (qcoar_na /dens_na_aer ) + (qcoar_ca /dens_ca_aer ) + &
1237 (qcoar_oin/dens_oin_aer) + (qcoar_oc /dens_oc_aer ) + &
1238 (qcoar_bc /dens_bc_aer ) )
1239
1240 qfine = -1.0e30
1241 qcoar = -1.0e30
1242
1243 ! identify the species by checking "nch" against the "lptr_xxx_a_amode(m)"
1244 do 2900 m = 1, nsize_aer(itype)
1245
1246 if (nch .eq. lptr_so4_aer(m,itype,iphase)) then
1247 qfine = qfine_so4
1248 qcoar = qcoar_so4
1249 else if (nch .eq. lptr_no3_aer(m,itype,iphase)) then
1250 qfine = qfine_no3
1251 qcoar = qcoar_no3
1252 else if (nch .eq. lptr_cl_aer(m,itype,iphase)) then
1253 qfine = qfine_cl
1254 qcoar = qcoar_cl
1255 else if (nch .eq. lptr_msa_aer(m,itype,iphase)) then
1256 qfine = qfine_msa
1257 qcoar = qcoar_msa
1258 else if (nch .eq. lptr_co3_aer(m,itype,iphase)) then
1259 qfine = qfine_co3
1260 qcoar = qcoar_co3
1261 else if (nch .eq. lptr_nh4_aer(m,itype,iphase)) then
1262 qfine = qfine_nh4
1263 qcoar = qcoar_nh4
1264 else if (nch .eq. lptr_na_aer(m,itype,iphase)) then
1265 qfine = qfine_na
1266 qcoar = qcoar_na
1267 else if (nch .eq. lptr_ca_aer(m,itype,iphase)) then
1268 qfine = qfine_ca
1269 qcoar = qcoar_ca
1270 else if (nch .eq. lptr_oin_aer(m,itype,iphase)) then
1271 qfine = qfine_oin
1272 qcoar = qcoar_oin
1273 else if (nch .eq. lptr_oc_aer(m,itype,iphase)) then
1274 qfine = qfine_oc
1275 qcoar = qcoar_oc
1276 else if (nch .eq. lptr_bc_aer(m,itype,iphase)) then
1277 qfine = qfine_bc
1278 qcoar = qcoar_bc
1279 else if (nch .eq. hyswptr_aer(m,itype)) then
1280 qfine = qfine_hysw
1281 qcoar = qcoar_hysw
1282 else if (nch .eq. waterptr_aer(m,itype)) then
1283 qfine = qfine_watr
1284 qcoar = qcoar_watr
1285 else if (nch .eq. numptr_aer(m,itype,iphase)) then
1286 dumvol1p = sqrt(volumlo_sect(m,itype)*volumhi_sect(m,itype))
1287 qfine = qfine_vol/dumvol1p
1288 qcoar = qcoar_vol/dumvol1p
1289 end if
1290
1291 if ((qfine >= 0.0) .and. (qcoar >= 0.0)) then
1292 qval = qfine*fr_fine(m) + qcoar*fr_coar(m)
1293 chem_bv = qval*alt
1294 goto 2910
1295 end if
1296
1297 2900 continue
1298 2910 continue
1299
1300 return
1301 end subroutine bdy_chem_value_mosaic
1302
1303