Topics:

• Error analysis of the Euler algorithm: Local error proportional to h^2, global error is proportional to h. (ODE).
• Runge-Kutta (4) algorithm: Local error is proportional to h^5, global error is proportional to h^4.
  Wikipedia: Runge-Kutta
  
• Dormand-Prince algorithm: Read the entry in PHASER Help and also see
  Wikipedia: Dormand-Prince_method
• How to compute e: comparison of several algorithms.
• What to do in practice: In real life applications it impossible to be certain if numerically generated solutions of differential equations are reliable, or how large the errors might be. The golden rule of practice is to compute the same solution using several algorithms with various step sizes. If the variations in these computations are large (more than you can accept for a particular application) your problem is most likely a dangerous one. If the variations in these computations are small enough, most likely you get a reliable answer out of the machine.

Assignment:

1. Use PHASER for solving the initial-value problem x' = x, x(0) = 1 at t=1, that is x(1)= e = 2.7182818284590452354 with several algorithms as follows. You should open the Console in the Numerics Editor to get the stats for your computations:
   • First use Runge-Kutta with step h = 0.01. Note that you get about 9 digits correct to the right of the decimal. How long does your computation take?
   • Now determine the "largest" step size in the Improved Euler algorithm to get the same number of correct digits. How long does it take to compute?
   • Now determine the "largest" step size in the Euler algorithm to get the same number of correct digits. How LOooooNG does it take to compute?
   • Explain your results in terms of the orders of the algorithms you use.

2. Consider the "explosion" problem x'= x^2, x(0)=1, and try to compute x(1). As we saw this value is infinite. Now in Phaser, compute this number using steps h = 0.01 and h = 0.005 and algorithms Euler, Improved Euler, Runge-Kutta4, Dormand-Prince5(4), and Dormand-Prince8(5,3). Following the golden rule above, conclude that this problem is indeed dangerous to compute. Next, In the Dormand-Prince 5(4), set the step h = 0.01 and decrease the Absolute error and Relative error from 1.0E-7 until Phaser refuses to compute the value x(1). What is the smallest error tolerance for which PHASER refuses to compute infinity? Monitoring the Console in the Numericas Editor, report at which t value PHASER stops computations for your set tolerances.

3. Again, consider the "explosion" problem x'= x^2, x(0)=1. Verify that the solution of this initial-value problem is x(t) = 1/(1-t). Now, we will compute x(0.9) using Euler's algorithm (staying away from the danger zone) and analyze the error.
   1. Note that the exact value of x(0.9) = 1/(1- 0.9).
   2. Recall that in lectures we showed that the global error bound for Euler is M*t_final*h where M is maximum of x''(t) for t between 0 and t_final. Determine the value of the proportionality constant M.
   3. Set h = 0.1 and compute Euler approximation for x(0.9). Determine the actual error in this approximation. Using the global error bound formula above, determine the the value of the theoretical error bound. How does this theoretical error bound compare to actual error?
   4. Repeat the calculations above for h = 0.001. Is the theoretical error bound now more realistic bound of the actual error?

4. Gompertz model of cancer growth: The differential equation
   x' = a*(exp(-b*t))*x
   is used to describe the growth of a tumor, where x(t) is a measure of its size (e.g. wieght or number of cells), and a and b are parameters specific to a particular tumor. To get started, let us take a = 3 and b = 2, and x(0) = 5; the solution in the XivsTime view is shown in the image below.

   

   Load the image into your PHASER by clicking on the picture. Explain what this model predicts for the tumor growth. Experiment with several different initial tumor sizes. How does the future size of tumor depend on its initial size? Experiment with several values of the parameters a and b. Explain the roles of these parameters in the tumor growth.

5. Install R: Please go to the website http://cran.opensourceresources.org/ and download and install the program R on your laptop. For more information on R, see http://www.r-project.org/ We will start using R next week.