Физика из Копенгагена

Python simulation

2006-10-01T11:09:00.000+02:00

The paper on this weeks simulation can be found here (Hurray for TeX!)

Original program with inserted gauss distributions can be found here (credit: Ian Bearden)

2006-09-28T15:15:00.000+02:00

Kanondata.xls

LIFTOFF!

2006-09-25T22:06:00.002+02:00

This week was somehow different from the last three, since our task this time was to investigate the use of physical simulations in a particular field of interest. I want to talk a little about computing and using a physical model for a given problem in the field of astronautics, which is defined as "the theory and practice of navigation through air or space".

I have a very deep passion for this field and want to make use of some calculations and models i have been working on with my current mathematical and physical knowledge.
The art of rocket propulsion has hardly changed for the last 50 years, the most popular rockets is still cheap and simple russian IBCM derived rockets. Utilizing modern technology has not changed this fact, the physics behind accelerating a particle into orbit has not changed. The main reason for the rocketing (literally speaking) prize tag of placing a payload in orbit is the need to carry all the energy needed to reach orbit velocity. Todays rocket design consists of two main designs, liquid fueled rockets (as the soyuz, SSME, Ariane5 main engine) and solid(most newly developed rocket boosters). These two designs, although very different in structure, has alot of things in common, including the need for an oxidizing agent to burn their fuel.
Carrying this oxidizer throughout the atmospheric trip adds an enormous mass to the rocket, by utilizing the air in the atmosphere as oxidizer over 50%! can be saved from the original takeoff mass.

Hot kerosene anyone?

Engines capable of extracting oxygen from the atmosphere and burning it with rocket fuel has been manufactured with very limited success. This situation however could very well change very rapidly, considering these engines both have applications within aeroplanes and rockets.
I have taken the liberty to demonstrate how a physical model can describe the ascent for this engine, and which applications and limitations the model might have, instead of discussing more generally how the physical simulations might be applied.

The model
Before we continue we must define what our model describes. The efficiency of a ramjet engine depends on the availability of oxygen. Our model therefore needs to describe the availability of oxygen as a function of height and speed. This must be a logical choice, since the density of the atmosphere decreases with its height, and likewise the amount of available oxygen increases if the engine goes faster through the atmosphere at a constant height.

We therefore have a function a function with two parameters, height and speed.

First off we look at the amount of availible oxygen at ground level with a speed of 1m/s. The first tool we need is the ideel gas law equation, that we recognise as:

pV = nRT

Isolating the parameter n(moles) we get:

n = (pV)/(RT)

From this equation we can see that the amount of moles are related to the temperature and pressure, the volume is the speed of the rocket times the area of the intake. If we define the area of intake as 2m², and use the paramter xfor the speed of the rocket, we get:

V = 2x

The temperature is not constant along the ascent of the rocket, so the best we can do here, is to use a mean value for the intended operational trajectory. The best solution would be to assign it a varying parameter, this however gets a bit more complicated.
The other parameter is the pressure, and we know from empirically developed mathematical models that the pressure decreases 100% for each 5km you ascent through the atmosphere.

We can therefore express our pressure at a given height:

p = (1.013*10^5 N/m²)/(2^y)

This function satisfy our demand that the pressure must decrease 100% for each 5km interval, where y= Z(integer).

We know that oxygen only makes up approximately 20% of the atmosphere, and we are only interested in oxygen since nitrogen cant replace the oxidizer function. We must therefore revise our pressure expression to:

p(O2) = (2.026*10^4 N/m²)/(2^y)

We can now use our modules to fabricate a complete function:

n = (pV)/(RT)

V = 2x

p(O2) = ((2.026*10^4 N/m²)2x)/(2^y)

n = ((2.026*10^4 N/m²)2x)/(2^y)(RT)

c = RT

R = 8.314 J/mol K

T = 273ºK

c = 2269.7 J/mol

n = ((2.026*10^4 N/m²)2x)/(2^y)(c)

We simply divide our c into the expression for the pressure and get:

n = ((8.926 N mol/J m²)2x)/(2^y)

We now have a complete function that describes the amount of available oxygen with the two parameters speed and height. But instead of having the function giving us values in moles, it would be more convinient to let it return values in kilograms.
We therefore use the linear relation:

n = m/M
m = nM
M(O2) = 32g/mol

Just by multiplying the last formula with 32g/mol, we can make the function more practical.

m(x,y) = ((285.64 N mol/J m²)2x)/(2^y)

Instead of having the function returning values in grams, we divide it by 10^3, this way, we will have our values in kilograms, which is a bit more convinient. At the same time, we can also multiply our only constant number with 2, to remove it from our parameter x.

m(x,y) = (0.57x)/(2^y)

Now we're done with building our formula, we must now look at it and see how much information we can extract from it.
Basicly it tells us how many kilograms of oxygen that is availible in a 2m^2 frame moving with a speed x and at a height y = Z(5km), where Z is any integer. eg when y = 3, the height corresponds to 15km.

Let's take a look at it:

The x axis describes the speed of the particle(rocket) in m/s, the interval ranges from 0 to Mach 10. The y axis decribes the height of the particle from 0 to 50km. The z axis is our function and decribes the amount of oxygen availible under the given circumstances.

The graph actually tells us some very interesting things. The amount of available oxygen at altitudes above y ~ 4 (20km) is very scarce and represents a challenging engineering problem because the intended operational range for this type of engine only goes to max 20km (approx.) as one can see from the graph. One must utilize the dense atmosphere for gaining maximum speed before ascending vertically upwards for a orbit injection.
So just by designing this function we have gained some valuable information about how such an enging behaves with respect to atmospheric pressure and speed.
This was a very simplified example and is of course only valid for fictional applications and for showing the use of physical simulations within the astronautic branch.

Right.......so tomorrow i will probably update my blog with a comment or two on climate models.

Projectile motion

2006-09-11T17:00:00.000+02:00

Purpose
The purpose of the second laboratory exercise was to work with uncertainties and statistics in relation to an experiment determining the start velocity of a ball fired from a cannon as shown on figure 1.

Figure 1

Theory and formulas
The primary task was to determine a value for the start velocity as a function of a measurable paramter. Several parameters were suggested, but the easiest to measure was agreed upon to be the length the ball travelled from the cannon until it hit the ground. By using the theory about a projectile movement in a gravitational field, the group managed to isolate the parameter start velocity as a function of time. First of all we need to look at the theory behind the projectile move in a gravitational field.

We consider a particle in a gravitational field, by using Newtons laws of motion we know that

x(t) = v0 cosΘt
y(t) = v0 sinΘ - gt

This vector describes a projectile trajectory in a gravitational field, and by exploiting our knowledge about vectors we can determine an expression for the start velocity as a function of shooting length by forcing the y(t) vector to zero.

v0 sinΘ - gt = 0
t = (v0 sinΘ)/(g)

We can now substitute this expression for the time, when the ball hits the ground (y vector equals zero) into the x vector.

x = v0 cosΘt
x = v0 cosΘ(v0 sinΘ)/(g)
x = (v0^2)(sin2Θ)/(g)

(I) v0 = sqrt(xg/sin2Θ)

This vector describes the position of the particle as a function of time, by forcing the y vector equal to zero, we define the position 0 relative to the y axis, which is the same as the point where the particle either takes off or hits the ground.
By using formula (I) we see that we can calculate the start velocity as a function of the length it has travelled projected onto the x-axis. We can also see that the relation between the length and the start velocity is linear.

Thoughts behind the experiment
Determining the best possible setup was an important factor in minimizing the errors associated with the exercise. Before the experiment could begin, the group discussed different setups. The single most important factor was that of the table height being equal to the firing height. By approximating the two heigths to eachother one could look at the system as a parabola movement without having to calculate in the negative y axis, which simplifies things a bit. Furthermore the importance of a stable setup was discussed, for instance, it was questioned whether or not the impulse change for the ball cannon would result in increasingly inaccurate results. The solution to this problem was to place a heavy brick on top of the cannon, but even though this might have reduced the effect it is still present.
It was discussed whether or not one should calculate the total length of the trajectory as a single dimension problem or expand it and look at the deviation from the center axis of the throw and calculate the length by using an algebraic approach. It was determined after the experiment that the only way forward was to calculate the problem as a one-dimensional system since no center axis had been marked on the experiment paper.

Statistical approach
For a comprehensive overview of the combined data from the experiment please see (45degree excel sheet and 75degree excel sheet). (windows file extension compatible)

FTP transfer

45degree excel sheet @ intro-soma
75degree excel sheet @ intro-soma

Our original goal was to determine the start velocity of the ball. We need to find the value indirect by measuring the one-dimensional trajectory length, and calculate the start velocity from this value. The first point in processing the data statistically is to determine the uncertainties in the direct measurements and afterwards apply these uncertainties to the interesting value, in this case the start velocity.

Uncertainty in measuring the length
The length from the cannon to the impact table was measured by a ruler with 0.001m precision but since we needed to measure the length in two revolutions for the 45degree experiment, the measuring uncertainty rises to twice of the original value.

Measuring the length for the 45degree experiment
( 1.000m +/- 0.001m ) + ( 0.033m +/- 0.001m ) = 1.033m +/- 0.002m

Measuring the length for the 75degree experiment
( 0.375m +/- 0.001m ) = 0.375m +/- 0.001m

This length however does not represent the total length of the trajectory, since it only reflects the length between the impact tables bezel and the canon. The next measurement is done on the paper and added to the original length, which gives us the following uncertainties.

Total length for the 45degree experiment
( 1.033m +/- 0.002m ) + ( impact length +/- 0.001m )

Total length for the 75degree experiment
( 0.375m +/- 0.001m ) + ( impact length +/- 0.001m )

Since the start velocity is a direct product of the measured length (linear relation), the uncertainty used to compute the end result is the same as the uncertainty used to describe the length.

I will return later on to the uncertainty but first we need to calculate the standard deviation in our samples.

Statistical dispersion
Being able to measure the deviation in a sample is an important factor when one wants to conclude anything about the statistical behaviour of the variable factor in a isolated system. In this exercise we wanted to define the ball canon by describing its ability to fire projectiles with a start velocity combined with a measurement uncertainty.
We saw earlier on that the relation between length and start velocity is linear and therefore we can directly apply the uncertainty from measuring the length to the start velocity, this however would not give a clear picture of the situation and this is where our standard deviation comes into force.
Besides describing the ball canon with an uncertainty taken from the length measurement we also need to look at the deviation from the measured data to the mean value.

The formual above calculates the standard deviation, which has been calculated in our experiment to the following values.
(Please see the excel sheets for an overview in the calculations. I will only state the resulting figures here)

Standard deviation 45 degree: 0.090
Standard deviation 75 degree: 0.156

After finding the standard deviation it would be interesting to state the standard deviation of mean for the two measurement series. SDOM gives us the uncertainty for our final result.

SDOM 45 degree: 0.019
SDOM 75 degree: 0.029

We can now state the start velocity with a known uncertainty taken from the SDOM, which gives us the following start velocities.

v0(45) = 3.58msˉ¹ (+/-) 0.019
v0(75) = 3.41msˉ¹ (+/-) 0.029

The page will be updated as soon as i figure out how to apply the uncertainties from the length measurement, since the above result is pretty optimistic......

Projectile motion (notes)

2006-09-11T16:45:00.000+02:00

As i don't seem to be the only one struggling with understanding the principles of statistical data processing i've included a couple of useful links to sites explaining about the terms in varying scientific complexity.

Standard deviation @ wikipedia
Standard deviation @ St. Louis Batesville public school
Root Mean Square (RMS) @ wikipedia
A comprehensive excel guide for understanding the terms graphically @ Harvey Mudd College Physics Department

If anyone else has any useful informations please feel free to comment and i will add them to the list.

Measuring g

1990-09-18T19:07:00.000+02:00

This week it was time to get our hands dirty with some more statistics related experiments. I was so lucky to end up in a group with two guys that apparently seemed to have gone through the theory before showing up in the lab, so the planning stage was quickly executed. The purpose of the exercise was to determine the gravitational acceleration a mass feels at the surface of the earth and to discuss our result in relation to the uncertainties connected to the experiment. We had different bobs and strings at our disposal, and by using these two things we wanted to determine g. At our initial lecture we had been introduced to the pendulum as a way of determining g and it also seemed like the most precise and easiest way. Before commenting on the experiment and the results i will use some lines to show the underlying mathematics.

At first we need to make three assumptions about the pendulum.

The string on which the bob is attached remains taut and massless.
The bobs mass is treated as a point mass (all mass is concentrated in one point)
The motion only takes place in a two dimensional system.

Figure 1

Figure 1 shows a pendulum with two forces acting upon it. The gravitational force can be resolved into components parallel and perpendicular in relation to the instantaneous velocity of the bob. The tension only has a radial component which is cancelled by the radial component of the gravity. The only force remaining in the system is the radial component of the gravity. The radial gravity force is not constant, but a function of the angle theta that the swing forms with the vertical axis. Knowing Newton's second law we can now describe the forces acting on our body.

F = ma
F = -mgsin(theta)

By looking at very small angle intervals we know that an angle approximates to the sin value of the angle, in other words: sin(theta) ~ theta

F = -mg(theta)

By looking at figure 1 we kan express theta as (x/l)

F = -(mgx/l)

By inserting the known values (mg/l) in a given constant k, we can rewrite the expression

F = -kx

This expression is a well known physical model defined as a harmonic oscillator and has the specific property that it will oscillate around an equilibrium with the period

T = 2pi sqrt(m/k)
T = 2pi sqrt(ml/mg)
T = 2pi sqrt(l/g)

We can now isolate g and thereby calculate the gravitational acceleration by known parameters, such as l and T.

g = l(2pi/T)²

However, this model only works with small angles since we have made the assumption that the arc length can be expressed as radius times the angle (in radians). Figure 2 shows the deviation for the model in relation to the real value, as one can notice, the uncertainty for angles below 30º approximates zero. In our experiment we only worked with angles <>

Figure 2 (credit: wikipedia)

Experiment setup
To measure g we used a pendulum suspended by two strings. By using this method we could minimize the uncontrolable string drag that arises from the fact that the string isnt completely frictionless in the spot where it's fixed. We then measured a period for the pendulum, 100 times with both a stopwatch and the smart timer system. We assumed that there was no uncertainty in measuring the time and therefore the only factors contributing with uncertainty is the string length and the standard deviation of mean in the period measurements.

Calculating the propagation of uncertainties
Figure 3 shows our setup with the imaginary tension force being the arithmetic sum of the two strings.

Figure 3 (credit: my amazing paint skills)

The measuring uncertainty was estimated to be (+/-) 0.0005m for a, and (+/-) 0.01m for b respectively:

a = 1.25m (+/-) 0.01m
b = 0.208m (+/-) 0.002m

Knowing that the fractional uncertainty in a² is just two times that of a, we can write:

a² = 1.56m² (+/-) 0.02m
b² = 0.043m² (+/-) 0.004m

Now we can use the quadratic sum of the previous uncertainties to calculate our best estimation for c with uncertainties.

∆c = sqrt( (0.02m)² + (0.004m)² )
∆c = (+/-) 0.02m

The numerical value of c from our measurements can be expressed as:

c = sqrt(a² + b²)
c = sqrt(1.56m² + 0.043m²)
c = 1.23m (+/-) 0.02m

We now know the uncertainty by which we can define our string length. This however, is not the only uncertainty in our experiment. We assumed that the stopwatches used for our experiment was without uncertainty, but other factors in the system might have contributed to varying period measurements and these interactions can be viewed in the different values of T. We define T as the best estimate (+/-) SDOM.

For a comprehensive overview of the statistical data, please see the excel sheets.

measuring_g.xls @ intro-soma
measuring_g.xls @ fys.ku.dk

Our best estimation for the period.

T(stopwatch) = 2.18s (+/-) 0.01s
T(smart timer) = 2.2296s

Figure 4

Figure 5

Figure 4 and 5 shows the histogram for both the stopwatch measurement and the sts. As one might notice, the interval in the sts histogram is muss less than that of the stopwatch histogram. Both curves approximate the shape of a parabola and have a considerable large amount of data points within their standard deviation.

Calculating g
By defining the pendulum length and the period as the only uncertainties in our experiment we can now calculate g as:

g = l(2pi/T)²
g(stopwatch) = 1.233m(2pi/2.18s)²
g(stopwatch) = 10.25m/s²

g(sts) = 1.233m(2pi/2.23s)²
g(sts) = 9.79m/s²

The two measurements has a natural uncertainty given by the uncertainty from the length measurement + the uncertainty in the period measurement. By applying the sum of the two uncertainties we can now state that our experiment showed that g lies in the interval of;

g(sts) = 9.79m/s² (+/-) 0.16m/s²
g(stopwatch) = 10.25m/s² (+/-) 0.21m/s²

The value of g from the sts experiment clearly lies within the normally accepted value of g. This however is not the case with our stopwatch experiment that has a value outside the generally accepted g interval.