Demolition Dust Model

by **Buddy McCloskey **** **
4/02/2008 **Contact****
**

**
1. Purpose:
**

**
The goals are mentioned in the preceding article World Trade Center Dust Bowl Analysis.
We want to extend the results discovered in the article, to apply to any building, if certain parameters are known,
then extend to any major city seen in 9 below. As a summary, the WTC Dust Bowl Model was
to predict quantity and region of possible dust settlement. Another more important function
will be to extend this model to a multitude of buildings in metropolitan area to estimate
loose dust generated that will be sucked into a vortex of a mushroom cloud during a nuclear
blast so more accurate estimation of downwind radioactive fallout can be achieved in case of
a nuclear attack on our major cities. City profiles of buildings heights and spacing density
will be approximated and a total dust weight and volume from each city's sampled spectrum of
buildings will be used. These tons of dust will be added to the ground crater's dust blown
out by the mechanical energy from the nuclear blast as input into the radioactive
fallout model RADFO as dust inventory.
**

**
2. Building Parameters:
**

**
For the nuclear, general conversion constants will be estimated, many of which will be the WTC
findings. For the dust fall from a planned building demolition, more detailed and actual values
should be used. In any case, the parameters derived or actual are listed below:
**

Densities | g/cm^3 |
---|---|

Accepted of Cement | 2.4 |

Structural Light Weight Concrete | 1.84 (WTC) |

My Cement Block Rub | 2.128 |

Accepted Average Cement | 2.12 |

Reinforced Steel | 7.85 |

% - Ratios | |
---|---|

Cement/Building by Volume | 3.461% |

Steel/Cement by Volume | 25.0021% |

Actual/Gravity Rate of Collapse | 1.6284 |

Max Fluff Factor | 1.873 |

Building Dimensions | Units |
---|---|

Height | ft |

Length | ft |

Width | ft |

Floor Slab Thickness ? | 2.5" in WTC |

Wall Thickness ? | in |

No of Floors ? | # |

Cement Weight * | Tons |

Cement Volume * | ft^3 |

Special Derived Constants | Value |
---|---|

Block Rub # 3 Pressure | .5236 PSI |

Block Rub # 3 Width | 4.25 in |

Block Rub # 3 Rub Length | 98.4 ft |

Block Rub # 3 Tons/Sec | 1.767e-8 |

Block Rub # 4/# 3 Pressure Ratio | 2.38861 |

Block Rub # 4/# 3 Grams Ratio | 6.01928 |

**
3. Dust by Physical Processes: - B%D.c
**

**
In demolition, we assume the base of the building is destroyed and falls upon
itself. It might be best to look at this general process here in trap3.
If these displays are shifted and jumbled when you open them up with your web
browser processor, open Word Pad or some other word processor, open them in MS
Accessories - Notepad. If by default you can't, you might right click and save
the display as a text file, then use Notepad.
The general idea is that the building walls are crushing the walls below them,
cumulating rubble that falls equally on the floor slab within the walls and outside
the walls - cumulating a large pile on the outside that rests against the accumulated
slab heights - the "pancake" effect often mentioned. The shape is a rectangular prism
with base radius Rc, height Hw in trap4.
The R's for radius rubble and to represent another term later. Notice, the building
and the next floor falling down.
Note the more detail in trap4, we note two processes. We see the chunks of cement
and the airborne dust is divided equally to either side of the crushing wall. This is indicated
by Dust # 1. Here, we have the bottom of the wall sliding in both directions < - >, causing dust
by the rubbing surfaces under great pressures, but assuming here that is not enough to greatly exceed the compression
limit of 6000 psi because the rubbing changes it to powder dust instantly. We calculate the dust by the wall edges rubbing,
and the remaining chunks are divided equally between the outside and inside piles. The first results of compression is the
crumbling of chunks caused by > 6000 psi on jagged edges.
In the lower diagram, we see the secondary processes at work indicated by Dust
#2 & #3. Now Dust #2 is caused by rubbing against the outside cone pile or lateral
surface area down side (hypotenuse L) of the triangle. We show spaces here but actually,
the chunks are continuous extension of the wall with the total pressure
still on the chunks as they traverse the outside wall, because the time between slabs
falls is less than .15 sec or 62 mph.
To the right, we have chunks accumulating above each floor slab (tbs) until the floor slab
above meets it at height Ht (height temporarily trapped). When the slab falls on each
of it's respective pile, it's the crushing action by the crush power formula in use
to produce dust #3 to combine with dust #1. Both of these are extruded with compressed
air puffs explained below. If there were no smaller pile, the slabs would stack like pancakes, and land flat,
producing no side rubbing. If they did slide, like the walls, the area is immensely
larger but the psi is less than 50 initially so not as much dust produced as one might
think.
**

**
4. Deriving Other Parameters from the Building Parameters:
A. Dust Puffs
**

**
Many derived parameters were developed in a forerunner program before the final
model here, not knowing if the data was needed but interesting to shed more light
on the physical process. We will present most of them here and tell if was used in
the final or not. Actually, we can view the data in the summary stat sheet ddmdat.
As a standard, we expand the WTC study again **

**
B. Separating Slabs and Walls
**

**
We treat the height of the rubble pile as the accumulated slab heights = 22.92'.
The volume of each slab is 2.5" times length x width = 9.915e+5 ft^3. Subtracting
this from the total cement volume leaves us the total volume of the walls. Pretend,
we unravel the enclosed building area by stacking the 4 walls end to end (perimeter)
at building height. Now to find the wall thickness, we divide the volume by the area
perimeter x height to find the thickness of this new shaped volume. However, a minor
adjustment, we must subtract off the heights of each slab or the stacked slab heights
here for the new height - walls with no slabs. This gives us a wall thickness of 11.3".
Multiplying the cement density by volume of each slab = 517.686 Tons and each floor's
wall volume = 9.383 Tons.
**

**
C. Rub Areas
**

**
Also needed is the rub areas of the wall. Since the slabs are accumulating vertically,
we assume no horizontal rubbing of the slabs. The wall rubs on itself and the residual
goes inside and outside the wall. Simply, the area is the wall thickness x perimeter
= 784.18'^2 (wtcaw). The rub 3 block rub length is 98.4' x it's width 4.25" from the above
table gives us an rub3a area of 34.85'^2 so the top of the wall area ratio is rarb=wtcaw/rub3a =
784.18/34.85 = 22.50 times as much dust as the rub 3 block area (rub3a) which by hand rubbing
alone would be minute.
**

**
D. Rub Time
**

**
We loop up the building height by the slab thickness. This gives us a better approximation
for the wall dimensions but allows us to evaluate the whole slab when it's turn comes. So for the WTC,
1365'/(2.5/12) = 7872 evaluations. If it takes 15 sec for the building to fall, then we only
evaluate for 15 sec/7872 = 0.00229 sec for the 2.5". This is multiplied by the Block Rub # 3
Tons/Sec (1.767e-8) to scale it down to the actual time the 2.5" of wall was rubbed. Thus now
this appears very minor.
**

**
E. PSI Application
**

**
As we loop up the building, actually building descending, it is better said that we are evaluating every 2.5"
as it collapses, so we examine every 2.5" section increments up the building when it arrives at the ground base.
As each 2.5" collapse, the weight of this 2.5" section is subtracted from the building weight for the next
2.5" evaluation. This reduces the PSI wall rub application. For the 1st 2.5", we calculate a PSI on the
bottom of the wall perimeter as 2070 (press) psi to the Block Rub # 3 = .5236 psi from the table above =
3953.42. Now the Block Rub # 4/# 3 Pressure Ratio from above = 2.38861 (prat0) so for prat0^ex = prat,
ex = log prat/log prat0 = 9.51215. From the above table, the Block Rub # 4/# 3 Grams Ratio is 6.01928,
so we have 6.01928^9.51215 = 2.601e7 (doratb).
**

**
F. Combining Factors C., D. and E.
**

**
Finally, we have 2.601e7 x Rub Area Ratios (22.50) x Block # 3 Tons/sec (1.767e-8) x Rub Time (0.00229 sec)
= .0237 tons. For the next n (2.5") of the loop, press = 2069.84 psi. All the other constants are the same,
only the new dorat = 6.01928^ex, finding ex as we did previously.
In this first case, .0237 tons is just the maximum potential of dust that could be produced, denoted by
the term tondb (tons of dust by building - wall rubbing on itself). Recall, in B. above, the tons/wall
per slab thickness height (2.5") is 9.383 tons (tpw). So, we can't rub off more than what's available
but we have .0237 << 9.383 so OK here. It doesn't appear that we would have a building where the potential
tons rub is greater than what was available if the wall thickness is large enough influencing the psi's
and tpw. In any case, we use the code check "if(tondb > tpw) tondb=tpw".
**

**
G. Budgeting Dust and Residue
**

**
We accumulate that tons of dust from the building source by the term 'cumtd'. Now the solid chunks that weren't turned
to dust is left as tonw = tpw-tondb. As we mention before, we split this to accumulate outside the wall and
inside the wall so each receives tonw/2. Redefine tonw now as original tonw/2. The inside wall
accumulating term is tons between slabs (tbs) incremented by tonw. It is reset to zero once a slab has fallen
on the pile.
**

**
H. Working With Volumes of Tonw
(1) Outside Wall Rubble Volume
**

**
We need to convert tonw to volume (tonwv) '^3 by dividing by density. We dump these solid chunks to the pile
outside the wall. Note in trap1, trap3
and trap4, that the rubble pile against the height of the accumulated
slabs (Hw). Now note in trap1, a vertical cross sectional mid length cut shape of the outside pile while
trap2 shows the top view. When we dump tonwv/2 on this pile, it is distributed all along the perimeter
2 Lr lengths and 2 Wr widths. However, since the perimeter wraps around at four 90 degree angles, what happens
at the corners? Here we view four 1/4 cone bases indicated by @. **

Thus we have the volume sum of 4 trapezoid prisms + 1 cone. For the cone, Vc=pi*Rc^2*Hw/3 and for the trapezoid Vt=4(1/2Rc*Hw)*2*(length+width). So total Vol = Vt+Vc = pi*Rc^2*Hw/3 + 4(1/2Rc*Hw)*2*(length+width) = pi*Rc^2*Hw/3 + 4*Rc*Hw*(length+width). Multiplying by 3, 3*vol=pi*Rc^2*Hw + 12*Rc*Hw*(length+width). We finally have a quadratic equation Rc^2*Hw*pi + Rc*12*Hw*(length+width) - 3Vol = 0. Solving for a quadratic ax^2 + bx + c= 0, x = -b +- sqrt(b^2-4ac)/2a. So here we have if x=Rc, a=Hw*pi, b=12*Hw*(length+width) and c=-3Vol, we solve for the positive root Rc = (-b + sqrt(b*b-4*a*c))/(2*a).

**
(2) Inside Wall Rubble Volume
**

**
Before we can apply Hw in (1) above here, we need to know the height Hw of the rubble from inside the wall.
We could define Hw as the accumulation of slab thicknesses if all the inside wall rubble (tbs) is turned to dust.
This dust has to be calculated. For the loop, we start out (n=0), and accumulate the 1st floor slab. Also, we
later make a cone pile out of tbs but initially, we have to set a tiny tbs=.1 ton so the math routines won't have
to work with 0's and terminate the program. Assume, we have a few floors collapsed already and we have some outside
wall height Hw. Now we accumulate the wall rubbing tonw to the tbs volume or tbsv '^3 by dividing by density.
Referring to trap4, we want to find Ht and Rb of the triangle rubble against the collapsing inside wall. This
triangle shape extends all around the inside perimeter of the building wall. So a cross section area cut (cac) =
tbsv/2(length+width).
By examining crushed stone piles (1" diameter), we seem to have a consistent Rb/Ht, or the inside cone Ht, has the ratio
near 9/10 call it tana, the tangent of it's angle of 42.0 degrees looking down from it's peak. Actually, these chunks are probably greater than 2-3" but the roughness of the
cement will provide more friction with one another to produce a possibly steeper height.
Actually, we tried to determine this angle with 3 types of sand. With a moist sand, we determined this angle to be 45 degrees.
With a wet sand, it first flows like a liquid and makes a pancake type base. As the water keeps dripping, the dryer sand makes
more consistency and the height builds up. Actually, this type of sand makes those high peaked dripping sand castles. We tried
to keep a consistent slope angle so settled with a peak (smallest angle) of 35 degrees. For dry sand,
the grains don't stick like being moist but flow like liquid wet sand in the beginning. Actually, this was controllable and
stabilized with an angle of 59 degrees. Then, continued pouring the dry sand until the dimensions where increased by over 2.5 times,
the angle was stabilized at 57 degrees. So averaging, we have 58 degrees dry, 45 degrees moist, and 35 degrees wet. So the stone pile
angle of 42 degrees is very similar to the moist sand.
This is all due to the material trying to pour like a liquid but
the increased horizontal friction caused by adhering of multiple sand grains when semi moist or the rough irregular
shape of stone particles digging vertically, causing higher psi in small columns inhibiting horizontal movement of adjacent stones.
Looking from the base, the slope upward is the 90 degree compliment or 32, 45 and 55 degrees for dry -> wet sand and 48 degrees for
crushed stone. Visualize this base angle (Ht/Rb) instead of looking down from the angle of the peak then dividing into halves, indicates
the material stacking is better. So generally, the slope is 1 near 45 degrees (Ht=Rb)
for most material, depending on it's roughness for larger chunks and less than 45 degrees for smoother dry material.
We want this cross section triangle to be in this ratio. So area of triangle (top view) = Rb/2 x Ht = cac.
Substitute Rb = Ht*tana, we have cac=Ht*tana/2 x Ht = Ht^2*tana/2, Ht^2=2cac/tana, Ht=sqrt(2cac/tana).
**

**
(A) Crushing the Inside Pile Theory
**

**
We use the energy of the slab (sle) to crush this inside pile of rubble. The weight on the slab is the weight
remaining from the building on this slab which travels a distance Ht for crushing. No work is done until the
slab is in contact with the pile. Thus sle is in ft-lbs, converting to kwh by multiplying by 3.76616e-7 kwh/ft-lb.
Now recall in WTC Dustbowl Model,
we had a crushp=1.15 KWH/Ton converting all the particles to the particle distribution
profile produced from the block rubs and considering 1/2 this value with cement on cement instead of metal on cement.
So potential tons by crushing is then tonbcp=sle/crushp. Is this crushp in KWH/Ton based on collecting the amount
of dust during some KWH of electricity by noting meter readings? If this is the case, is the machine efficiency
taken into account. With just a wild guess, for now, we might use me=.25, a 75% loss of efficiency by friction.
A Mr. John Dochner from Rohrer's Quarry in Lititz, Pa. has forwarded me some limestone crushing parameters
seen in **

**
(B) Limits to Crushed Dust
**

**
We can't produce more dust than tonbcp and we have the pile (tbs) and the slab (tps) to work with. We set a variable
ton = to the minimum of tbs and tps. With ton set to the minimum of slabs or pile, the crushing is effecting both pile
and slab so as with the rub, dust is being generated from both so the maximum dust will be double the minimum between
them or 2*ton. We may have to reduce tonbcp from it's potential to an actual by the following code. We have
if(tonbcp > 2*ton) tonbcp=2*ton, bringing the max potential tonbcp down to the actual 2*ton. As an example, at the
onset of the 2nd slab, tonbcp=945.7 with tbs=276.1 and tps=517.7. The min ton=tbs=276.1 which will usually be in all cases.
Now 2*ton = 552.2 < potential (945.7) so we reduced the potential tonbcp down to 552.2 tons of dust produced by the slab crush power.
Later, when there is less weight remaining from the building, the potential (tonbcp) is less than 2*ton. In this case, the actual tonbcp
is limited to it's original potential (sle/crushp). This happens after 46 stories collapse so this procedure is crucial.
**

**
(C) Budgeting Dust and Residual
**

**
Again, we immediately accumulate this crushed dust to variable cumtd, where the wall dust accumulates also.
Now with the residual chunks, we add it to the other pile inside the wall (cumti). How much debris
do we add? Well we took 1/2 from each pile and slab source so what remains is
(tps-tonbcp/2)+(tbs-tonbcp/2) or tps+tbs-tonbcp. The first expression explains what is really
happening where the 2nd is mathematically equivalent and easier - or could say, add slab and pile
then take away the dust they produced but it doesn't happen in that order - we really add what is
left in the slab after this dust and also what is left in the pile after an equal amount of dust (276.1 tons)
or 517.7-276.1=241.6 tons left in the slab and none in the pile between slabs in our example. We then reset the tons between slab (tbs) to 0 for the next slab.
**

**
(3) Outside Wall Rubble Dust Process
**

**
Procedure and values are different than the rub wall application.
**

**
(A) Changing Rub Areas
**

**
Now that we have a wall of debris from (C) above, we have it's height Hw, where we apply the quadratic
roots formula in H.(1) above. Here we have the base or cone radius Rc in trap4 by solving the quadratic.
With Hw and Rc known, we have L = sqrt(Rc^2 + Hw^2), the length the residual tonw slides
down the side of the outside pile (Airborne Dust # 2). The lateral surface area of the outside rubble trapezoid is L x perimeter
or L*2*(length+width). Now L is also set to the lateral length of the four quarter cones whose total lateral surface
area is the 1 cone pi x Rc x L. Both these areas are summed to represent the rub working area wca. Now just
as we had in 4.C. Rub Areas above for the wall, we now have the debris trapezoid - cone area to the rub3 area
as rarc=wca/rub3a.
**

**
(B) Changing Rub Time Factors
**

**
With Hw and Rc, we also have the 1/2 angle of the cone peak trap-cone as angc=Arc Tan(Rc/Hw). For doratc,
the psi component normal towards the side of the trap-cone is sin(angc) so prat is reduced, so is ex and
the resulting doratc. However, the velocity component may be slightly faster because of less friction from
less pressure so we use a factor of doratb/doratc to approximate the speed up the velocity. However, the
velocity component is now avs*cos(angc) where avs is the average speed of fall @ bh/sec from Actual/Gravity
Rate of Collapse = 1.6284 from table above = 91.0'/sec or 62 mph done above in first paragraph of 3.
Dust by Physical Processes:. So for the travel time secc along L, secc=L/(avs*cos(angc)*doratb/doratc).
**

**
(C) Combining All Factors
**

**
As with the wall rub above (4.F. Combining Factors), we still have the same physical terms but different values
so different variable names using c's as in tondc=rubtps3*doratc*rarc*secc.
**

**
(D) Limits to Rubbed Dust
**

**
Again, we are careful to account for all the tons - dust and pile residual. This is similar to the dust caused
by crushing in 4.H.(2)(C) above but important differences. Here the potential max dust is represented by tondc.
Again, we set ton to the minimum of the two rub sources this time, tonw again - the residual from the wall dust rub
and cumto - the cumulative tons of outside rubble against the wall. Recall, tonw also went to inside the wall to be
crushed by the slab, this tonw will be rubbed by sliding down the trap-cone lateral side or rubbed the 2nd time.
Again we have "if(tondc > 2*ton) tondc=2*ton" bringing the max potential tonbdc down to the actual 2*ton.
**

**
(E) Budgeting Dust and Residual
**

**
Now here, we want to add what residual is left to the outside wall cumto pile. This residual is tonw-tondc/2.
The other half (tondc/2) of the dust was from the cumto pile itself. So we have another tondc/2
to subtract from cumto. Altogether, we add (tonw-tondc/2) - tondc/2 or tonw-tondc. Again, the simpler
last expression conceals the actual process. We add tonw less half the dust that was rubbed off from
it to the pile, which also contributed half the dust. The difference here is that the pile that we add
the residual to is growing and provides a larger rubbing surface area itself rather than in the crushing, where
we crush the whole pile between slabs each time. Also, we add the generated dust tondc to total dust term cumtd.
**

**
5. Results:
**

**
Reviewing the bottom of ddmdat, we see in Rubble:, we listed the tons of residual from within and outside
the wall and total dust from all sources combined (cumtd). This should total the tons of the original building
before it was demolished. Corresponding and underneath is the % of each contribution of tons. It is seen that
for the WTC, this method produced 62.6% to dust compared 72.5% in the wtcdbm model.
This lower value is interesting but acceptable. In the previous, we took a % of the total cone volumes for each foot of building fall.
Since we used the resulting volume to help predict the total volume, there were some bias built in. However,
there was some merit to quantify the physical process that we could use here. In this study, we used actual psi on
walls when they occurred and treated the slabs separately as a crushing estimate. The pile and the rubble were
accumulated but had no bias built in. We will be examining some variables more accurately as mentioned above
in 4. Deriving Other Parameters from the Building Parameters:. **

**
6. Adding Steel & Fluff to Pile for Final Dimensions: (Not Used in Final Model)
**

**
We will add steel (vols) to the inside of wall only as that is most likely. Now vols=psv*volc where psv is
Steel/Cement by % Vol in table above. Thus, the new inside volume iwv=cumti*2.e3/den + vols. Now the new
hw=iwv/iwsa where iwsa = inside wall slab area as used above. This means we have a new height for the total
complex outside wall volume owv. Again, we apply the quadratic root equation in 4. H. (1) above for the
radius/base Rb but substitute owv for Vol.
For the fluff factor, let's just boost owv directly because of the crushing slabs, more air spaces would be
jammed up and filled. Also, since the WTC pile was maxed at 3 stories high, lets not increase that anymore. So
we solved the quadratic root equation for Rb again and substitute owv for Vol again. Only the Rb extends
outward but still short of the observed. The real wide umbrella cone in the WTC really brings up questions and
a mystery when it should have just fallen vertically. Could it be that the added energy of the plane and fuel
burn, primed the WTC for greater efficiency for when the gravity took over? It should appear that the fall
velocity would be greater but the extreme lateral spreading indicates some other energy assistance was supplied
horizontally internally?
**

**
7. Particle % Distribution from Crushing:
**

**
Above, we calculated % dust but assumed the particle distribution
produced was from the rub blocks @ 12,800 dpi plot. This time, the dust produced
from inside the wall by crushing will have a different particle distribution.
In rubbing, the finest individual particles are removed and become airborne. The
rubbing of larger particles breaking down to smaller is assumed minimal. In the
crushing process, larger particles are breaking down to smaller particles, a different
process than rubbing.
From the quarry, we refer to **

**
A. Split Rub and Crush Outputs in Program
**

**
We now modify the program DUSTBM.c (DUST Box Model), to output separate files of total dust from the outside wall
rub and inside wall crush components. One output file is the total dust as before. Every time
a slab crushes the inside wall residual pile, another output file contains variables rcaf and tonbcp. To refer to these terms, we view the diagram and explanation in crush.dia.
Here we assume a slab travels 2 thicknesses of itself then breaks up. Before breaking up, the slab is attached to
the building taking all the tons weight above it (ctr). So for this distance (2*slabt), x ctr*2e3 = ft-lb =
sle (slab energy). For the rest of the distance, (hpi-2*slabt), we have the remaining tons/slab (tps) weight applied
which is rat = 4*stl*(length+width)/iwsa where stl=cac/slabt, the cross sectional area cut of the pile (cac) divided
by a new height of 1 slab thickness.
That is because the rubble triangle base is spread by being squashed like a pancake with height =
1 slab, and length stl (slab thickness length). The 4*(length+width) is the perimeter as
used before so 4*stl*(length+width) is the area along the 4 sided building perimeter.
This area to the inside wall area (iwsa)=(length-wallt)*(width-wallt), is the ratio rat.
This detached energy contribution is minimal compared to building weight original sle but
added to it by sle = sle + (hpi-2*slabt)*rat*tps. Now the total sle is converted to KWH
by sle/ft-lb*3.76616e-7 kwh/ft-lb. Now then tonbcp=sle/crushp where crushp=.004315 KWH/Ton
is Rohrer's Quarry Vertical Shaft Crusher from table QuarryCrush.dat
and tonbcp = Tons of dust by crush power available. Actually, this is the potential (poton) for that
production depending how much cement is also available as explained in 4.H.(2)(b)
In this update, we calculate the remaining tons per slab (tps) or rat*tps now instead of just
tps because of the breakup of the slab crushing the residual pile.
As before, we choose the minimum (ton) of the two (rat*tps or tbs= tons of rubble between
slabs), then each (crusher and pile) supply ton or 2 ton total. Now before we used
if(tonbcp > 2*ton) tonbcp = 2*ton which is still good but also, how many times can the
actual tonbcp be crushed if available? Well that simply is rcaf=poton/tonbcp but of what use is this?
Now we have how many rubbles this size could we crush to produce the quarry's dust
profile as seen in the graph sievecpd.gif. Since there is only one rubble
pile here, we use the remaining energy sle to keep crushing this pile to smaller particle sizes
because we ran out available cement to crush.
**

**
B. Keep on Crushing - CRUSH%D.c Update 10/22/08
**

**
In another program, we read the Quarry's Dust Profile values interpolated from the graph
sievecpd.gif as we did before with a few programs in
World Trade Center Dust Bowl Analysis under 2. % Particle Distribution:.
We call this file "crush.%d" and will be used in all subsequent programs.
Suppose we want to extend the crushing beyond this profile to produce a new particle % distribution
profile reduced to smaller particles. Well we could determine how much energy was available and
keep applying the sqrt(d/id) for each d until the crushing energy is matched. Rather than doing
this over and over, lets determine this beforehand by crushing this original file crush.%d in
increments to produce a different profile for each corresponding energy crush increment.
We determine how much crushing energy is produced from each slab and find it's corresponding
calculated compressed energy - % profile increment. However, if crushed less than that
required to produce crush.%d, we have larger particles, rather than smaller. We won't extrapolate
this backwards here but remember that more larger dust particles could be available with less crushing.
In the expression crush=sqrt(60/d), lets change the constant 60 to any diameter d and d to a
new standard constant id. To crush from d microns (u) to id microns (u), requires crush=sqrt(d/id).
Or in other words to increase the energy factor of sqrt(60/d), we divide by sqrt(60/id).
That is sqrt(60/d)/sqrt(60/id) would also be sqrt[60/id)/(60/d)] = sqrt(60/60 * d/id) = sqrt(d/id) again.
In this program, we read the quarry's crush.%d for d and % particle distribution pd. We start and reduce
all d's to 1u by crush=sqrt(d/id), then times pd (weighted) then cumulated by ccrush = ccrush + pd*crush.
Since the sum of pd's = 100%, then the sum of crush's (ccrush) is the final weighted energy required
to crush all particles in crush.%d down to id=1u. We won't take all the particles to that extreme and may not
have that required energy. Also, the WTC photos of 911 show dust particles settling, where 1u particles would
essentially stay afloat in the observed downwind distance.
Starting from 2"-3" cement chunks, we have crushed to a profile as seen in
sievecpd.gif again, it appears as a double normal curve.
Other normal curves are seen in NormDis.png. Here
u is the mean=median=mode at u=0, and u=-2. The standard deviation is sigma
but shows the plots of various sig^2's. The formula for a normal curve is given
by pdf.png. A sig represents one standard deviation of the
data. Now the area of the peak to 1 sig on either side of any normal curve
represents 68.26% of the data. If the area under the whole curve = 1.0, then -1 sig
to + 1 sig = .6826. The exp( ) part formula above can be rewritten as exp[-(x-u)/sig]^2/2
where (x-u)/sig is the number of sig units we call (ri). So if x=10, u=4 and sig=2, then
ri=(10-4)/2 = 3 sigs. The area of the curve between + - 3 sig is 99.74%. Actually,
most use 2.5 sigs to represent all the data (98.76%). For our case, we require not to
lose any dust so we choose + - 5 sigs = 100.00%, 100.000% or more decimal placed zeros.
What do we plan to do with this normal curve? We want to break particles in half, or try
to, that even though most won't be exactly half, determined by what precision you judge by,
they will distribute about the mean=mode=median of 1/2 the size to crush. If we consider a
particle size of 504u, which we won't use that large in our model, then 1/2
that or 252u will be the mean=mode=median as indicated by pb2.gif.
Here we don't have the sig but estimating 68%, we might eyeball guess sig = 65u to get
the feel of a normal distribution. If we want ri = 5 sigs to represent all the data, then
sig must be 252u/5 = 50.4 u, less than our eyeball estimate.
In our program, let's use 1u intervals, determine the size to crush (max), and
ri= -5 sig -> 0 sig -> +5 sig or 0 u -> max/2 u -> max u or 0 u -> 252 u -> 504u.
So we have sr2pi=sqrt(2*pi) and sig=(max/2)/5, and ri= (max/2 - i)/sig where
i=particle diameter index [i] and u's. We loop i=0 -> max-1 and ri2=ri^2. We have
y or erf[i]=exp(-ri2/2)/(sig*sr2pi). Two examples are seen in erf.dat.
The simple 1st case is sig=1u centered around diameter = 6u. Notice the symmetry of the
erf x 100 values and the cumulative area of 100.0%. In the 2nd example, we still have the
cumulative area of 100.0%, symmetry around 26 u but many more values with sig= 5.00 making a
wider bell curve similar to the blue plot in NormDis.png.
We decide to use diameter (max). What decides the size max? Well, we could just start with the
max diameter size, crush a bit, then the next largest size down to 10u or less. The most likely of
sequence particles to be crushed will be all the volumes occupied by particle diameter i. Simply,
this is the % particle distribution pd x volume (vol) where vol=4/3*pi*r^3. However, vol x npar
(%) is also termed particle radio activity so pact[i]=pd*vol for each particle i to np. For our 1st
max, we read diameter or index i & pd from "crush.%d" and calculate vol x np, for pact[i].
For our 1st case, max pact[i] (or pd*vol) was np=346u.
Here mean=mode=median max/2 = 173u. So the 1st crushing will be on the np max = 346u particles. How much of this
volume should we crush? If we crush all of it, that is too much and we might show a rough % particle distribution. Too little,
and the program and data output files will become too large. We can take a percentage, so we did and the best
choice is 20%. Why? Well with the 6000 loops, the np=346u repeated about 5 times as the max? This is what probably
happens in the real world. Just 20% of the volume is taken off the max particle volume and then after 2 or 3 times, there was
another size to replace 346u as max. Then 346u comes back a few more times as max again until most of the 346u particles have
been depleted by being crushed smaller. So even though vol is small with smaller particles, the np has accumulated immensely
and hardly any more larger particles are to be crushed.
We find the total volume occupied by 346u's as vols=tvol*pd, correcting pd(%)/100 where tvol is total crushed dust volume.
Lets use delta volume (dvol) of vols. For dvol, we use a factor 20% = .2. So we have dvol=.2*vols. We take this dvol value
from vols. What we do now is to reduce the % particle distribution of 346u by pd[346]=(vols-dvol)/tvol then * 100 for %.
Next we distribute dvol by the normal curve's erf values as just seen by the two examples in erf.dat.
To see this probability or % distribution normal (pdn) applied to all particles < 346u, see crush.346.
Note, at the bottom, it is subtracted from 346u, while for the other particles it is added. The added effect on pd x 100 appears nil
sometimes showing no change within the 4 decimal % near the tails. However, as you approach the max/2 = 173u, it's added effect (->)
by the erf's pdn is the greatest and diminishes towards the beginning (1u) and end (345u) of the profile tails.
We call this the crushing loops l. This is repeated 6000 times so there are few larger particles to be crushed into smaller particles.
If we go more than 6000 crush increments, they jam into the 1-10 u range and becomes a blocky distribution. Also, 1-10u particles are
considered to stay aloft for many miles so no need to distinguish between these small sizes and neither did the earlier models mentioned.
For 6000 crushes, we end up with the profile as seen in loop.6e3. Here we see the jamming tendency
already in the 1 - 10u range (cpd=90%) and we have our cumulated % distribution of 99.9% within 1 - 100u range.
For each of the 6000 crushes, we calculate the crush energy required for each new change of % particle distribution.
So we have crush=sqrt(max/id) for each diameter max crushed to diameter id as explained above where
diameter max replaces 60u. We have id in the denominator because we crush to the small sizes being in
the denominator < numerator so crush > 1. Now we accumulate crush by ccrush. However, crush is really
weighted by the number of particles np or pd, really pd/100 since pd is given by %'s. So the cumulated
crush would be ccrush = ccrush + crush*(pd/100), but we're applying the % of the normal erf curve (pdn)
representing how the crushed particles are distributed. Here, pdn @ max/2 would be the highest.
So now, ccrush = ccrush + crush*(pd/100)*pdn when pdn is in decimal form and pd in %. For each of the 6000 cumulative
crush intervals (ccrush), we have a unique particle % distribution. Both the ccrush and associated new % (pd) are written
as a 3 character binary file with ccrush for each of the 449 % pd. This file (cumcaf.pd) is over 8 mb in size.
**

**
C. Matching Crush Energy Between Files
**

**
So we read the total dust from one file, then the WTC rcaf,tonbcp - dust produced by each slab tonbcp (cton)
from another file. Now real caf (rcaf) above is defined as poton/tonbcp, or a factor of how many times tonbcp could be
crushed to the quarry profile % dist (crush.%d). Now in review from above, poton (potential tons) = sle/crushp = tonbcp from
ton (lesser of tbs=tons between slabs, tps=tons/slab). if(tonbcp > 2*ton) tonbcp=2*ton. It's 2*ton because tons from slab
& pile are crushing each other equally. In topic above, we just calculated the cumulated crushing adjustment factor (ccrush).
So now we match these energy factors, after we rewrite the WTC rcaf,tonbcp file in ascending order ( crushr.wtc)
to match the ascending order of the 8 mb cumcaf.pd file. We read the crushr.wtc file first, then ccrush in the cumcaf.pd
- %pd file. As soon as we match the real rcaf by the code if(ccrush > rcaf), we continue the id loop to read the np particles % distribution (pd). An example of this is seen in
slabs.dat. The search for the ccrush > slab's rcaf is given for the 1st and last 11 slabs. Note, we start the loop
1 back from the previous matched slab in the loop.
Now, we matched ccrush with slab's rcaf so we have a particles % distribution stored for each slab in an
array of particles dimension pact[id] = pact[id] + pd*tonbcp. We used pact[] again, which originally represented particle activity
but now used here as tons of dust for each particle size diameter
accumulated for each slab. Actually, tons - volumes can be used interchangeably by a constant of density so it actually is still
particle activity pact[] again but expressed in terms of tons this time but really means ton[i]. As a check, we now totaled the
final pact[]'s, from all collapsed slabs, and the crushed tons cton = cton + tonbcp from the crushr.wtc file = Crush Tons = 59455.6.
The sum of pact[i]'s was 59525.0 or .12% too high. With rounding off and accuracy, this is excellent agreement, however
all pact[i]'s were scaled by 1/1.0012 = .9988 correction factor anyway.
**

**
D. Recombining Rub and Crush Dust % Dist
**

**
Finally, we have the total crushed tonbcp's (cton) as dust from slabs. We also read the total tons of dust (tton) from
a file earlier, then rton=tton-cton is the total rubbed tons of dust. From our previous model runs, we have
Total Tons 73788.063 = Crush Tons 59455.637 + Rub Tons 14332.426.
In the rubbed % part dist file (rub.%d),
the pd read after each particle diameter i produces rubt=pd/100*rton of dust. This added to the
crushed dust ton[i] or denoted by pact[i], the new % of dust for each particle i is dustp=(pact[i]+rubt)/tton*100. No need to
work with an array [i] any more, so i and dustp are dumped to the final dust output file
dust.%d. However, we did not have any rubbed particles > 399u, so
the rest of the tons/particle i diameter is the crushed tons only from 400 -> 449u.
To see a comparison of rubbed and crushed tons tabulated for each particle i = id diameter, see
crush.rub. Note, the crushed tons < rub tons from 1-13u then
crushed tons > rub tons from >= 14u.
Using dust.%d in the DUSTBM.c again, and removing the fluff factor because of rethinking the packing of
tiny particles layered to a few inches, compared to pulverized limestone, we have the updated output file
DUSTBM.dat. The rubble pile has not changed nor the amount of dust but because
we have overall larger particles in crushing than rubbing. As seen in crush.dia, we
have a best fit with the observed WTC street level dust, so we repeat here, experimenting, we found that
3.6 x slab thickness rather than 2, and the occurrences of no lateral movement is at 1/1.4 times the slab
travel distance of 3.6 times it's thickness. So 1/3 --> 1/1.4 and 2 --> 3.6 correlates better with the
observed WTC street dust photos from the web. So finally for the WTC's,
The max dust was 2.57" @ .147 mi and .45" @ 2 mi with 52.5% of dust deposited within the 2 mile range.
62.9% of WTC buildings went to dust. - End Update 10/22/08
**

**
8. Kinetic Energy Effects of a Nuclear Blast:
**

**
**

**
A. Pressure Wave
**

**
After viewing 1 Megaton PSI Effects, we chose
within a radius of 5 psi from a any nuclear blast, to collapse all buildings, at least the cement. All the calculated dust
will be combined with the ground crater dust from a surface blast. Each component has it's own particle distribution to be
inhaled by the updraft velocity of the stem in the nuclear mushroom cloud.
We plotted each of their radial distance x vs. psi. For quick reference, it was as follows:
**

R mi | PSI |
---|---|

1.7 | 12 |

2.7 | 5 |

4.7 | 2 |

7.4 | 1 |

**
**

**
B. Crater Dimensions
**

**
For a surface burst of 1 megaton, the PBS website indicates a hole of 200' deep, 1000' wide.
This is seen but inverted in crater.dia. The crater is the top segment
of the sphere (hole in the ground) with a radius of curvature R. The width of 1000' is
broken into the smaller radii r, which is the radius of the hole at the ground surface - remember, this is inverted.
This dirt segment volume of the sphere is also called the cap of a sphere.
We can compute this dirt volume by either of the 2 formulas shown v=pi/3*h^2*(3*R-h) that requires the radius of curvature
R=(r^2+h^2)/2h or v=pi/6*h*(3*r^2 + h^2) using r and h directly. The radius of curvature here has it's center 725' above the
bottom of the hole or 525' above the ground. In both cases, v = 8.273e7'^3 and the volume V of the whole sphere, of which
the segment is a part is, V = 1.596e9'^3. This means the dirt blown out the hole is only 5.2% of the energy used as if the
bomb were implanted below the surface deep enough before the explosion, but near the surface enough to just blow the hole out
in the ground above it. This 5.2% will hold for all bomb sizes so we simply multiply v = 8.273e7'^3/mega x megatons for any
other surface burst bombs. However, we don't have the particle size % distribution yet. Now we could take dirt, shake it up and
blow it onto the plate glass for scanning as we did before with the limestone from the quarry, however, this is certainly what
a nuclear blast does. We may not have to do this as we explain later.
**

**
9. Extrapolate WTC Dust Generation to Major Cities: - CITY%D.c
**

**
**

**
**

**
A. City Building Profiles
**

**
We now have the amount of dust for each particle size from 1-449u from the combination of rubbing and crushing. We will apply
this method to different cities as we estimate each city's building density distribution from photos. If any reader has done
this with better information, please contribute to this model development as we will acknowledge all those that contribute.
All the buildings within the 5 psi ring will use some of the other building parameters from the table above rather than the WTC's data.
As a first case, since we did WTC, we might as well start with New York City (NYC). Looking at a map,
we classify Manhattan. We will have approximate city profile by building estimates on each city in a file called city.dat.
Such a line of data is as follows: NYC cityl= 9 cityw= 2 cityd= 65 floors1= 65 floors2= 25 length= 240 width= 225.
Now cityl=city length in mi, cityw=city width in mi, cityd = city building density in %, floors1 = highest
average # of floors (stories) in center of city, floors2 = average # of floors (stories) at end of city (cityl),
length = average building length in ft throughout city, and width = average building width in ft throughout city.
Past the length of the city, we have floors=4 stories.
However, using this city's building height profile and size, we need to estimate the height for every .1 mile from ground
zero to the city's outer limit shown by the comma symbol with city's equivalent area radius Rc in citydust.dia. So for each dr=.1 mi, we round off to the nearest
number of buildings floors, rather than the building height because the distance between floors is assumed to be divisible by the
slab thickness. For the WTC, it was 12.4 ft, and rule of thumb is 10' floor-ceiling. We will center between 10' & 12.4' to near 11' including the slab.
We can run the model and observe the effects of different slab thicknesses. For whole # slabs, nsc=11/slabt+.5. Then the actual ceiling height ch=slabt*nsc.
We loop up the building by the slab thicknesses, so the number of times to loop n is ln=nsc*floors. With every nsc, mod = remainder
n/nsc = 0, we evaluate the slab crushing routine. We linearly interpolate and round off the # of floors between floors1 -> floors2
to distance x5psi. Any distance past Rc, floors=4 stories. We have lx=Rc/dr+1 buildings to evaluate along the radial, since we start at center of city
ix=0, to ix<=lx at x5psi.
**

**
B. How Many Building Fall?
**

**
There are 3 possible cases where we must examine the 5 psi ring for city and metro-suburbs. These are presented in
citydust.dia. In the top diagram, this is the most likely case, the 5 psi ring is greater than
the city limits and encloses it. All the city buildings (#) are gone, as well as some of the suburbs (^) buildings out to
the 5 psi ring. In the 2nd case, we reverse the 1st case, the city limits are bigger than the 5 psi ring. The 5 psi ring
never reaches to the end of the city buildings or any of the suburbs. Many city buildings and all the suburbs buildings
are still standing. In both of these cases, we can use the city building density cityd without modification because
all the city lies within the 5 psi radius or the 5 psi radius is a circle that includes only a part of the city.
Our 3rd case is complex in that the 5 psi ring criss crosses the city's boundary. If we use an equivalent city radius
Rc, this won't criss cross but will be either in or out of the circle boundary as either of the 2 cases above. This situation is
likened to NYC's Manhattan's Island, where we have a 2 mi x 9 mi dimension. It will be best with such an elongated dimension to perform an
alternative method. By the diagram, we see in reality, the length of the city is not reached by the 5 psi ring so that
forces the buildings in the city's length into it's width which over estimates because for NYC, much of it is water, and
should not be included. If the city's dimensions are almost a square, use of the equivalent radius Rc is fine.
We need to know, the city area inside the 5 psi ring. So by integration, we first use the equation of a circle and solve
for y in x^2 + y^2 = a^2, y = (a^2 - x^2)^1/2 where a = radius 5 psi ring (x5psi). Now the x
limits are from x1=0, to x2=cityw/2 for right upper quad the of circle. So INT y dx =
[x*sqrt(a^2 - x^2) + a^2*asin(x/a)]/2, remove /2 and we have right semicircle since % in area of
quad = % in area of semicircle and circle. Substituting limits x1 and x2 for x, we have NYC's 5psi
area mi^2 cpsi5a = x2*sqrt(a^2 - x2^2) + a^2*asin(x2/a) - (x1*sqrt(a^2 - x1^2) + a^2*asin(x1/a)).
Subbing with x1=0, we are left with cpsi5a = x2*sqrt(a^2 - x2^2) + a^2*asin(x2/a).
Now the area of the 5 psi ring is psi5a=pi*x5psi^2 so the ratio of the integral area to it is
rat=cpsi5a/(psi5a/2), using /2 for a semicircle. Now we can reduce the cityd to rat*cityd.
Back to the 5 psi ring (psi5a), we convert mi^2 to feet^2
and divide by the area of the building (width*length), we would cram all the buildings into the 5
psi ring or 100% and have the max # of buildings (nbm). When we multiply this, nba = nbm*cityd/100,
we now have the actual # of buildings (nba) within the 5 psi ring. In the each building model run loop,
we have the tton of dust but that was for lx = 28 buildings along the radial. So the average for 1 building is tton/lx. Now times
nba actual buildings, we have total tons of dust for the city. However, this is only a fair estimate.
**

**
**

**
(1). Weighting by Concentric Rings Areas - Update 11/20/08
**

**
We use dr=.1 mi out to the 5 psi ring, we have 28 buildings dr's. There is a building dust calculation for each dr.
For any 2 radii, the area of a ring is the difference between their areas or pi (r2^2-r1^2) or pi (r2+r1) (r2-r1).
If dr is the thickness of a ring, and r is the inner radius, then r1=r, r2=r+dr. So ring area = pi (r+dr+r) (r+dr-r) =
pi (2r+dr) dr. Or if r2=r and r1=r-dr, area = pi (r2+r1) (r2-r1) = pi (r+r-dr) (r-(r-dr) = pi (2r-dr) dr. In the 1st case, r can
start the loop at r=0, but not in the 2nd case so we would start at r2 = dr, then r1=0 for the start of the loop.
We actually took this case as we didn't want to start with r=0 since it may have caused problems later in using r later in the program loops.
For the 28 dr ring buildings, we see the results in city%.dm. Here we see the Dust Ring Area
and it's accumulated columns. The decimal or Fraction column is the cumulated areas ratio to the 5psi area (psi5a)
up to 1.000.
**

**
(2). Another Look at PSI's
**

**
After some subsequent runs of the later programs, to cut off the city's building dust at the 5 psi ring seemed too short and abrupt.
In the real world it, would rather be a diminishing of dust . Reviewing
1 Megaton PSI Effects again,
less buildings will fall at 2psi so let's cut off any buildings falling at 1 psi instead.
Using our formula above to work with in the programs, instead of the table of points, we use
x5psi=2.76, x2psi=4.69 and x1psi=7.00 mi. Now for our new radial distance r, how do we reduce the # of buildings falling
between x5psi and x1psi? We could linear interpolate 100% at x5psi and 0% at x1psi. However above, we saw it as
a power reduction formula as seen in rpsi.gif as
psi=mega*28.5 mi^-1.653 or a slightly better fit as psi=mega*29 mi^-1.73.
Here mega=1 and we just solved for mi for x5psi,x2psi and x1psi above. We introduce a new ratio, pr
as the Power Ratio or PSI Ratio as pr=(x2psi/x5psi)^-1.73 = 0.4000 as an example. So the ratio
of these distances is still that inverse power law, but not quite the inverse power law squared (-2).
Subbing for the rest for mi, x1psi/x2psi = 0.5000 and x1psi/x5psi = 0.2000.
If just the ratios of distances, we have 2.76/4.69=.588 instead of .4,
4.69/7=.67 instead of .5, and 2.76/7= .394 instead of .2 indicating the power law as a better and
faster reduction of the pressure wave with distance. Notice that the simpler form is
just the ratios of the psi's, psi2/psi5 = 2/5 = .4 yields same results. However, we
need the psi with each dr so pr = psi/5 = mega*29/5 mi^-1.73 where mi=r. We will use
the pr to represent the reduction factor for the # of buildings to fall, but for x=r>x1psi,
we have 0.0. We rather drop pr from .2 to .0 past x1psi rather than 1.0 to .0 past x5psi.
For r <= x5psi, we keep pr=1.
**

**
(3). A 70 Ring Circus
**

**
Where do we apply pr? Rather than put it off and apply it's correction in
some of the subsequent programs and maybe make a mistake or lose it's meaning
in the later complexity, we thought it best to apply it directly as a correction factor now
to the Dust Ring Area that is stored in an array called perimeter area pa[] for use
in the following major nested loop. Also, since we now use x1psi instead of x5psi, with dr=.1 mi out to the
1 psi ring (7.00), we now have 70 dr rings or 70 buildings. Now instead of working with the max # of
buildings (nbm) for the whole 1psi area, we compute nbm and adjust to nba (actual) by the city's density for each
ring area pa[ix] and now further adjust pa[ix] by the factor pr.
The dust components and new ring areas with pr applied is seen in city.%pd.
Notice here, that the Dust Ring Area = pa[ix] until pa[ix] is corrected by the pr < 1 factor.
Also notice the dust components and the crater dust volume comparison. Not only do we have a 70 concentric ring circus,
but we have more smaller performers in the larger outer rings.
**

**
C. Super Duper Looper I - CITY%PD.c
**

**
The 1st file to read is city.dat - the city's building profile.
The 2nd file to read is rub1.NYC which is the Rubbed Tons (rubt) for each of the 70 ix (dr=.1mi) ring
buildings from CITY%D.c . Since we have assumed, but not exactly, to distinguish between
the two dust generation methods, we have the same particle distributions from all rubbed slabs. That is,
more dust, not more smaller particles as crushing does. We open the 3rd file rub.%d
which is the Rubbed % Particle Distribution (pd) up to 399u. With these last two files, we initialize the tons of dust for each
particle size id as ton[id]=rubt*pd/100*nba*pa[ix]. The pd is the % distribution of the particle size id and the nba*pa[ix] is
the real # of buildings per ring [ix] as explained in the previous paragraph. Still reading the 1st ix ring data, we now open the
4th file crush1r.NYC, which is the dust potential from the
crushing slabs program CITY%D.c
reading each slab's rcaf (column 2) and tonbcp (column 3) from all 70 ix ring buildings (column 1).
In crush1r.NYC, we have an indicator record (1.000 0.001) when we are reading a new
stack of slabs in the next ix building ring. We now open the 5th file,
cumcaf.pd ,
the binary file from CRUSH%D.c
and run another nested loop to match 'ccrush' @ (id=0) with rcaf as explained above. Since there is no id=0, it's a perfect
place to identify 'ccrush'. With pd = particle size (id) decimal distribution for each slab in building ix, we now add this to the
rubbed tons id array ton[id] = ton[id] + tonbcp*pd/100*nba*pa[ix], which is the same procedure as for rubbed tons except tonbcp replaces rubt. **

As the first ring of building's (ix=0 -> .1mi) is complete, we write it out to ctydst1.NYC, as the last ix (lix) and the ton[id], and reset ton[id]=0 for the next ix. We then repeat all these loops by continuing to read the next ix ring and must rewind files rub.%d and cumcaf.pd again for particle size distributions. Finally, when all the buildings are read, we reread all the ix and tons from ctydst1.NYC for a check to see if it matches the totals before all the pd's and pa[ix]'s sub divided the total against the 6th file totald1.NYC . This difference was < .05% error so no errors in the loops or the summated pd factors were evident.

**
**

**
D. Super Duper Looper II CDRCRD.c & PD%R.c Update 03/08/12
**

**
Particle Distribution % by Regression replaces a few loops in CITY%PD.c but is now further complicated by
the list of discrete buildings with their random radii (ir) from the LAS format instead of generated consecutive looped building previously.
Since we have discrete buildings, we would like to store the % pd x tons 'tons' by radial distance 'ir' & particle diameter 'id' depicting a 3D surface,
for a final output file but is too large for our CPU ram here. We output to a file then read it by file positions generated by id,ir. Then we have to do
the same over again for the crushed dust after the rubbed dust. Most all large files in these RADFO programs use this procedure.
We finally approach this problem because Super Duper Looper I above was complex and was developed to combine all these files in
CITY%PD.c rather than all the effort in CITY%D.c instead. A two-program approach was to make this easier.
However, we are forced to a new method in CDRCRD.c because we now evaluate different building dimensions in the same psi ring (ir) and
the number of buildings create their own actual density spacing and height instead of estimating it.
The matching of all the files in
CITY%PD.c was done by sequencing through until a match was found. Matching 'rcaf' with 'crush'
could apply to 6000 records for every building vertically for every slab thickness. In the estimated building height profile, linear decreasing
with distance (ir), we could estimate where to start searching in the files since it is a few feet lower than the last building. We could also
reverse the file order and start with the smaller buildings and the next larger building's matches should be the next slab group into the file.
Here, any building height appears in any 'ir' ring. Therefore, we may have excess CPU times if we continue this method here.
Once the match is found, we extract all 1->449 particle id's pd's (% distribution) in the binary file's (displayed here as ASCII) seen to the right of id=0
column in cumcaf.pd which was calculated from CRUSH%D.c .
PD%R.c is just a modified CITY%D.c from above. Instead of incrementing
each ir psi ring and approximate a common building height vs. distance ir with fixed dimension for each building, then estimating the number of these buildings
in this ir psi ring by a density guess, we calculated each ir psi ring by vector range from each building to the city's ground zero building as well as their
real dimensions and even their orientation for psi impacts. All this will be developed later but now seen in 920515dr.lwh
as a later input file.
**

**
(1) New Files & Order
**

**
Now the new type of city data (1st file) to read is Miami.dat - the Ground Zero's building name & location
by 3 different coordinate systems (Lat/Lon, UTM & State Plane).
The 2nd file to read is rub1.MIA which is the Rubbed Tons (rubt) for each of the 70-ir (dr=.1mi) ring
buildings from CDRCRD.c . Since we have assumed, but not exactly, to distinguish between
the two dust generation methods, we have the same particle distributions from all wall-pile rubs. That is, more dust, not more smaller
particles as crushing does. We open the 3rd file rub.%d which is the Rubbed % Particle Distribution (pd)
With these last two rub files (tons/ir & % dist/id), we have the % tons of dust for each particle size id up to 399u for each
given radial increment 'ir' to complete our 3D surface area.
Still reading the 1st 'ir' ring data, we now open the revised 4th file loop crush1.MIA which is the dust potential
from crushing slabs program CDRCRD.c
reading each slab's rcaf (column 2) and tonbcp (column 3) from all 70-ir ring buildings (column 1). Notice each building's 'ir' is
in random order but remains the same until vertically profiled by slab thicknesses.
We now open up that same 5th file cumcaf.pd ,
the binary file from, CRUSH%D.c but run in
our new program PD%R.c that replaces the end
portion of CITY%PD.c . Now we develop a regression
formula to match it's 'crush' with 'rcaf' directly.
**

**
** E. Quick Match Reference?

As seen by the first id=0, which is a diameter of 0 in ** CRUSH%D.c **, is the ccrush term.
Ignoring the rest of the file record, which is the 449 particles % distribution, we converted the 3-char ccrush term to ASCII
as seen in ** rcumcaf.b2 **.
Note that the interval for 'ccrush' changes and sometimes erratic near the bottom. The test for this is to plot ccrush
for each 'n' and try a linear regression analysis. Note that we are trying to match the 6000 record #'s (y) from the match of
ccrush <-> rcaf (x). We ignore the first record and handle this later. Our first plot is seen in
RCAF # 1.

** (1) Regression Analysis **

The 6000 rcaf's (x) plot vs. the record # (y) is seen in blue and the best fit linear line is seen in yellow. Actually, it's non-linear but appears linear. Since the raw data's blue line is non-linear, we modify the regression equation's to handle the non-linear form y=ax^b instead of the linear y=bx+a. We need a working formula, that is, we don't have to go back and read each x,y again after we have the mean of x & y. This is developed in reg.dia.

** (2) Regression Analysis Results **

The results are b=0.914 & a=144.421 with correlation coefficient r=.976. Actually, the term r^2 = .953 is more important meaning that 95.3% of the observed y variation can be predicted by the regression equation, or 4.7% is random noise? Well this seems pretty good to use, however, suppose we try to tweak and fudge it for better results. Actually, we could have written a program to loop a & b until the 'r' becomes better.

** (3) A Bad Surprise **

Before this, we want to test the fit for other known equations so we generate x,y by a formula and then try to see if the non-linear regression discovers a,b exactly. Generating 150 points for y=10 x^b, where we have b from 1.0 -> .2, we see in axeb.txt the diminishing results of r with b. Not bad but thought we would have an exact match instead. So the conclusion is that with a set of related data, we may have a best fit line but not perfect on exact (r=1.000) correlated data.

** (4) Play Graphics **

What happens if we force the a,b and plot the curves to see if the equation matches the raw data. It may be better to see what's
happening rather than just watch the end result of r. Actually, it produced stunning results. This is time consuming but we only
need a,b for one time. Now with the new a,b, we have the best visual fit in RCAF # 2. This is almost a
perfect fit but we didn't calculate r but b=.678 & a=325. However, another statistical term called the standard of estimate (S.E.E)
went from 538.8 down to 118.3.

** (5) Play Calculus **

Well the black gaps between the yellow and blue lines in RCAF # 2 indicate where the data isn't an exact fit.
Can we improve on this? Well we have to match rcaf with ccrush. Using y=a*x^b with b=0.678 & a=325, we find y (rec #) to open to read ccrush.
If it isn't a match, we find another ccrush, say 10 records further by file position filp2=filp1+10*rl where we just read ccrush1
with filp1 using record length 'rl'. So the delta ccrush is dc=(ccrush2-ccrush1)/10 so then the number of records to match rcaf is
delta records dr=(rcaf-ccrush1)/dc. These results are seen in rcaf2.b1 , rcaf2.b2
and rcaf2.b3. The last column tells what % error the final ccrush from rcaf occurs. The real absolute average
difference is given by 'dif avg' on the last line. The last ( rcaf2.b3) is probably the most important range of
rcaf's used which average less than a 1% error.

Could we still improve on this? We have calculated a dy/dx slope near the ccrush by using 10 records above to find a closer approximation to rcaf.
Since the data has some noise and we may have jumpy slopes along the raw data, suppose we take the slope of the non-linear curve y=a*x^b. This is a classic
basic derivative theorem in Calculus. So with y' = dy/dx = dy(ax^b)/dx = a*dy(x^b)/dx = ab*x^(b-1), we have the slope of the best fitted curve at ccrush1.
Note this is displayed in rcaf.dia. The choice of this method
also smoothes some erratic data. So we have y'= 325*.678*x^(.678-1). Applying this techniques, we have rcaf4.b1 , rcaf4.b2
and rcaf4.b3. Even though the 1st 2 appear worse than the other method above,
the last is the important range used since the higher values with small rcaf values
calculate to small tons of dust. We are satisfied by using this method which should produce < .4% error in the most weighted used rcaf range.

With pd = particle size (id) % by decimal distribution for each slab thickness in each building's ir, we now add this to the rubbed tons' diameter ir array
ton[id]. This means we will initialize the final output file 'ctydst1.NYC' for NYC's example with 'ir' range vertically 0 -> 70 (dr=.1)'s and tons/id diameter's
horizontally 0-449u. Here, we don't know what building will be the first and last ring 'ir' so we just keep summating the file array with it's % pr until complete.
We have no more loops by continuing to read the next 'ir' ring and rewind the files
rub.%d and ** cumcaf.pd ** for matching then write out particle size distributions.
The final tons/id/ir/ is seen in ctydst1.NYC. Notice, the variable tons value length from Super Dooper Looper I.
However, in our new file ** PD%R.c **, we need a fixed value length to produce a constant record length
to use file positions 'filp' to read, summate and write again - such as ctydst1.MIA. We also changed a few variable
names and defined them in a glossary with a better pic in **Trap4 Modified.dia **.
- End Update 03/08/12

E. PSI Wave => Buildings: - 02/16/11

F. Crater vs. Building Dust Results

We convert the tons of city's building to volume and compare it to the volume of the crater.
Usually the density of dirt is considered at 2.5 g/cm^3, and we are using cement as
2.4, since we have a mixture of road and sidewalk, they are so close already, let's consider
them equal, which we should do in the radioactive fallout model to really save CPU time.
Running the model for NYC using 1,5 & 10 megaton bombs, we have the tally below:

1 Mega

City Buildings: 1.197e+8 Crushed Tons of Dust = 1.597e+9'^3 = 19.31 x Crater Vol (8.273e+7'^3)

City Buildings: 1.415e+8 Crushed + Rubbed Tons of Dust = 1.889e+9'^3 = 22.84 x Crater Vol (8.273e+7'^3)

5 Mega

City Buildings: 5.461e+8 Crushed Tons of Dust = 7.289e+9'^3 = 17.62 x Crater Vol (4.136e+8'^3)

City Buildings: 1.275e+9 Crushed + Rubbed Tons of Dust = 1.702e+10'^3 = 41.16 x Crater Vol (4.136e+8'^3)

10 Mega

City Buildings: 9.469e+8 Crushed Tons of Dust = 1.264e+10'^3 = 15.28 x Crater Vol (8.273e+8'^3)

City Buildings: 3.606e+9 Crushed + Rubbed Tons of Dust = 4.814e+10'^3 = 58.19 x Crater Vol (8.273e+8'^3)

From 1-5 megatons, if we don't include all the city buildings, the ratio is smaller.
If we include all the city buildings, the ratio is higher but if we include many of the city's outskirt smaller buildings, it is lower again
because we weighted more smaller buildings while the dust from the crater hole is proportional to the size of the megaton blast within a
large city's linear building profile. We noticed that the larger the megaton blast, the building dust increases faster than the crater hole
dust increasing the ratio to 50+ times as much for 10 megatons.
- End Update 11/20/08

Experimenting with the different slab and wall thicknesses produces somewhat complex offsetting effects. It is
difficult to pinpoint what effects each have on the total dust. Finally, we decided to go with the 4" slab thickness,
as many have quoted for older buildings combining with even thicker slabs with pillar supports. We averaged this with
the lightweight strong cement of 2.5" The wall was considered a constant of 12" thick. These are the ground rules we will apply to all the cities.
The % to dust for a 65 story building was 38.0%, and for a 20 story building was 35.6%.

10. Combining Crater and Building Dust into Mushroom Cloud:

For the crater dust, it is obliterated and immediately sucked into the nuclear cloud plasma and vaporized
into the near smallest observed dust particle of 1u so the % particle distribution of the crater dust
doesn't matter as it will all be 1u when passing through the plasma. Now, also for the building dust, much of it will
also be sucked into the stem of the nuclear mushroom cloud and pass through the plasma reducing it's
size to 1u also. The question is, how much of the building dust will be swept up?

In order to answer this, let's assume a minimum horizontal wind speed for the largest of particles (450u). The derived fall velocity formula is complex and will be referenced more when developing the RADFO model components. For reference, the calculated terminal velocity of 1.3564e+1 fps or 9.25 mph for <= 450u to 2.5082e-2 fps or .02 mph for <= 10u. How far out from ground zero is the inrush of air?