September 9, 2014 Class 5---Error Analysis for the Euler Method (II),
Introduction to Newton's 2nd Law
Activities:
Error estimates
Notes in error estimate may be more readable in the pdf file
~sg/euler_error_analysis.pdf
There is also a link on the class web page.
Most of of this class was devoted to understanding the
error in the Euler algorithm for integrating a differential
equation. We did a careful estimate of the error for the
simple case
dy/dx = f(x).
This is simple because the right hand side depends on x,
not y(x). So this is just doing an integral.
If we use a step size dx, the error in a bin with left hand
edge at x is:
(1/2!) (dx)^2 f'(x) + (1/3!) (dx)^3 f''(x) + higher order terms
In our first example of the Euler program, f(x)=2x. Thus,
f'(x)=2 and f''(x)=0, so the second term above is 0. Why
is the error order dx, not (dx)^2? That is because the
above expression is for each bin. In our example, we
integrate x from 1 to 2, so the number of bins n=(1/dx).
Thus, the total error is dx for this case.
In our second example, f(x)=3x^2. f'(x)=6x and f''(x)=6,
so the second error term does not vanish. If you very
carefully evalute the two error terms you will find that
the first term is
3 (dx)^2 * n (1 + (n-1)*dx/2 )
where n=1/dx. Using the fact that n*dx=1, this first
term can be written as
3 dx (3/2 -dx/2) or 9/2 *dx - 3/2 * (dx)^2.
Actually, only this term requires care. The f'' term is
easily seen to be
(dx)^2.
Adding the two terms, our final result is
9/2 *dx - 1/2 * (dx)^2.
Before we did the careful estimate, we used a program polyfit
to fit the result as a function of step size with a second
order polynomial. That result agrees with the error calculation
above. It is also possible to plot the error
divided by dt and find that the coefficient of the second
term is negative, not positive as we would get from just
looking at the second term in the Taylor expansion.
Solving Newton's Law of Motion
We discussed how to turn a 2nd order
differential equation into a series of coupled first order
differential equations. We also discussed why position,
velocity and acceleration are so important, and the next
derivative, called the jerk, is not well known.
The program fall_euler_ansi.c in ~sg contains code to
integrate the equations for a falling body.
It is set up for the Euler algorithm. We will looked at the
code in some detail. One key issue of discussion was pointer
dereferencing as when a pointer variable x_ptr is used with the
* operator in front of it.
*x_ptr = 6.3
places the value 6.3 in the memory location that x_prt points to.
If desired, we can talk more about this later.
In the next class, we will talk about energy conservation.