Mathematics is the door and the key to science.
Lecture #5
Numerical Methods
1) Mathematical Model
geometry
equations
2) Discrete Algebratic Approximation
3) Simulation Algorithm
4) Technical Plan
5) Code Validation
6) Timing and memory tests
7) Simulation
8) Interpretation and visualization
9) Analysis and Conclusions
How do I know this works?
syntax errors
does the darn thing compile?
Fix minor typographical errors
logic errors
are there problems with the coding?
algorithm errors
does the code solve the mathematical model?
model errors
is the mathematical model correct?
symmetry
does geometry affect the results
convergence
do the algorithms actually converge
conservation
are conserved quantities conserved?
consistency
is it consistent with previously published results or analytic solutions
test your codes!
Check for conservation!
Check for convergence!
Compare with analytic approximations!
simulations which rely on random numbers
good random numbers are essential
wide variety of simulations can be done using random numbers
1. Drop a large number of points at random locations inside the box.
2. Determine what fraction of these points is within the circle.
3. Determine the relative area using this fraction.
1. Generate a number of test states
2. Determine the best of the test cases.
3. Mutate the best test cases look.
4. Repeat steps 2 and 3 until the problem is solved.
First set of paths.
Mutated versions of the best path.
Simulations which have large numbers of causally unrelated events.
1. Generate a random initial state
2. Iterate state until it is stable.
3. Bin or track the final state of this system.
4. Repeat steps 1-3 until enough statistically significant results have been generated.
1. Throw a rock off of Mars in a random direction at a random time.
2. Determine if the rock hits Earth.
3. Repeat 1-2 until enough statistics have been created for meaningful results.
Chains of discrete events.
The ultimate success or failure of a particular state is determined by a random number generator and the probability from the previous state.
Previous states do not directly influence the next state.
1. Atoms decay from their original state to a final state.
2. The decay of a given atom is independent of the other atoms in the material.
3. The chances of a specific atom decaying in any given period of time is always the same.
4. All atoms have the same chance of decay.
What type of Monte Carlo Simulation is this?
Genetic Algorithm?
Discrete Event?
Markov Chain?
1. Set the probability of atoms surviving in a given time period.
2. Generate an atom in the first time bin.
3. Determine if the atom survives in this time bin.
4. If it survives, move to the next time bin and repeat 3 until you have reached the final time bin.
5. If it does not survive, increment the number of decays in that particular time bin.
6. Repeat steps 2-5 until done.
Conrad Hilton
Michael Wilding
Mike Todd
Eddie Fisher
Richard Burton
Richard Burton
John Warner
Larry Fortensky
f(x) = f(x0) + f'(x0)(x - x0) + 1/2 f''(x0)(x -x0)2 + ... + 1/n! f(n)(x0) (x-x0)n
+ 1/(n+1)! f(n+1)(z) (x-x0)n+1
where z is some point between x and
x0
We can drop the remainer term and approximate
this with a polynomial of degree n.
where does it come from?
first order approximation
The largest error in the Taylor series is related to the final term of the series.
If :
max (abs ( f(n+1)(z) )) <= M
The largest error will be approximately :
M/(n+1)! (x-x0)n+1
This is NOT an exact relationship.
The true value depends on the function.
given a data set, how do you find the values at points between date elements?
You must interpolate the values
Given by piecewise polynomial of the form:
p(x) = f1 + (f1 - f0)(x - x1) / (x1 -x0)
Where x = [x0, x1]
Given:
p(x) = f1 + (f1 - f0)(x - x1) / (x1 -x0)
we can also write:
p(x) = a0 + a1 x
where we require
p(xi) = fi
so we can form a set linear equations to determine a0 and a1
we can extend the idea of linear interpolation unto a general polynomial of degree n
Then
p(x) = a0 + a1 x + a2 x + ... + an x
with the constraint equations:
a0 + a1 x1 + a2 x2 + ... + an xni = fi
where i = 0... N
With the constraint equations at each point, we can write the set of linear equations to 0be solved as:
1 x0 x02 ... x0na0f0
1 x1 x12 ... x1na1f1
... ... ... ... = ...
1 xn xn2 ... xnn anfn
This matrix must be inverted for each interval that is solved for.
for n unique points, we can define
(x- x0) (x-x1) ... (x- xn)
lj(x) = ----------------------
(xj - x0)(xj -x1) ... (xj - xn)
( x - xk)
= (product) -------
k = 0,n ( xj - xk)
k != j
With this definition, you can see that
lj( xi) = 1, if i =j
0, if i != j
with this property
lj( xi) = 1, if i =j
0, if i != j
of the Lagrange polynomial, we can create a polynomial for interpolation.
P(x) = (sum) lj(x) fj
j= 0,n
the higher order polynomials are great for interpolation if you have a small number of points OR lots of cpu time.
However, it is sometimes useful to break the problem into a series of piecewise functions.
Much lower CPU time, and still may do a good job with interpolation.
one type of piecewise linear function is a spline.
Splines are piecewise linear functions which have continuous derivatives up to the order of the spline.
Linear spline- continuous 0th derivative
Quadratic spline - continuous 1st
derivative
Cubic spline - continous 2nd derivative
the cubic basis spline is based on a set of basis equations
Bi(x) = the value of this function at each point in the interval
if we let:
c(x) = sum ( alpha(i) Bi(x) )
i=1, n
we can use the constraint equations
c(xi) = fi over I= 1,...,n
to find alpha(I).
With the properties of the B-spline,
this becomes a tridiagonal system and can be easily solved.
Splines are continuous, but not continuous over all derivatives. It is computationally cheap to calculate a spline, since only (2m+1) data values are involved where m is the order of the spline.
It is expensive to create an nth order polynomial for a data set of size n. Large matrixes must be created or at the very least, large multiplications must be calculated. However, the polynomial is more continuous and MAY be a better fit to the data.
the least squares method does NOT attempt to exactly fit the data at each point
the 'best" fit of the data is found by finding the minimum residual in a given set of data
sum wi ( fi - p(xi))2
I=1,m
must be minimized
least squares minimizes the residuals when data is compared to a model
interpolation forces data to be exactly fit with a polynomial
Learn about shell scripts
sed awk sort if case
sed - editing file streams
awk - editing file streams
sort - sorting file streams
shell script
removing old files with *.o extensions
Simple way to do changes in a file.
sed 'list of ed commands' filenames...
The output will be standard out unless you pipe it to a file.
Simple ed commands:
s/old / new / f
where old is the old string
new is the new string
f is either
p = print
g = replace all occurances
w file = write to 'file"
/ pattern /command
executes the command every time a pattern is found
/ pattern /!command
executes the command every time a pattern is NOT found
sed -n 'command' file
will only print lines explicitly listed in the command string
^beginning of line
$end of line
.Any single character
[...] any set of characters in range
[g-l] for example
[^...] any set of characters NOT in range
*anything
r*zero or more occurences of r
a\append lines to output until one not ending in \
b label branch to label
c\change lines to following text
ddelete lines
i\insert following text before output
llist line- make all characters visible
r fileread file
t labelbranch to label if substitution is made
on current line
= print line number
: labelprint label
sed '/dog/p' test
will take all the lines from the file 'test" which have the string 'dog" in them and print them.
These lines will appear to be double printed.
sed -n '/dog/p' test
will print only the lines with the string 'dog" from the file test. Only these lines will appear.
Sed 's/dog/cat/' test
will print all lines of the file 'test" and substitute the FIRST occurrence of 'dog" on each line into the string 'cat"
sed 's/^dog/cat/' test
will replace all the occurrences where 'dog" appears at the beginning of a line into 'cat".
sed 's/.og/cat/g' test
will replace any occurance of the string '.og" to 'cat" anywhere in the line.
sed '/dog/s/g/ggie/g' file
will replace find any line with 'dog" on it and then replace all the 'g"'s with 'ggie".
general syntax:
awk 'program' filenames
where a 'program' has a syntax of:
pattern {action}
The pattern is just a regular pattern expression similar to sed. For example,
/dog/ will try to match the string 'dog"
/.dog/ will match all the dogs at the beginning of a line
the actions include any set of functions.
Functions include
print - print a set of variables or message
printf - print a formatted set of variables - similar to C
Math functions in awk
cos(expr) - take the cosine
sin(expr) - take the sin
exp(expr) - take the exponent
int(expr) - take the integer part
log(expr) - take the natural log
getline() - get the next input line
length(s) - length of string s
substr(s,m,n) - take the substring
split(s,a,c) - split string s on character c into a[1]...a[n] and return the value n
index(s1,s2) - return postion of s2 in s1
Variables in awk
FILENAME - current filename
FS - field separator (default=blank+tab)
NR = number of input records
NF - number of fields in input record
1...NF - variable for the value of each record
0 - entire current line
operators include
=
+=
-=
*=
/=
'=
||
&&
!
>
>=
<
<=
==
!=
~
!~match and not match
ls -l | awk '{ print NR, $0, $1}'
ls -l | awk '/dog/ { print NR, $1 }'
ls -l | awk '{ print length($1), $1}'