module_mosaic_movesect.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_movesect
10
11
12 use module_data_mosaic_asect
13 use module_data_mosaic_other
14 use module_peg_util
15
16
17 !
18 ! if apply_n1_inflow = 1, then subr move_sections_apply_n1_inflow
19 ! applies an ad_hoc correction to bin 1 to account for growth of
20 ! smaller particles (which are not simulated) growing into bin 1
21 ! if apply_n1_inflow = other, this correction is not applied
22 !
23 integer, parameter :: apply_n1_inflow = 1
24
25 contains
26
27
28
29 !-----------------------------------------------------------------------
30 !
31 ! zz08movesect.f - created 24-nov-01 from scm movebin code
32 ! 04-dec-01 rce - added code to treat bins that had state="no_aerosol"
33 ! 10-dec-01 rce - diagnostic writes to fort.96 deleted
34 ! 08-aug-02 rce - this is latest version from freshair scaqs-aerosol code
35 ! 03-aug-02 rce - bypass this routine when msectional=701
36 ! output nnewsave values to lunout when msectional=702
37 ! 29-aug-03 rce - use nspec_amode_nontracer in first "do ll" loop
38 ! 12-nov-03 rce - two approaches now available
39 ! jacobson moving-center (old) - when msectional=10
40 ! tzivion mass-number advection (new) - when msectional=20
41 ! 28-jan-04 rce - eliminated the move_sections_old code
42 ! (no reason to keep it)
43 !
44 ! rce 2004-dec-03 - many changes associated with the new aerosol "pointer"
45 ! variables in module_data_mosaic_asect
46 ! rce 2005-feb-17 - fixes involving drydens_pre/aftgrow, drymass_...,
47 ! & mw_aer. All are dimensioned (isize,itype) but were being used
48 ! as (itype). In old compiler, this was OK, and treated as (itype,1).
49 ! In new compiler, this caused an error.
50 ! Also, in subr test_move_sections, set iphase based on iflag.
51 !-----------------------------------------------------------------------
52
53
54 !-----------------------------------------------------------------------
55 subroutine move_sections( iflag, iclm, jclm, k, m)
56 !
57 ! routine transfers aerosol number and mass between sections
58 ! to account for continuous aerosol growth
59 ! this routine is called after the gas condensation module (MOSAIC) or
60 ! aqueous chemistry module has increased the mass within sections
61 !
62 ! moving-center algorithm or mass-number advection algorithm is used,
63 ! depending on value of mod(msectional,100)
64 ! section mean diameter is given by
65 ! vtot = ntot * (pi/6) * (dmean**3)
66 ! where vtot and ntot are total dry-volume and number for the section
67 ! if dmean is outside the section boundaries (dlo_sect & dhi_sect), then
68 ! all the mass and number in the section are transfered to the
69 ! section with dlo_sect(nnew) < dmean < dhi_sect(nnew)
70 !
71 ! mass mixing ratios (mole/mole-air or g/mole-air) are in
72 ! rsub(massptr_aer(ll,n,itype,iphase),k,m)
73 ! number mixing ratios (particles/mole-air) are in
74 ! rsub(numptr_aer(n,itype,iphase),k,m)
75 ! these values are over-written with new values
76 ! the following are also updated: adrydens_sub(n,k,m), admeandry_sub(n,k,m),
77 ! awetdens_sub(n,k,m), admeanwet_sub(n,k,m)
78 !
79 ! particle sizes are in cm
80 !
81 ! input parameters
82 ! iflag = 1 - do transfer after trace-gas condensation/evaporation
83 ! = 2 - do transfer after aqueous chemistry
84 ! = -1/-2 - do some "first entry" tasks for the iflag=+1/+2 cases
85 ! iclm, jclm, k = current i,j,k indices
86 ! m = current subarea index
87 ! drymass_pregrow(n,itype) = dry-mass (g/mole-air) for section n
88 ! before the growth
89 ! drymass_aftgrow(n,itype) = dry-mass (g/mole-air) for section n
90 ! after the growth but before inter-section transfer
91 ! drydens_pregrow(n,itype) = dry-density (g/cm3) for section n
92 ! before the growth
93 ! drydens_aftgrow(n,itype) = dry-density (g/cm3) for section n
94 ! after the growth but before inter-section transfer
95 ! (drymass_pregrow and drydens_pregrow are not used by the moving-center
96 ! algorithm, but would be needed for other sectional algorithms)
97 !
98 ! this routine is the top level module for the post-october-2003 code
99 ! and will do either moving-center or mass-number advection
100 !
101 implicit none
102
103 ! include 'v33com'
104 ! include 'v33com2'
105 ! include 'v33com3'
106 ! include 'v33com9a'
107 ! include 'v33com9b'
108
109 ! subr arguments
110 integer iflag, iclm, jclm, k, m
111
112 ! local variables
113 integer idiag_movesect, iphase, itype, &
114 l, ll, llhysw, llwater, lnew, lold, l3, &
115 method_movesect, n, nnew, nold
116 integer nnewsave(2,maxd_asize)
117
118 real densdefault, densh2o, smallmassaa, smallmassbb
119 real delta_water_conform1, delta_numb_conform1
120
121 real drydenspp(maxd_asize), drydensxx0(maxd_asize), &
122 drydensxx(maxd_asize), drydensyy(maxd_asize)
123 real drymasspp(maxd_asize), drymassxx0(maxd_asize), &
124 drymassxx(maxd_asize), drymassyy(maxd_asize)
125 real dryvolxx(maxd_asize), dryvolyy(maxd_asize)
126 real rmassxx(maxd_acomp+2,maxd_asize), &
127 rmassyy(maxd_acomp+2,maxd_asize)
128 real rnumbpp(maxd_asize), rnumbxx0(maxd_asize), &
129 rnumbxx(maxd_asize), rnumbyy(maxd_asize)
130 real specdensxx(maxd_acomp), specmwxx(maxd_acomp)
131 real xferfracvol(2,maxd_asize), xferfracnum(2,maxd_asize)
132 real wetvolxx(maxd_asize), wetvolyy(maxd_asize)
133 real wetmassxx(maxd_asize), wetmassyy(maxd_asize)
134
135 character*160 msg
136
137
138 !
139 ! check for valid inputs
140 !
141 if (msectional .le. 0) return
142 if (ntype_aer .le. 0) return
143 if (nphase_aer .le. 0) return
144
145 !
146 ! run diagnostic tests
147 ! (these will only be run for certain values of idiag_movesect
148 ! rce 2004-nov-29 - to avoid recursion, test_move_sections
149 ! is now called from mosaic_driver
150 !
151 ! call test_move_sections( iflag, iclm, jclm, k, m )
152
153
154 ! get "method_movesect" from digits 1-2 of msectional (treat 1-9 as 10)
155 method_movesect = mod( msectional, 100 )
156 if (method_movesect .le. 10) method_movesect = 10
157
158 ! get "idiag_movesect" from digits 3-4 of msectional
159 idiag_movesect = mod( msectional, 10000 )/100
160
161 if ((method_movesect .eq. 10) .or. &
162 (method_movesect .eq. 11) .or. &
163 (method_movesect .eq. 20) .or. &
164 (method_movesect .eq. 21) .or. &
165 (method_movesect .eq. 30) .or. &
166 (method_movesect .eq. 31)) then
167 continue
168 else if ((method_movesect .eq. 19) .or. &
169 (method_movesect .eq. 29) .or. &
170 (method_movesect .eq. 39)) then
171 return
172 else
173 msg = '*** subr move_sections error - ' // &
174 'illegal value for msectional'
175 call peg_error_fatal( lunerr, msg )
176 end if
177
178 if (iabs(iflag) .eq. 1) then
179 iphase = ai_phase
180 else if (iabs(iflag) .eq. 2) then
181 iphase = cw_phase
182 msg = '*** subr move_sections error - ' // &
183 'iflag=2 (after aqueous chemistry) not implemented'
184 call peg_error_fatal( lunerr, msg )
185 else
186 msg = '*** subr move_sections error - ' // &
187 'iabs(iflag) must be 1 or 2'
188 call peg_error_fatal( lunerr, msg )
189 end if
190
191
192 ! when iflag=-1/-2, call move_sections_checkptrs then return
193 ! if ((ncorecnt .le. 0) .and. (k .le. 1)) then
194 if (iflag .le. 0) then
195 write(msg,9040) 'method', method_movesect
196 call peg_message( lunout, msg )
197 write(msg,9040) 'idiag ', idiag_movesect
198 call peg_message( lunout, msg )
199 call move_sections_checkptrs( iflag, iclm, jclm, k, m )
200 return
201 end if
202 9040 format( '*** subr move_sections - ', a, ' =', i6 )
203
204
205 densdefault = 2.0
206 densh2o = 1.0
207
208 ! if bin mass mixrat < smallmassaa (1.e-20 g/mol), then assume no growth
209 ! AND no water AND conform number so that size is within bin limits
210 smallmassaa = 1.0e-20
211 ! if bin mass OR volume mixrat < smallmassbb (1.e-30 g/mol), then
212 ! assume default density to avoid divide by zero
213 smallmassbb = 1.0e-30
214
215
216 ! process each type, one at a time
217 do 1900 itype = 1, ntype_aer
218
219 if (nsize_aer(itype) .le. 0) goto 1900
220
221 call move_sections_initial_conform( &
222 iflag, iclm, jclm, k, m, iphase, itype, &
223 method_movesect, idiag_movesect, llhysw, llwater, &
224 densdefault, densh2o, smallmassaa, smallmassbb, &
225 delta_water_conform1, delta_numb_conform1, &
226 drydenspp, drymasspp, rnumbpp, &
227 drydensxx0, drymassxx0, rnumbxx0, &
228 drydensxx, drymassxx, dryvolxx, rmassxx, rnumbxx, &
229 wetmassxx, wetvolxx, &
230 specmwxx, specdensxx )
231
232 if (method_movesect .le. 19) then
233 call move_sections_calc_movingcenter( &
234 iflag, iclm, jclm, k, m, iphase, itype, &
235 method_movesect, idiag_movesect, llhysw, llwater, nnewsave, &
236 densdefault, densh2o, smallmassaa, smallmassbb, &
237 drydensxx, drymassxx, dryvolxx, rmassxx, rnumbxx, &
238 wetmassxx, wetvolxx, &
239 xferfracvol, xferfracnum )
240 else
241 call move_sections_calc_masnumadvect( &
242 iflag, iclm, jclm, k, m, iphase, itype, &
243 method_movesect, idiag_movesect, llhysw, llwater, nnewsave, &
244 densdefault, densh2o, smallmassaa, smallmassbb, &
245 drydenspp, drymasspp, rnumbpp, &
246 drydensxx0, drymassxx0, rnumbxx0, &
247 drydensxx, drymassxx, dryvolxx, rmassxx, rnumbxx, &
248 wetmassxx, wetvolxx, &
249 xferfracvol, xferfracnum )
250 end if
251
252 call move_sections_apply_moves( &
253 iflag, iclm, jclm, k, m, iphase, itype, &
254 method_movesect, idiag_movesect, llhysw, llwater, nnewsave, &
255 densdefault, densh2o, smallmassaa, smallmassbb, &
256 delta_water_conform1, delta_numb_conform1, &
257 drymassxx, dryvolxx, rmassxx, rnumbxx, wetmassxx, wetvolxx, &
258 drymassyy, dryvolyy, rmassyy, rnumbyy, wetmassyy, wetvolyy, &
259 xferfracvol, xferfracnum )
260
261 ! rce 08-nov-2004 - this call is new (and perhaps temporary)
262 if (apply_n1_inflow .eq. 1) then
263 call move_sections_apply_n1_inflow( &
264 iflag, iclm, jclm, k, m, iphase, itype, &
265 method_movesect, idiag_movesect, llhysw, llwater, nnewsave, &
266 densdefault, densh2o, smallmassaa, smallmassbb, &
267 delta_water_conform1, delta_numb_conform1, &
268 drymassxx, dryvolxx, rmassxx, rnumbxx, wetmassxx, wetvolxx, &
269 drymassyy, dryvolyy, rmassyy, rnumbyy, wetmassyy, wetvolyy, &
270 xferfracvol, xferfracnum, &
271 specmwxx, specdensxx )
272 end if
273
274 ! call move_sections_final_conform( &
275 ! iflag, iclm, jclm, k, m, iphase, itype )
276
277 1900 continue
278
279 return
280 end subroutine move_sections
281
282
283 !-----------------------------------------------------------------------
284 subroutine move_sections_initial_conform( &
285 iflag, iclm, jclm, k, m, iphase, itype, &
286 method_movesect, idiag_movesect, llhysw, llwater, &
287 densdefault, densh2o, smallmassaa, smallmassbb, &
288 delta_water_conform1, delta_numb_conform1, &
289 drydenspp, drymasspp, rnumbpp, &
290 drydensxx0, drymassxx0, rnumbxx0, &
291 drydensxx, drymassxx, dryvolxx, rmassxx, rnumbxx, &
292 wetmassxx, wetvolxx, &
293 specmwxx, specdensxx )
294
295 !
296 ! routine does some initial tasks for the section movement
297 ! load rmassxx & rnumbxx from rsub
298 ! load specmwxx & specdensxx
299 ! set drymassxx & dryvolxx from drymass_aftgrow & drydens_aftgrow,
300 ! OR compute them from rmassxx, specmwxx, specdensxx if need be
301 ! set wetmassxx & wetvolxx from dry values & water mass
302 ! conform rnumbxx so that the mean particle size of each section
303 ! (= dryvolxx/rnumbxx) is within the section limits
304 !
305 implicit none
306
307 ! include 'v33com'
308 ! include 'v33com2'
309 ! include 'v33com3'
310 ! include 'v33com9a'
311 ! include 'v33com9b'
312
313 ! subr arguments
314 integer iflag, iclm, jclm, iphase, itype, k, &
315 m, method_movesect, idiag_movesect, llhysw, llwater
316 real densdefault, densh2o, smallmassaa, smallmassbb
317 real delta_water_conform1, delta_numb_conform1
318 real drydenspp(maxd_asize), drydensxx0(maxd_asize), &
319 drydensxx(maxd_asize)
320 real drymasspp(maxd_asize), drymassxx0(maxd_asize), &
321 drymassxx(maxd_asize)
322 real dryvolxx(maxd_asize)
323 real rmassxx(maxd_acomp+2,maxd_asize)
324 real rnumbpp(maxd_asize), rnumbxx0(maxd_asize), &
325 rnumbxx(maxd_asize)
326 real specdensxx(maxd_acomp), specmwxx(maxd_acomp)
327 real wetvolxx(maxd_asize)
328 real wetmassxx(maxd_asize)
329
330
331 ! local variables
332 integer l, ll, lnew, lold, l3, n, nnew, nold
333
334 real dummass, dumnum, dumnum_at_dhi, dumnum_at_dlo, dumr, &
335 dumvol, dumvol1p, dumwatrmass
336
337
338 ! assure positive definite
339 do l = 1, ltot2
340 rsub(l,k,m) = max( 0., rsub(l,k,m) )
341 end do
342
343 ! load mixrats into working arrays and assure positive definite
344 llhysw = ncomp_plustracer_aer(itype) + 1
345 llwater = ncomp_plustracer_aer(itype) + 2
346 do n = 1, nsize_aer(itype)
347 do ll = 1, ncomp_plustracer_aer(itype)
348 l = massptr_aer(ll,n,itype,iphase)
349 rmassxx(ll,n) = rsub(l,k,m)
350 end do
351 rmassxx(llhysw,n) = 0.
352 l = 0
353 if (iphase .eq. ai_phase) l = hyswptr_aer(n,itype)
354 if (l .gt. 0) rmassxx(llhysw,n) = rsub(l,k,m)
355 rmassxx(llwater,n) = 0.
356 l = 0
357 if (iphase .eq. ai_phase) l = waterptr_aer(n,itype)
358 if (l .gt. 0) rmassxx(llwater,n) = rsub(l,k,m)
359
360 rnumbxx(n) = rsub(numptr_aer(n,itype,iphase),k,m)
361 rnumbxx0(n) = rnumbxx(n)
362 rnumbpp(n) = rnumbxx(n)
363
364 drydenspp(n) = drydens_pregrow(n,itype)
365 drymasspp(n) = drymass_pregrow(n,itype)
366
367 drydensxx(n) = drydens_aftgrow(n,itype)
368 drymassxx(n) = drymass_aftgrow(n,itype)
369 drydensxx0(n) = drydensxx(n)
370 drymassxx0(n) = drymassxx(n)
371 end do
372
373 ! load specmw and specdens also
374 do ll = 1, ncomp_plustracer_aer(itype)
375 specmwxx(ll) = mw_aer(ll,itype)
376 specdensxx(ll) = dens_aer(ll,itype)
377 end do
378
379 delta_water_conform1 = 0.0
380 delta_numb_conform1 = 0.0
381
382
383 do 1390 n = 1, nsize_aer(itype)
384
385 !
386 ! if drydens_aftgrow < 0.5, then bin had state="no_aerosol"
387 ! compute volume using default dry-densities, set water=0,
388 ! and conform the number
389 ! also do this if mass is extremely small (below smallmassaa)
390 !
391 if ( (drydensxx(n) .le. 0.5) .or. &
392 (drymassxx(n) .le. smallmassaa) ) then
393 dummass = 0.
394 dumvol = 0.
395 do ll = 1, ncomp_aer(itype)
396 dumr = rmassxx(ll,n)*specmwxx(ll)
397 dummass = dummass + dumr
398 dumvol = dumvol + dumr/specdensxx(ll)
399 end do
400 drymassxx(n) = dummass
401 if (min(dummass,dumvol) .le. smallmassbb) then
402 drydensxx(n) = densdefault
403 dumvol = dummass/densdefault
404 dumnum = dumvol/volumcen_sect(n,itype)
405 else
406 drydensxx(n) = dummass/dumvol
407 dumnum = rnumbxx(n)
408 dumnum_at_dhi = dumvol/volumhi_sect(n,itype)
409 dumnum_at_dlo = dumvol/volumlo_sect(n,itype)
410 dumnum = max( dumnum_at_dhi, min( dumnum_at_dlo, dumnum ) )
411 end if
412 delta_numb_conform1 = delta_numb_conform1 + dumnum - rnumbxx(n)
413 rnumbxx(n) = dumnum
414 rnumbpp(n) = rnumbxx(n)
415 delta_water_conform1 = delta_water_conform1 - rmassxx(llwater,n)
416 rmassxx(llwater,n) = 0.
417 end if
418
419 ! load dry/wet mass and volume into "xx" arrays
420 ! which hold values before inter-mode transferring
421 dryvolxx(n) = drymassxx(n)/drydensxx(n)
422 dumwatrmass = rmassxx(llwater,n)*mw_water_aer
423 wetmassxx(n) = drymassxx(n) + dumwatrmass
424 wetvolxx(n) = dryvolxx(n) + dumwatrmass/densh2o
425
426 1390 continue
427
428 return
429 end subroutine move_sections_initial_conform
430
431
432 !-----------------------------------------------------------------------
433 subroutine move_sections_calc_movingcenter( &
434 iflag, iclm, jclm, k, m, iphase, itype, &
435 method_movesect, idiag_movesect, llhysw, llwater, nnewsave, &
436 densdefault, densh2o, smallmassaa, smallmassbb, &
437 drydensxx, drymassxx, dryvolxx, rmassxx, rnumbxx, &
438 wetmassxx, wetvolxx, &
439 xferfracvol, xferfracnum )
440 !
441 ! routine calculates section movements for the moving-center approach
442 !
443 ! material in section n will be transfered to section nnewsave(1,n)
444 !
445 ! the nnewsave are calculated here
446 ! the actual transfer is done in another routine
447 !
448 implicit none
449
450 ! include 'v33com'
451 ! include 'v33com2'
452 ! include 'v33com3'
453 ! include 'v33com9a'
454 ! include 'v33com9b'
455
456 ! subr arguments
457 integer iflag, iclm, jclm, iphase, itype, k, &
458 m, method_movesect, idiag_movesect, llhysw, llwater
459 integer nnewsave(2,maxd_asize)
460 real densdefault, densh2o, smallmassaa, smallmassbb
461 real drydensxx(maxd_asize)
462 real drymassxx(maxd_asize)
463 real dryvolxx(maxd_asize)
464 real rmassxx(maxd_acomp+2,maxd_asize)
465 real rnumbxx(maxd_asize)
466 real xferfracvol(2,maxd_asize), xferfracnum(2,maxd_asize)
467 real wetmassxx(maxd_asize)
468 real wetvolxx(maxd_asize)
469
470 ! local variables
471 integer ll, n, ndum, nnew, nold
472 real dumnum, dumvol, dumvol1p, sixoverpi, third
473 character*160 msg
474
475
476 sixoverpi = 6.0/pi
477 third = 1.0/3.0
478
479 !
480 ! compute mean size after growth (and corresponding section)
481 ! particles in section n will be transferred to section nnewsave(1,n)
482 !
483 do 1390 n = 1, nsize_aer(itype)
484
485 nnew = n
486
487 ! don't bother to transfer bins whose mass is extremely small
488 if (drymassxx(n) .le. smallmassaa) goto 1290
489
490 dumvol = dryvolxx(n)
491 dumnum = rnumbxx(n)
492
493 ! check for number so small that particle volume is
494 ! above that of largest section
495 if (dumnum .le. dumvol/volumhi_sect(nsize_aer(itype),itype)) then
496 nnew = nsize_aer(itype)
497 goto 1290
498 ! or below that of smallest section
499 else if (dumnum .ge. dumvol/volumlo_sect(1,itype)) then
500 nnew = 1
501 goto 1290
502 end if
503
504 ! dumvol1p is mean particle volume (cm3) for the section
505 dumvol1p = dumvol/dumnum
506 if (dumvol1p .gt. volumhi_sect(n,itype)) then
507 do while ( (nnew .lt. nsize_aer(itype)) .and. &
508 (dumvol1p .gt. volumhi_sect(nnew,itype)) )
509 nnew = nnew + 1
510 end do
511
512 else if (dumvol1p .lt. volumlo_sect(n,itype)) then
513 do while ( (nnew .gt. 1) .and. &
514 (dumvol1p .lt. volumlo_sect(nnew,itype)) )
515 nnew = nnew - 1
516 end do
517
518 end if
519
520 1290 nnewsave(1,n) = nnew
521 nnewsave(2,n) = 0
522
523 xferfracvol(1,n) = 1.0
524 xferfracvol(2,n) = 0.0
525 xferfracnum(1,n) = 1.0
526 xferfracnum(2,n) = 0.0
527
528 1390 continue
529
530
531 ! diagnostic output
532 ! output nnewsave values when msectional=7xxx
533 if (idiag_movesect .eq. 70) then
534 ndum = 0
535 do n = 1, nsize_aer(itype)
536 if (nnewsave(1,n) .ne. n) ndum = ndum + 1
537 end do
538 if (ndum .gt. 0) then
539 write(msg,97751) 'YES', iclm, jclm, k, m, &
540 ndum, (nnewsave(1,n), n=1,nsize_aer(itype))
541 call peg_message( lunout, msg )
542 else
543 write(msg,97751) 'NO ', iclm, jclm, k, m, &
544 ndum, (nnewsave(1,n), n=1,nsize_aer(itype))
545 call peg_message( lunout, msg )
546 end if
547 end if
548 97751 format( 'movesect', a, 4i3, 3x, i3, 3x, 14i3 )
549
550 return
551 end subroutine move_sections_calc_movingcenter
552
553
554 !-----------------------------------------------------------------------
555 subroutine move_sections_calc_masnumadvect( &
556 iflag, iclm, jclm, k, m, iphase, itype, &
557 method_movesect, idiag_movesect, llhysw, llwater, nnewsave, &
558 densdefault, densh2o, smallmassaa, smallmassbb, &
559 drydenspp, drymasspp, rnumbpp, &
560 drydensxx0, drymassxx0, rnumbxx0, &
561 drydensxx, drymassxx, dryvolxx, rmassxx, rnumbxx, &
562 wetmassxx, wetvolxx, &
563 xferfracvol, xferfracnum )
564 !
565 ! routine calculates section movements for the mass-number-advection approach
566 !
567 ! material in section n will be transfered to sections
568 ! nnewsave(1,n) and nnewsave(2,n)
569 ! the fractions of mass/volume transfered to each are
570 ! xferfracvol(1,n) and xferfracvol(2,n)
571 ! the fractions of number transfered to each are
572 ! xferfracnum(1,n) and xferfracnum(2,n)
573 !
574 ! the nnewsave, xferfracvol, and xferfracnum are calculated here
575 ! the actual transfer is done in another routine
576 !
577 implicit none
578
579 ! include 'v33com'
580 ! include 'v33com2'
581 ! include 'v33com3'
582 ! include 'v33com9a'
583 ! include 'v33com9b'
584
585 ! subr arguments
586 integer iflag, iclm, jclm, iphase, itype, k, &
587 m, method_movesect, idiag_movesect, llhysw, llwater
588 integer nnewsave(2,maxd_asize)
589
590 real densdefault, densh2o, smallmassaa, smallmassbb
591 real drydenspp(maxd_asize), drydensxx0(maxd_asize), &
592 drydensxx(maxd_asize)
593 real drymasspp(maxd_asize), drymassxx0(maxd_asize), &
594 drymassxx(maxd_asize)
595 real dryvolxx(maxd_asize)
596 real rmassxx(maxd_acomp+2,maxd_asize)
597 real rnumbpp(maxd_asize), rnumbxx0(maxd_asize), &
598 rnumbxx(maxd_asize)
599 real xferfracvol(2,maxd_asize), xferfracnum(2,maxd_asize)
600 real wetvolxx(maxd_asize)
601 real wetmassxx(maxd_asize)
602
603 ! local variables
604 integer ierr, n, nnew, nnew2
605 integer iforce_movecenter(maxd_asize)
606
607 real dum1, dum2, dum3
608 real dumaa, dumbb, dumgamma, dumratio
609 real dumfracnum, dumfracvol
610 real dumntot
611 real dumv
612 real dumvbar_aft, dumvbar_pre
613 real dumvcutlo_nnew_pre, dumvcuthi_nnew_pre
614 real dumvlo_pre, dumvhi_pre, dumvdel_pre
615 real dumvtot_aft, dumvtot_pre
616 real dumzlo, dumzhi
617 real sixoverpi, third
618
619 character*4 dumch4
620 character*1 dumch1
621 character*160 msg
622
623
624 sixoverpi = 6.0/pi
625 third = 1.0/3.0
626
627 !
628 ! compute mean size after growth (and corresponding section)
629 ! some of the particles in section n will be transferred to section nnewsave(1,n)
630 !
631 ! if the aftgrow mass is extremely small,
632 ! OR if the aftgrow mean size is outside of
633 ! [dlo_sect(1,itype), dhi_sect(nsize_aer(itype),itype)]
634 ! then use the moving-center method_movesect for this bin
635 ! (don't try to compute the pregrow within-bin distribution)
636 !
637 do 3900 n = 1, nsize_aer(itype)
638
639 nnew = n
640 iforce_movecenter(n) = 0
641
642 xferfracvol(1,n) = 1.0
643 xferfracvol(2,n) = 0.0
644 xferfracnum(1,n) = 1.0
645 xferfracnum(2,n) = 0.0
646
647 dumvtot_aft = -1.0
648 dumvtot_pre = -1.0
649 dumvbar_aft = -1.0
650 dumvbar_pre = -1.0
651 dumvlo_pre = -1.0
652 dumvhi_pre = -1.0
653 dumgamma = -1.0
654 dumratio = -1.0
655 dumvcutlo_nnew_pre = volumlo_sect(nnew,itype)*(dumvbar_pre/dumvbar_aft)
656 dumvcuthi_nnew_pre = volumhi_sect(nnew,itype)*(dumvbar_pre/dumvbar_aft)
657 dumfracvol = -1.0
658 dumfracnum = -1.0
659
660 ! don't bother to transfer bins whose mass is extremely small
661 if (drymassxx(n) .le. smallmassaa) then
662 iforce_movecenter(n) = 1
663 goto 1290
664 end if
665
666 dumvtot_aft = dryvolxx(n)
667 dumntot = rnumbxx(n)
668
669 ! check for particle volume above that of largest section
670 ! or below that of smallest section
671 if (dumntot .le. dumvtot_aft/volumhi_sect(nsize_aer(itype),itype)) then
672 nnew = nsize_aer(itype)
673 iforce_movecenter(n) = 2
674 goto 1290
675 else if (dumntot .ge. dumvtot_aft/volumlo_sect(1,itype)) then
676 nnew = 1
677 iforce_movecenter(n) = 3
678 goto 1290
679 end if
680
681 ! dumvbar_aft is mean particle volume (cm3) for the section
682 ! find the section that encloses this volume
683 dumvbar_aft = dumvtot_aft/dumntot
684 if (dumvbar_aft .gt. volumhi_sect(n,itype)) then
685 do while ( (nnew .lt. nsize_aer(itype)) .and. &
686 (dumvbar_aft .gt. volumhi_sect(nnew,itype)) )
687 nnew = nnew + 1
688 end do
689
690 else if (dumvbar_aft .lt. volumlo_sect(n,itype)) then
691 do while ( (nnew .gt. 1) .and. &
692 (dumvbar_aft .lt. volumlo_sect(nnew,itype)) )
693 nnew = nnew - 1
694 end do
695
696 end if
697
698 1290 nnewsave(1,n) = nnew
699 nnewsave(2,n) = 0
700
701 if (iforce_movecenter(n) .gt. 0) goto 3700
702
703
704 ! if drydenspp (pregrow) < 0.5 (because bin had state="no_aerosol" before
705 ! growth was computed, so its initial mass was very small)
706 ! then use the moving-center method_movesect for this bin
707 ! (don't try to compute the pregrow within-bin distribution)
708 ! also do this if pregrow mass is extremely small (below smallmassaa)
709 if ( (drydenspp(n) .le. 0.5) .or. &
710 (drymasspp(n) .le. smallmassaa) ) then
711 iforce_movecenter(n) = 11
712 goto 3700
713 end if
714
715 dumvtot_pre = drymasspp(n)/drydenspp(n)
716
717 dumvlo_pre = volumlo_sect(n,itype)
718 dumvhi_pre = volumhi_sect(n,itype)
719 dumvdel_pre = dumvhi_pre - dumvlo_pre
720
721 ! if the pregrow mean size is outside of OR very close to the bin limits,
722 ! then use moving-center approach for this bin
723 dumv = dumvhi_pre - 0.01*dumvdel_pre
724 if (dumntot .le. dumvtot_pre/dumv) then
725 iforce_movecenter(n) = 12
726 goto 3700
727 end if
728 dumv = dumvlo_pre + 0.01*dumvdel_pre
729 if (dumntot .ge. dumvtot_pre/dumv) then
730 iforce_movecenter(n) = 13
731 goto 3700
732 end if
733
734 ! calculate the pregrow within-section size distribution
735 dumvbar_pre = dumvtot_pre/dumntot
736 dumgamma = (dumvhi_pre/dumvlo_pre) - 1.0
737 dumratio = dumvbar_pre/dumvlo_pre
738
739 if (dumratio .le. (1.0001 + dumgamma/3.0)) then
740 dumv = dumvlo_pre + 3.0*(dumvbar_pre-dumvlo_pre)
741 dumvhi_pre = min( dumvhi_pre, dumv )
742 dumvdel_pre = dumvhi_pre - dumvlo_pre
743 dumgamma = (dumvhi_pre/dumvlo_pre) - 1.0
744 dumratio = dumvbar_pre/dumvlo_pre
745 else if (dumratio .ge. (0.9999 + dumgamma*2.0/3.0)) then
746 dumv = dumvhi_pre + 3.0*(dumvbar_pre-dumvhi_pre)
747 dumvlo_pre = max( dumvlo_pre, dumv )
748 dumvdel_pre = dumvhi_pre - dumvlo_pre
749 dumgamma = (dumvhi_pre/dumvlo_pre) - 1.0
750 dumratio = dumvbar_pre/dumvlo_pre
751 end if
752
753 dumbb = (dumratio - 1.0 - 0.5*dumgamma)*12.0/dumgamma
754 dumaa = 1.0 - 0.5*dumbb
755
756 ! calculate pregrow volumes corresponding to the nnew
757 ! section boundaries
758 dumvcutlo_nnew_pre = volumlo_sect(nnew,itype)*(dumvbar_pre/dumvbar_aft)
759 dumvcuthi_nnew_pre = volumhi_sect(nnew,itype)*(dumvbar_pre/dumvbar_aft)
760
761 ! if the [dumvlo_pre, dumvhi_pre] falls completely within
762 ! the [dumvcutlo_nnew_pre, dumvcuthi_nnew_pre] interval,
763 ! then all mass and number go to nnew
764 if (nnew .eq. 1) then
765 if (dumvhi_pre .le. dumvcuthi_nnew_pre) then
766 iforce_movecenter(n) = 21
767 else
768 nnew2 = nnew + 1
769 end if
770 else if (nnew .eq. nsize_aer(itype)) then
771 if (dumvlo_pre .ge. dumvcutlo_nnew_pre) then
772 iforce_movecenter(n) = 22
773 else
774 nnew2 = nnew - 1
775 end if
776 else
777 if ((dumvlo_pre .ge. dumvcutlo_nnew_pre) .and. &
778 (dumvhi_pre .le. dumvcuthi_nnew_pre)) then
779 iforce_movecenter(n) = 23
780 else if (dumvlo_pre .lt. dumvcutlo_nnew_pre) then
781 nnew2 = nnew - 1
782 else
783 nnew2 = nnew + 1
784 end if
785 end if
786 if (iforce_movecenter(n) .gt. 0) goto 3700
787
788 ! calculate the fraction of ntot and vtot that are within
789 ! the [dumvcutlo_nnew_pre, dumvcuthi_nnew_pre] interval
790 dumzlo = (dumvcutlo_nnew_pre - dumvlo_pre)/dumvdel_pre
791 dumzhi = (dumvcuthi_nnew_pre - dumvlo_pre)/dumvdel_pre
792 dumzlo = max( dumzlo, 0.0 )
793 dumzhi = min( dumzhi, 1.0 )
794 dum1 = dumzhi - dumzlo
795 dum2 = (dumzhi**2 - dumzlo**2)*0.5
796 dum3 = (dumzhi**3 - dumzlo**3)/3.0
797 dumfracnum = dumaa*dum1 + dumbb*dum2
798 dumfracvol = (dumvlo_pre/dumvbar_pre) * (dumaa*dum1 + &
799 (dumaa*dumgamma + dumbb)*dum2 + (dumbb*dumgamma)*dum3)
800
801 if ((dumfracnum .le. 0.0) .or. (dumfracvol .le. 0.0)) then
802 iforce_movecenter(n) = 31
803 nnewsave(1,n) = nnew2
804 else if ((dumfracnum .ge. 1.0) .or. (dumfracvol .ge. 1.0)) then
805 iforce_movecenter(n) = 32
806 end if
807 if (iforce_movecenter(n) .gt. 0) goto 3700
808
809 nnewsave(2,n) = nnew2
810
811 xferfracvol(1,n) = dumfracvol
812 xferfracvol(2,n) = 1.0 - dumfracvol
813 xferfracnum(1,n) = dumfracnum
814 xferfracnum(2,n) = 1.0 - dumfracnum
815
816 3700 continue
817
818 ! output nnewsave values when msectional=7xxx
819 if (idiag_movesect .ne. 70) goto 3800
820
821 if (nnewsave(2,n) .eq. 0) then
822 if (nnewsave(1,n) .eq. 0) then
823 dumch4 = 'NO X'
824 else if (nnewsave(1,n) .eq. n) then
825 dumch4 = 'NO A'
826 else
827 dumch4 = 'YESA'
828 end if
829 else if (nnewsave(1,n) .eq. 0) then
830 if (nnewsave(2,n) .eq. n) then
831 dumch4 = 'NO B'
832 else
833 dumch4 = 'YESB'
834 end if
835 else if (nnewsave(2,n) .eq. n) then
836 if (nnewsave(1,n) .eq. n) then
837 dumch4 = 'NO Y'
838 else
839 dumch4 = 'YESC'
840 end if
841 else if (nnewsave(1,n) .eq. n) then
842 dumch4 = 'YESD'
843 else
844 dumch4 = 'YESE'
845 end if
846
847 dumch1 = '+'
848 if (drymasspp(n) .gt. drymassxx(n)) dumch1 = '-'
849
850 msg = ' '
851 call peg_message( lunout, msg )
852 write(msg,97010) dumch1, dumch4, iclm, jclm, k, m, &
853 n, nnewsave(1,n), nnewsave(2,n), iforce_movecenter(n)
854 call peg_message( lunout, msg )
855 write(msg,97020) 'pre mass, dens ', &
856 drymasspp(n), drydenspp(n)
857 call peg_message( lunout, msg )
858 write(msg,97020) 'aft mass, dens, numb', &
859 drymassxx(n), drydensxx(n), rnumbxx(n)
860 call peg_message( lunout, msg )
861 if ((drydensxx(n) .ne. drydensxx0(n)) .or. &
862 (drymassxx(n) .ne. drymassxx0(n)) .or. &
863 (rnumbxx(n) .ne. rnumbxx0(n) )) then
864 write(msg,97020) 'aft0 mas, dens, numb', &
865 drymassxx0(n), drydensxx0(n), rnumbxx0(n)
866 call peg_message( lunout, msg )
867 end if
868 write(msg,97020) 'vlop0, vbarp, vhip0', &
869 volumlo_sect(n,itype), dumvbar_pre, volumhi_sect(n,itype)
870 call peg_message( lunout, msg )
871 write(msg,97020) 'vlop , vbarp, vhip ', &
872 dumvlo_pre, dumvbar_pre, dumvhi_pre
873 call peg_message( lunout, msg )
874 write(msg,97020) 'vloax, vbarax, vhiax', &
875 dumvcutlo_nnew_pre, dumvbar_pre, dumvcuthi_nnew_pre
876 call peg_message( lunout, msg )
877 write(msg,97020) 'vloa0, vbara, vhia0', &
878 volumlo_sect(nnew,itype), dumvbar_aft, volumhi_sect(nnew,itype)
879 call peg_message( lunout, msg )
880 write(msg,97020) 'dumfrvol, num, ratio', &
881 dumfracvol, dumfracnum, dumratio
882 call peg_message( lunout, msg )
883 write(msg,97020) 'frvol,num1; vol,num2', &
884 xferfracvol(1,n), xferfracnum(1,n), &
885 xferfracvol(2,n), xferfracnum(2,n)
886 call peg_message( lunout, msg )
887
888 97010 format( 'movesect', 2a, 7x, 4i3, 4x, &
889 'n,nnews', 3i3, 4x, 'iforce', i3.2 )
890 97020 format( a, 1p, 4e13.4 )
891
892 3800 continue
893
894 !
895 ! check for legal combinations of nnewsave(1,n) & nnewsave(2,n)
896 ! error if
897 ! nnew1 == nnew2
898 ! both are non-zero AND iabs(nnew1-nnew2) != 1
899 ierr = 0
900 if (nnewsave(1,n) .eq. nnewsave(2,n)) then
901 ierr = 1
902 else if (nnewsave(1,n)*nnewsave(2,n) .ne. 0) then
903 if (iabs(nnewsave(1,n)-nnewsave(2,n)) .ne. 1) ierr = 1
904 end if
905 if (ierr .gt. 0) then
906 write(msg,97010) 'E', 'RROR', iclm, jclm, k, m, &
907 n, nnewsave(1,n), nnewsave(2,n), iforce_movecenter(n)
908 call peg_message( lunout, msg )
909 end if
910
911
912 ! if method_movesect == 30-31 then force moving center
913 ! this is just for testing purposes
914 if ((method_movesect .ge. 30) .and. (method_movesect .le. 39)) then
915 nnewsave(1,n) = nnew
916 nnewsave(2,n) = 0
917 xferfracvol(1,n) = 1.0
918 xferfracvol(2,n) = 0.0
919 xferfracnum(1,n) = 1.0
920 xferfracnum(2,n) = 0.0
921 end if
922
923 3900 continue
924
925 return
926 end subroutine move_sections_calc_masnumadvect
927
928
929 !-----------------------------------------------------------------------
930 subroutine move_sections_apply_moves( &
931 iflag, iclm, jclm, k, m, iphase, itype, &
932 method_movesect, idiag_movesect, llhysw, llwater, nnewsave, &
933 densdefault, densh2o, smallmassaa, smallmassbb, &
934 delta_water_conform1, delta_numb_conform1, &
935 drymassxx, dryvolxx, rmassxx, rnumbxx, wetmassxx, wetvolxx, &
936 drymassyy, dryvolyy, rmassyy, rnumbyy, wetmassyy, wetvolyy, &
937 xferfracvol, xferfracnum )
938 !
939 ! routine performs the actual transfer of aerosol number and mass
940 ! between sections
941 !
942 implicit none
943
944 ! include 'v33com'
945 ! include 'v33com2'
946 ! include 'v33com3'
947 ! include 'v33com9a'
948 ! include 'v33com9b'
949
950 ! subr arguments
951 integer iflag, iclm, jclm, iphase, itype, k, &
952 m, method_movesect, idiag_movesect, llhysw, llwater
953 integer nnewsave(2,maxd_asize)
954
955 real densdefault, densh2o, smallmassaa, smallmassbb
956 real delta_water_conform1, delta_numb_conform1
957 real drymassxx(maxd_asize), drymassyy(maxd_asize)
958 real dryvolxx(maxd_asize), dryvolyy(maxd_asize)
959 real rmassxx(maxd_acomp+2,maxd_asize), &
960 rmassyy(maxd_acomp+2,maxd_asize)
961 real rnumbxx(maxd_asize), rnumbyy(maxd_asize)
962 real xferfracvol(2,maxd_asize), xferfracnum(2,maxd_asize)
963 real wetvolxx(maxd_asize), wetvolyy(maxd_asize)
964 real wetmassxx(maxd_asize), wetmassyy(maxd_asize)
965
966 ! local variables
967 integer jj, l, ll, n, ndum, nnew, nold
968 integer jja, jjb, jjc
969
970 real dumbot, dumnum, dumnum_at_dhi, dumnum_at_dlo, &
971 dumvol, dumvol1p, dumxfvol, dumxfnum, sixoverpi, third
972 real dumxx(maxd_asize), dumyy(maxd_asize), dumout(maxd_asize)
973
974 character*160 msg
975 character*8 dumch8
976 character*4 dumch4
977
978
979 sixoverpi = 6.0/pi
980 third = 1.0/3.0
981
982 !
983 ! initialize "yy" arrays that hold values after inter-mode transferring
984 ! "yy" = "xx" for sections that do not move at all
985 ! "yy" = 0.0 for sections that do move (partially or completely)
986 !
987 do 1900 n = 1, nsize_aer(itype)
988
989 if ( (nnewsave(1,n) .eq. n) .and. &
990 (nnewsave(2,n) .eq. 0) ) then
991 ! if nnew == n, then material in section n will not be transferred, and
992 ! section n will contain its initial material plus any material
993 ! transferred from other sections
994 ! so initialize "yy" arrays with "xx" values
995 drymassyy(n) = drymassxx(n)
996 dryvolyy(n) = dryvolxx(n)
997 wetmassyy(n) = wetmassxx(n)
998 wetvolyy(n) = wetvolxx(n)
999 rnumbyy(n) = rnumbxx(n)
1000 do ll = 1, ncomp_plustracer_aer(itype) + 2
1001 rmassyy(ll,n) = rmassxx(ll,n)
1002 end do
1003
1004 else
1005 ! if nnew .ne. n, then material in section n will be transferred, and
1006 ! section n will only contain material that is transferred from
1007 ! other sections
1008 ! so initialize "yy" arrays to zero
1009 drymassyy(n) = 0.0
1010 dryvolyy(n) = 0.0
1011 wetmassyy(n) = 0.0
1012 wetvolyy(n) = 0.0
1013 rnumbyy(n) = 0.0
1014 do ll = 1, ncomp_plustracer_aer(itype) + 2
1015 rmassyy(ll,n) = 0.0
1016 end do
1017
1018 end if
1019
1020 1900 continue
1021
1022 !
1023 ! do the transfer of mass and number
1024 !
1025 do 2900 nold = 1, nsize_aer(itype)
1026
1027 if ( (nnewsave(1,nold) .eq. nold) .and. &
1028 (nnewsave(2,nold) .eq. 0 ) ) goto 2900
1029
1030 do 2800 jj = 1, 2
1031
1032 nnew = nnewsave(jj,nold)
1033 if (nnew .le. 0) goto 2800
1034
1035 dumxfvol = xferfracvol(jj,nold)
1036 dumxfnum = xferfracnum(jj,nold)
1037
1038 do ll = 1, ncomp_plustracer_aer(itype) + 2
1039 rmassyy(ll,nnew) = rmassyy(ll,nnew) + rmassxx(ll,nold)*dumxfvol
1040 end do
1041 rnumbyy(nnew) = rnumbyy(nnew) + rnumbxx(nold)*dumxfnum
1042
1043 drymassyy(nnew) = drymassyy(nnew) + drymassxx(nold)*dumxfvol
1044 dryvolyy(nnew) = dryvolyy(nnew) + dryvolxx(nold)*dumxfvol
1045 wetmassyy(nnew) = wetmassyy(nnew) + wetmassxx(nold)*dumxfvol
1046 wetvolyy(nnew) = wetvolyy(nnew) + wetvolxx(nold)*dumxfvol
1047
1048 2800 continue
1049
1050 2900 continue
1051
1052 !
1053 ! transfer among sections is completed
1054 ! - check for conservation of mass/volume/number
1055 ! - conform number again
1056 ! - compute/store densities and mean sizes
1057 ! - if k=1, save values for use by dry deposition routine
1058 ! - copy new mixrats back to rsub array
1059 !
1060 call move_sections_conserve_check( &
1061 1, iflag, iclm, jclm, k, m, iphase, itype, &
1062 method_movesect, idiag_movesect, llhysw, llwater, nnewsave, &
1063 densdefault, densh2o, smallmassaa, smallmassbb, &
1064 delta_water_conform1, delta_numb_conform1, &
1065 drymassxx, dryvolxx, rmassxx, rnumbxx, wetmassxx, wetvolxx, &
1066 drymassyy, dryvolyy, rmassyy, rnumbyy, wetmassyy, wetvolyy )
1067
1068 do 3900 n = 1, nsize_aer(itype)
1069
1070 dumvol = dryvolyy(n)
1071 if (min(drymassyy(n),dumvol) .le. smallmassbb) then
1072 dumvol = drymassyy(n)/densdefault
1073 rnumbyy(n) = dumvol/volumlo_sect(n,itype)
1074 adrydens_sub(n,k,m) = densdefault
1075 awetdens_sub(n,k,m) = densdefault
1076 admeandry_sub(n,k,m) = dcen_sect(n,itype)
1077 admeanwet_sub(n,k,m) = dcen_sect(n,itype)
1078 else
1079 dumnum = rnumbyy(n)
1080 dumnum_at_dhi = dumvol/volumhi_sect(n,itype)
1081 dumnum_at_dlo = dumvol/volumlo_sect(n,itype)
1082 dumnum = max( dumnum_at_dhi, min( dumnum_at_dlo, dumnum ) )
1083 rnumbyy(n) = dumnum
1084 adrydens_sub(n,k,m) = drymassyy(n)/dumvol
1085 dumvol1p = dumvol/dumnum
1086 admeandry_sub(n,k,m) = (dumvol1p*sixoverpi)**third
1087 awetdens_sub(n,k,m) = wetmassyy(n)/wetvolyy(n)
1088 dumvol1p = wetvolyy(n)/dumnum
1089 admeanwet_sub(n,k,m) = min( 100.*dcen_sect(n,itype), &
1090 (dumvol1p*sixoverpi)**third )
1091 end if
1092 avoldry_sub(n,k,m) = dumvol*cairclm(k)
1093 amassdry_sub(n,k,m) = drymassyy(n)*cairclm(k)
1094
1095 ! if (k .eq. 1) then
1096 ! awetdens_sfc(n,iclm,jclm) = awetdens_sub(n,k,m)
1097 ! admeanwet_sfc(n,iclm,jclm) = admeanwet_sub(n,k,m)
1098 ! end if
1099
1100 do ll = 1, ncomp_plustracer_aer(itype)
1101 l = massptr_aer(ll,n,itype,iphase)
1102 rsub(l,k,m) = rmassyy(ll,n)
1103 end do
1104 l = 0
1105 if (iphase .eq. ai_phase) then
1106 l = waterptr_aer(n,itype)
1107 if (l .gt. 0) rsub(l,k,m) = rmassyy(llwater,n)
1108 l = hyswptr_aer(n,itype)
1109 if (l .gt. 0) rsub(l,k,m) = rmassyy(llhysw,n)
1110 end if
1111 rsub(numptr_aer(n,itype,iphase),k,m) = rnumbyy(n)
1112
1113 3900 continue
1114
1115 call move_sections_conserve_check( &
1116 2, iflag, iclm, jclm, k, m, iphase, itype, &
1117 method_movesect, idiag_movesect, llhysw, llwater, nnewsave, &
1118 densdefault, densh2o, smallmassaa, smallmassbb, &
1119 delta_water_conform1, delta_numb_conform1, &
1120 drymassxx, dryvolxx, rmassxx, rnumbxx, wetmassxx, wetvolxx, &
1121 drymassyy, dryvolyy, rmassyy, rnumbyy, wetmassyy, wetvolyy )
1122
1123
1124 ! diagnostic output
1125 ! output nnewsave values when msectional=7xxx
1126 if (idiag_movesect .ne. 70) goto 4900
1127
1128 ndum = 0
1129 do n = 1, nsize_aer(itype)
1130 if (nnewsave(1,n)+nnewsave(2,n) .ne. n) ndum = ndum + 1
1131 end do
1132 ! if (ndum .gt. 0) then
1133 ! write(msg,97010) 'SOME', iclm, jclm, k, m, &
1134 ! ndum, (nnewsave(1,n), nnewsave(2,n), n=1,nsize_aer(itype))
1135 ! call peg_message( lunout, msg )
1136 ! else
1137 ! write(msg,97010) 'NONE', iclm, jclm, k, m, &
1138 ! ndum, (nnewsave(1,n), nnewsave(2,n), n=1,nsize_aer(itype))
1139 ! call peg_message( lunout, msg )
1140 ! end if
1141
1142 dumch4 = 'NONE'
1143 if (ndum .gt. 0) dumch4 = 'SOME'
1144 msg = ' '
1145 call peg_message( lunout, msg )
1146 write(msg,97010) dumch4, iclm, jclm, k, m, ndum
1147 call peg_message( lunout, msg )
1148 do jjb = 1, nsize_aer(itype), 10
1149 jjc = min( jjb+9, nsize_aer(itype) )
1150 write(msg,97011) (nnewsave(1,n), nnewsave(2,n), n=jjb,jjc)
1151 call peg_message( lunout, msg )
1152 end do
1153
1154 ! write(msg,97020) 'rnumbold', (rnumbxx(n), n=1,nsize_aer(itype))
1155 ! call peg_message( lunout, msg )
1156 ! write(msg,97020) 'rnumbnew', (rnumbyy(n), n=1,nsize_aer(itype))
1157 ! call peg_message( lunout, msg )
1158
1159 ! write(msg,97020) 'drvolold', (dryvolxx(n), n=1,nsize_aer(itype))
1160 ! call peg_message( lunout, msg )
1161 ! write(msg,97020) 'drvolnew', (dryvolyy(n), n=1,nsize_aer(itype))
1162 ! call peg_message( lunout, msg )
1163
1164 dumbot = log( volumhi_sect(1,itype)/volumlo_sect(1,itype) )
1165 do n = 1, nsize_aer(itype)
1166 dumxx(n) = -9.99
1167 dumyy(n) = -9.99
1168 if ((rnumbxx(n) .ge. 1.0e-35) .and. &
1169 (dryvolxx(n) .ge. 1.0e-35)) then
1170 dumvol1p = dryvolxx(n)/rnumbxx(n)
1171 dumxx(n) = 1.0 + log(dumvol1p/volumlo_sect(1,itype))/dumbot
1172 end if
1173 if ((rnumbyy(n) .ge. 1.0e-35) .and. &
1174 (dryvolyy(n) .ge. 1.0e-35)) then
1175 dumvol1p = dryvolyy(n)/rnumbyy(n)
1176 dumyy(n) = 1.0 + log(dumvol1p/volumlo_sect(1,itype))/dumbot
1177 end if
1178 end do
1179
1180 ! write(msg,97030) 'lnvolold', (dumxx(n), n=1,nsize_aer(itype))
1181 ! call peg_message( lunout, msg )
1182 ! write(msg,97030) 'lnvolnew', (dumyy(n), n=1,nsize_aer(itype))
1183 ! call peg_message( lunout, msg )
1184
1185 do jja = 1, 6
1186 if (jja .eq. 1) then
1187 dumch8 = 'rnumbold'
1188 dumout(:) = rnumbxx(:)
1189 else if (jja .eq. 2) then
1190 dumch8 = 'rnumbnew'
1191 dumout(:) = rnumbyy(:)
1192 else if (jja .eq. 3) then
1193 dumch8 = 'drvolold'
1194 dumout(:) = dryvolxx(:)
1195 else if (jja .eq. 4) then
1196 dumch8 = 'drvolnew'
1197 dumout(:) = dryvolyy(:)
1198 else if (jja .eq. 5) then
1199 dumch8 = 'lnvolold'
1200 dumout(:) = dumxx(:)
1201 else if (jja .eq. 6) then
1202 dumch8 = 'lnvolnew'
1203 dumout(:) = dumyy(:)
1204 end if
1205 do jjb = 1, nsize_aer(itype), 10
1206 jjc = min( jjb+9, nsize_aer(itype) )
1207 if (jja .le. 4) then
1208 write(msg,97020) dumch8, (dumout(n), n=jjb,jjc)
1209 else
1210 write(msg,97030) dumch8, (dumout(n), n=jjb,jjc)
1211 end if
1212 call peg_message( lunout, msg )
1213 dumch8 = ' '
1214 end do
1215 end do
1216
1217 !97010 format( / 'movesectapply', a, 4i3, 3x, i3 / 5x, 10(3x,2i3) )
1218 !97020 format( a, 1p, 10e9.1 / (( 8x, 1p, 10e9.1 )) )
1219 !97030 format( a, 10f9.3 / (( 8x, 10f9.3 )) )
1220 97010 format( 'movesectapply', a, 4i3, 3x, i3 )
1221 97011 format( 5x, 10(3x,2i3) )
1222 97020 format( a, 1p, 10e9.1 )
1223 97030 format( a, 10f9.3 )
1224
1225 4900 continue
1226 return
1227 end subroutine move_sections_apply_moves
1228
1229
1230 !-----------------------------------------------------------------------
1231 ! rce 08-nov-2004 - this routine is new (and perhaps temporary)
1232 !
1233 subroutine move_sections_apply_n1_inflow( &
1234 iflag, iclm, jclm, k, m, iphase, itype, &
1235 method_movesect, idiag_movesect, llhysw, llwater, nnewsave, &
1236 densdefault, densh2o, smallmassaa, smallmassbb, &
1237 delta_water_conform1, delta_numb_conform1, &
1238 drymassxx, dryvolxx, rmassxx, rnumbxx, wetmassxx, wetvolxx, &
1239 drymassyy, dryvolyy, rmassyy, rnumbyy, wetmassyy, wetvolyy, &
1240 xferfracvol, xferfracnum, &
1241 specmwxx, specdensxx )
1242 !
1243 ! routine applies an ad_hoc correction to bin 1 to account for
1244 ! growth of smaller particles (which are not simulated) growing into
1245 ! bin 1
1246 ! the correction to particle number balances the loss of particles from
1247 ! bin 1 to larger bins by growth
1248 ! the correction to mass assumes that the particles coming into bin 1
1249 ! are slightly larger than dlo_sect(1)
1250 !
1251 implicit none
1252
1253 ! include 'v33com'
1254 ! include 'v33com2'
1255 ! include 'v33com3'
1256 ! include 'v33com9a'
1257 ! include 'v33com9b'
1258
1259 ! subr arguments
1260 integer iflag, iclm, jclm, iphase, itype, k, &
1261 m, method_movesect, idiag_movesect, llhysw, llwater
1262 integer nnewsave(2,maxd_asize)
1263
1264 real densdefault, densh2o, smallmassaa, smallmassbb
1265 real delta_water_conform1, delta_numb_conform1
1266 real drymassxx(maxd_asize), drymassyy(maxd_asize)
1267 real dryvolxx(maxd_asize), dryvolyy(maxd_asize)
1268 real rmassxx(maxd_acomp+2,maxd_asize), &
1269 rmassyy(maxd_acomp+2,maxd_asize)
1270 real rnumbxx(maxd_asize), rnumbyy(maxd_asize)
1271 real xferfracvol(2,maxd_asize), xferfracnum(2,maxd_asize)
1272 real wetvolxx(maxd_asize), wetvolyy(maxd_asize)
1273 real wetmassxx(maxd_asize), wetmassyy(maxd_asize)
1274 real specdensxx(maxd_acomp), specmwxx(maxd_acomp)
1275
1276
1277 ! local variables
1278 integer jj, l, ll, n, nnew, nold
1279
1280 real deltanum, deltavol, dumvol1p, dumxfvol, dumxfnum
1281
1282
1283 !
1284 ! compute fraction of number transferred out of bin 1 by growth
1285 !
1286 nold = 1
1287 n = nold
1288 if ( (nnewsave(1,nold) .eq. nold) .and. &
1289 (nnewsave(2,nold) .eq. 0 ) ) goto 3900
1290
1291 dumxfnum = 0.0
1292 do 2800 jj = 1, 2
1293 nnew = nnewsave(jj,nold)
1294 if (nnew .le. 0) goto 2800
1295 if (nnew .eq. nold) goto 2800
1296 dumxfnum = dumxfnum + xferfracnum(jj,nold)
1297 2800 continue
1298
1299 !
1300 ! compute "inflow" change to number and volume
1301 ! number change matches that lost by growth
1302 ! volume change assume inflow particles slightly bigger then dlo_sect
1303 !
1304 dumvol1p = 0.95*volumlo_sect(n,itype) + 0.05*volumhi_sect(n,itype)
1305 deltanum = dumxfnum*rnumbxx(n)
1306 deltavol = deltanum*dumvol1p
1307 if (dumxfnum .le. 0.0) goto 3900
1308 if (deltanum .le. 0.0) goto 3900
1309 if (deltavol .le. 0.0) goto 3900
1310
1311 !
1312 ! increment the number and masses
1313 ! if the old dryvol (dryvolxx) > smallmassbb, then compute mass increment for
1314 ! each species from the old masses (rmassxx)
1315 ! otherwise only increment the first mass species
1316 !
1317 if (dryvolxx(n) .gt. smallmassbb) then
1318 dumxfvol = deltavol/dryvolxx(n)
1319 do ll = 1, ncomp_plustracer_aer(itype) + 2
1320 rmassyy(ll,n) = rmassyy(ll,n) + dumxfvol*rmassxx(ll,n)
1321 end do
1322 else
1323 ll = 1
1324 rmassyy(ll,n) = rmassyy(ll,n) + deltavol*specdensxx(ll)/specmwxx(ll)
1325 end if
1326 rnumbyy(n) = rnumbyy(n) + deltanum
1327
1328 !
1329 ! transfer results into rsub
1330 !
1331 do ll = 1, ncomp_plustracer_aer(itype)
1332 l = massptr_aer(ll,n,itype,iphase)
1333 rsub(l,k,m) = rmassyy(ll,n)
1334 end do
1335 if (iphase .eq. ai_phase) then
1336 l = waterptr_aer(n,itype)
1337 if (l .gt. 0) rsub(l,k,m) = rmassyy(llwater,n)
1338 l = hyswptr_aer(n,itype)
1339 if (l .gt. 0) rsub(l,k,m) = rmassyy(llhysw,n)
1340 end if
1341 rsub(numptr_aer(n,itype,iphase),k,m) = rnumbyy(n)
1342
1343
1344 3900 continue
1345 return
1346 end subroutine move_sections_apply_n1_inflow
1347
1348
1349 !-----------------------------------------------------------------------
1350 subroutine move_sections_conserve_check( ipass, &
1351 iflag, iclm, jclm, k, m, iphase, itype, &
1352 method_movesect, idiag_movesect, llhysw, llwater, nnewsave, &
1353 densdefault, densh2o, smallmassaa, smallmassbb, &
1354 delta_water_conform1, delta_numb_conform1, &
1355 drymassxx, dryvolxx, rmassxx, rnumbxx, wetmassxx, wetvolxx, &
1356 drymassyy, dryvolyy, rmassyy, rnumbyy, wetmassyy, wetvolyy )
1357 !
1358 ! routine checks for conservation of number, mass, and volume
1359 ! by the move_sections algorithm
1360 !
1361 ! ipass = 1
1362 ! initialize all thesum(jj,ll) to zero
1363 ! computes thesum(1,ll) from rmassxx, rnumbxx, ...
1364 ! computes thesum(2,ll) from rmassyy, rnumbyy, ...
1365 ! compares thesum(1,ll) with thesum(2,ll)
1366 ! computes thesum(3,ll) from rsub before section movement
1367 ! ipass = 2
1368 ! computes thesum(4,ll) from rsub after section movement
1369 ! compares thesum(3,ll) with thesum(4,ll)
1370 !
1371 ! currently only implemented for condensational growth (iflag=1)
1372 !
1373 implicit none
1374
1375 ! include 'v33com'
1376 ! include 'v33com2'
1377 ! include 'v33com3'
1378 ! include 'v33com9a'
1379 ! include 'v33com9b'
1380
1381 ! subr arguments
1382 integer ipass, iflag, iclm, jclm, iphase, itype, k, &
1383 m, method_movesect, idiag_movesect, llhysw, llwater
1384 integer nnewsave(2,maxd_asize)
1385 real densdefault, densh2o, smallmassaa, smallmassbb
1386 real delta_water_conform1, delta_numb_conform1
1387 real drymassxx(maxd_asize), drymassyy(maxd_asize)
1388 real dryvolxx(maxd_asize), dryvolyy(maxd_asize)
1389 real rmassxx(maxd_acomp+2,maxd_asize), &
1390 rmassyy(maxd_acomp+2,maxd_asize)
1391 real rnumbxx(maxd_asize), rnumbyy(maxd_asize)
1392 real wetvolxx(maxd_asize), wetvolyy(maxd_asize)
1393 real wetmassxx(maxd_asize), wetmassyy(maxd_asize)
1394
1395 ! local variables
1396 integer jj, l, ll, llworst, llworstb, n
1397 integer nerr, nerrmax
1398 save nerr, nerrmax
1399 data nerr, nerrmax / 0, 999 /
1400
1401 real dumbot, dumtop, dumtoler, dumerr, dumworst, dumworstb
1402 real thesum(4,maxd_acomp+7)
1403 save thesum
1404
1405 character*8 dumname(maxd_acomp+7)
1406 character*160 msg
1407
1408
1409 if (ipass .eq. 2) goto 2000
1410
1411 do ll = 1, ncomp_plustracer_aer(itype)+7
1412 do jj = 1, 4
1413 thesum(jj,ll) = 0.0
1414 end do
1415 end do
1416
1417 do n = 1, nsize_aer(itype)
1418 do ll = 1, ncomp_plustracer_aer(itype)+2
1419 thesum(1,ll) = thesum(1,ll) + rmassxx(ll,n)
1420 thesum(2,ll) = thesum(2,ll) + rmassyy(ll,n)
1421 end do
1422 ll = ncomp_plustracer_aer(itype)+3
1423 thesum(1,ll) = thesum(1,ll) + rnumbxx(n)
1424 thesum(2,ll) = thesum(2,ll) + rnumbyy(n)
1425 ll = ncomp_plustracer_aer(itype)+4
1426 thesum(1,ll) = thesum(1,ll) + drymassxx(n)
1427 thesum(2,ll) = thesum(2,ll) + drymassyy(n)
1428 ll = ncomp_plustracer_aer(itype)+5
1429 thesum(1,ll) = thesum(1,ll) + dryvolxx(n)
1430 thesum(2,ll) = thesum(2,ll) + dryvolyy(n)
1431 ll = ncomp_plustracer_aer(itype)+6
1432 thesum(1,ll) = thesum(1,ll) + wetmassxx(n)
1433 thesum(2,ll) = thesum(2,ll) + wetmassyy(n)
1434 ll = ncomp_plustracer_aer(itype)+7
1435 thesum(1,ll) = thesum(1,ll) + wetvolxx(n)
1436 thesum(2,ll) = thesum(2,ll) + wetvolyy(n)
1437 end do
1438
1439
1440 2000 continue
1441 !
1442 ! calc sum over bins for each species
1443 ! for water, account for loss in initial conform (delta_water_conform1)
1444 !
1445 do n = 1, nsize_aer(itype)
1446 do ll = 1, ncomp_plustracer_aer(itype)+3
1447 if (ll .le. ncomp_plustracer_aer(itype)) then
1448 l = massptr_aer(ll,n,itype,iphase)
1449 else if (ll .eq. ncomp_plustracer_aer(itype)+1) then
1450 l = 0
1451 if (iphase .eq. ai_phase) l = hyswptr_aer(n,itype)
1452 else if (ll .eq. ncomp_plustracer_aer(itype)+2) then
1453 l = 0
1454 if (iphase .eq. ai_phase) l = waterptr_aer(n,itype)
1455 else
1456 l = numptr_aer(n,itype,iphase)
1457 end if
1458 if (l .gt. 0) &
1459 thesum(ipass+2,ll) = thesum(ipass+2,ll) + rsub(l,k,m)
1460 end do
1461 end do
1462 if (ipass .eq. 1) then
1463 ll = ncomp_plustracer_aer(itype)+2
1464 thesum(3,ll) = thesum(3,ll) + delta_water_conform1
1465 ll = ncomp_plustracer_aer(itype)+3
1466 thesum(3,ll) = thesum(3,ll) + delta_numb_conform1
1467 end if
1468
1469 !
1470 ! now compare either sum1-sum2 or sum3-sum4
1471 ! on ipass=1, jj=1, so compare sum1 & sum2
1472 ! on ipass=2, jj=3, so compare sum3 & sum4
1473 !
1474 do ll = 1, ncomp_plustracer_aer(itype)+7
1475 dumname(ll) = ' '
1476 write(dumname(ll),'(i4.4)') ll
1477 if (ll .le. ncomp_plustracer_aer(itype)) dumname(ll) = &
1478 name(massptr_aer(ll,1,itype,iphase))(1:4)
1479 if (ll .eq. ncomp_plustracer_aer(itype)+1) dumname(ll) = 'hysw'
1480 if (ll .eq. ncomp_plustracer_aer(itype)+2) dumname(ll) = 'watr'
1481 if (ll .eq. ncomp_plustracer_aer(itype)+3) dumname(ll) = 'numb'
1482 if (ll .eq. ncomp_plustracer_aer(itype)+4) dumname(ll) = 'drymass'
1483 if (ll .eq. ncomp_plustracer_aer(itype)+5) dumname(ll) = 'dryvol'
1484 if (ll .eq. ncomp_plustracer_aer(itype)+6) dumname(ll) = 'wetmass'
1485 if (ll .eq. ncomp_plustracer_aer(itype)+7) dumname(ll) = 'wetvol'
1486 end do
1487
1488 jj = 2*ipass - 1
1489 dumworst = 0.0
1490 dumworstb = 0.0
1491 llworst = 0
1492 llworstb = 0
1493 do ll = 1, ncomp_plustracer_aer(itype)+7
1494 dumtop = thesum(jj+1,ll) - thesum(jj,ll)
1495 dumbot = max( abs(thesum(jj,ll)), abs(thesum(jj+1,ll)), 1.0e-35 )
1496 dumerr = dumtop/dumbot
1497 if (abs(dumerr) .gt. abs(dumworst)) then
1498 llworstb = llworst
1499 dumworstb = dumworst
1500 llworst = ll
1501 dumworst = dumerr
1502 end if
1503 end do
1504
1505 dumtoler = 5.0e-7
1506 dumtoler = 9.0e-7
1507 if (abs(dumworst) .gt. dumtoler) then
1508 nerr = nerr + 1
1509 if (nerr .le. nerrmax) then
1510 msg = ' '
1511 call peg_message( lunout, msg )
1512 write(msg,97110) iclm, jclm, k, m, ipass, llworst
1513 call peg_message( lunout, msg )
1514 write(msg,97120) ' nnew(1,n)', &
1515 (nnewsave(1,n), n=1,nsize_aer(itype))
1516 call peg_message( lunout, msg )
1517 write(msg,97120) ' nnew(2,n)', &
1518 (nnewsave(2,n), n=1,nsize_aer(itype))
1519 call peg_message( lunout, msg )
1520
1521 ll = llworst
1522 if (ll .eq. 0) ll = ncomp_plustracer_aer(itype)+2
1523 write(msg,97130) 'name/relerr/thesum', jj, '/thesum', jj+1, &
1524 dumname(ll), dumworst, thesum(jj,ll), thesum(jj+1,ll)
1525 call peg_message( lunout, msg )
1526
1527 if ( (ll .eq. ncomp_plustracer_aer(itype)+3) .and. &
1528 (abs(dumworstb) .gt. dumtoler) ) then
1529 ll = max( 1, llworstb )
1530 dumtop = thesum(jj+1,ll) - thesum(jj,ll)
1531 dumbot = max( abs(thesum(jj,ll)), abs(thesum(jj+1,ll)), 1.0e-35 )
1532 dumerr = dumtop/dumbot
1533 write(msg,97130) 'name/relerr/thesum', jj, '/thesum', jj+1, &
1534 dumname(ll), dumerr, thesum(jj,ll), thesum(jj+1,ll)
1535 call peg_message( lunout, msg )
1536 end if
1537 end if
1538 end if
1539
1540 97110 format( 'movesect conserve ERROR - i/j/k/m/pass/llworst', &
1541 4i3, 2x, 2i3 )
1542 97120 format( a, 64i3 )
1543 97130 format( a, i1, a, i1, 2x, a, 1p, 3e16.7 )
1544
1545 return
1546 end subroutine move_sections_conserve_check
1547
1548
1549 !-----------------------------------------------------------------------
1550 subroutine test_move_sections( iflag_test, iclm, jclm, k, m )
1551 !
1552 ! routine runs tests on move_sections, using a matrix of
1553 ! pregrow and aftgrow masses for each section
1554 !
1555 implicit none
1556
1557 ! include 'v33com'
1558 ! include 'v33com2'
1559 ! include 'v33com3'
1560 ! include 'v33com9a'
1561 ! include 'v33com9b'
1562
1563 ! subr arguments
1564 integer iflag_test, iclm, jclm, k, m
1565
1566 ! local variables
1567 integer idiag_movesect, iflag, ii, iphase, itype, jj, &
1568 l, ll, n, nn
1569 integer ientryno
1570 save ientryno
1571 data ientryno / 0 /
1572
1573 real dumnumb, dumvolpre, dumvolaft
1574 real dumsv_rsub(l2maxd)
1575 real dumsv_drymass_pregrow(maxd_asize,maxd_atype)
1576 real dumsv_drymass_aftgrow(maxd_asize,maxd_atype)
1577 real dumsv_drydens_pregrow(maxd_asize,maxd_atype)
1578 real dumsv_drydens_aftgrow(maxd_asize,maxd_atype)
1579
1580 character*160 msg
1581
1582 integer maxvolfactpre, maxvolfactaft
1583 parameter (maxvolfactpre=15, maxvolfactaft=23)
1584
1585 real dumvolfactpre(maxvolfactpre)
1586 data dumvolfactpre / &
1587 2.0, 0.0, 1.0e-20, 0.5, 0.9, &
1588 1.0, 1.01, 1.1, 2.0, 4.0, 7.9, 7.99, 8.0, &
1589 8.1, 16.0 /
1590
1591 real dumvolfactaft(maxvolfactaft)
1592 data dumvolfactaft / &
1593 4.0, 0.0, 1.0e-20, 0.01, 0.02, 0.05, 0.1, 0.5, 0.9, &
1594 1.0, 1.01, 1.1, 2.0, 4.0, 7.9, 7.99, 8.0, &
1595 8.1, 16.0, 32., 64., 128., 256. /
1596
1597
1598 !
1599 ! check for valid inputs
1600 ! and first entry
1601 !
1602 if (msectional .le. 0) return
1603 if (nsize_aer(1) .le. 0) return
1604
1605 idiag_movesect = mod( msectional, 10000 )/100
1606 if (idiag_movesect .ne. 70) return
1607
1608 ientryno = ientryno + 1
1609 if (ientryno .gt. 2) return
1610
1611 !
1612 ! save rsub and drymass/dens_aft/pregrow
1613 !
1614 do l = 1, ltot2
1615 dumsv_rsub(l) = rsub(l,k,m)
1616 end do
1617 do itype = 1, ntype_aer
1618 do n = 1, nsize_aer(itype)
1619 dumsv_drymass_pregrow(n,itype) = drymass_pregrow(n,itype)
1620 dumsv_drymass_aftgrow(n,itype) = drymass_aftgrow(n,itype)
1621 dumsv_drydens_pregrow(n,itype) = drydens_pregrow(n,itype)
1622 dumsv_drydens_aftgrow(n,itype) = drydens_aftgrow(n,itype)
1623 end do
1624 end do
1625
1626 !
1627 ! make test calls to move_sections
1628 !
1629 do 3900 iflag = 1, 1
1630
1631 iphase = ai_phase
1632 if (iabs(iflag) .eq. 2) iphase = cw_phase
1633
1634 do 3800 itype = 1, ntype_aer
1635
1636 do 2900 nn = 1, nsize_aer(itype)
1637
1638 do 2800 ii = 1, maxvolfactpre
1639
1640 do 2700 jj = 1, maxvolfactaft
1641
1642 ! zero out rsub and dryxxx_yyygrow arrays
1643 do n = 1, nsize_aer(itype)
1644 do ll = 1, ncomp_plustracer_aer(itype)
1645 rsub(massptr_aer(ll,n,itype,iphase),k,m) = 0.0
1646 end do
1647 l = 0
1648 if (iphase .eq. ai_phase) l = waterptr_aer(n,itype)
1649 if (l .gt. 0) rsub(l,k,m) = 0.0
1650 l = numptr_aer(n,itype,iphase)
1651 if (l .gt. 0) rsub(l,k,m) = 0.0
1652 drymass_pregrow(n,itype) = 0.0
1653 drymass_aftgrow(n,itype) = 0.0
1654 drydens_pregrow(n,itype) = -1.0
1655 drydens_aftgrow(n,itype) = -1.0
1656 end do
1657
1658 ! fill in values for section nn
1659 n = nn
1660 dumnumb = 1.0e7
1661 rsub(numptr_aer(n,itype,iphase),k,m) = dumnumb
1662 ll = 1
1663 l = massptr_aer(ll,n,itype,iphase)
1664
1665 dumvolpre = volumlo_sect(n,itype)*dumvolfactpre(ii)*dumnumb
1666 drydens_pregrow(n,itype) = dens_aer(ll,itype)
1667 drymass_pregrow(n,itype) = dumvolpre*drydens_pregrow(n,itype)
1668 if (ii .eq. 1) drydens_pregrow(n,itype) = -1.0
1669
1670 dumvolaft = volumlo_sect(n,itype)*dumvolfactaft(jj)*dumnumb
1671 drydens_aftgrow(n,itype) = dens_aer(ll,itype)
1672 drymass_aftgrow(n,itype) = dumvolaft*drydens_aftgrow(n,itype)
1673 if (jj .eq. 1) drydens_aftgrow(n,itype) = -1.0
1674
1675 rsub(l,k,m) = drymass_aftgrow(n,itype)/mw_aer(ll,itype)
1676
1677 msg = ' '
1678 call peg_message( lunout, msg )
1679 write(msg,98010) nn, ii, jj
1680 call peg_message( lunout, msg )
1681 write(msg,98011) dumvolfactpre(ii), dumvolfactaft(jj)
1682 call peg_message( lunout, msg )
1683
1684 ! make test call to move_sections
1685 call move_sections( iflag, iclm, jclm, k, m )
1686
1687 msg = ' '
1688 call peg_message( lunout, msg )
1689 write(msg,98010) nn, ii, jj
1690 call peg_message( lunout, msg )
1691 write(msg,98011) dumvolfactpre(ii), dumvolfactaft(jj)
1692 call peg_message( lunout, msg )
1693
1694 2700 continue
1695 2800 continue
1696 2900 continue
1697
1698 3800 continue
1699 3900 continue
1700
1701 98010 format( 'test_move_sections output - nn, ii, jj =', 3i3 )
1702 98011 format( 'volfactpre, volfactaft =', 1p, 2e12.4 )
1703
1704 !
1705 ! restore rsub and drymass/dens_aft/pregrow
1706 !
1707 do l = 1, ltot2
1708 rsub(l,k,m) = dumsv_rsub(l)
1709 end do
1710 do n = 1, nsize_aer(itype)
1711 drymass_pregrow(n,itype) = dumsv_drymass_pregrow(n,itype)
1712 drymass_aftgrow(n,itype) = dumsv_drymass_aftgrow(n,itype)
1713 drydens_pregrow(n,itype) = dumsv_drydens_pregrow(n,itype)
1714 drydens_aftgrow(n,itype) = dumsv_drydens_aftgrow(n,itype)
1715 end do
1716
1717 return
1718 end subroutine test_move_sections
1719
1720
1721 !-----------------------------------------------------------------------
1722 subroutine move_sections_checkptrs( iflag, iclm, jclm, k, m )
1723 !
1724 ! checks for valid number and water pointers
1725 !
1726 implicit none
1727
1728 ! include 'v33com'
1729 ! include 'v33com2'
1730 ! include 'v33com3'
1731 ! include 'v33com9a'
1732 ! include 'v33com9b'
1733
1734 ! subr parameters
1735 integer iflag, iclm, jclm, k, m
1736
1737 ! local variables
1738 integer l, itype, iphase, n, ndum
1739 character*160 msg
1740
1741 do 1900 itype = 1, ntype_aer
1742 do 1800 iphase = 1, nphase_aer
1743
1744 ndum = 0
1745 do n = 1, nsize_aer(itype)
1746 l = numptr_aer(n,itype,iphase)
1747 if ((l .le. 0) .or. (l .gt. ltot2)) then
1748 msg = '*** subr move_sections error - ' // &
1749 'numptr_amode not defined'
1750 call peg_message( lunerr, msg )
1751 write(msg,9030) 'mode, numptr =', n, l
1752 call peg_message( lunerr, msg )
1753 write(msg,9030) 'iphase, itype =', iphase, itype
1754 call peg_message( lunerr, msg )
1755 call peg_error_fatal( lunerr, msg )
1756 end if
1757 ! checks involving nspec_amode and nspec_amode_nontracer
1758 ! being the same for all sections are no longer needed
1759 l = 0
1760 if (iphase .eq. ai_phase) l = waterptr_aer(n,itype)
1761 if ((l .gt. 0) .and. (l .le. ltot2)) ndum = ndum + 1
1762 end do
1763 if ((ndum .ne. 0) .and. (ndum .ne. nsize_aer(itype))) then
1764 msg = '*** subr move_sections error - ' // &
1765 'waterptr_aer must be on/off for all modes'
1766 call peg_message( lunerr, msg )
1767 write(msg,9030) 'iphase, itype =', iphase, itype
1768 call peg_message( lunerr, msg )
1769 call peg_error_fatal( lunerr, msg )
1770 end if
1771 9030 format( a, 2(1x,i6) )
1772
1773 1800 continue
1774 1900 continue
1775
1776 return
1777 end subroutine move_sections_checkptrs
1778
1779
1780
1781 end module module_mosaic_movesect