Name: Math 428
Exam 1
5 March 2002

Yes, mathematics has two faces; it is the rigorous science of Euclid but it is also something else. Mathematics presented in the Euclidean way appears as a systematic, deductive science; but mathematics in the making appears as an experimental, inductive science.
- from ``How to Solve It'' by G. Polya

Instructions: Show all work to receive full or partial credit. You may use a scientific calculator on this exam. All University rules and guidelines for student conduct are applicable.

1.

[15 pts] Circle the integrals for which you would expect the midpoint quadrature to converge to the exact integral as one increases the number of panels.

$$\int_0^{\pi^2} \frac{t^3-3e^t}{5-\cos(\sqrt{t})} dt, \ \ \ \ ... ...} dt, \ \ \ \ \ \int_0^{\pi^2} \frac{t^3-3e^t}{1-\sin(t)} dt.$$

Explain your reasoning.

All of our convergence analysis of quadratures so far has been based on the assumptions that one can expand the functions in a Taylor series. To do so, the functions must be continuous and have the right number of derivatives. In the case of the midpoint rule, we need three derivatives. By inspection, the first two functions, have no problems with asymptotes and plenty of derivatives, so we would expect the midpoint rule to converge. The latter two functions have asymptotes so they are not even continuous on the domain of integration, and we would not expect the midpoint rule to converge to the integral, assuming that the integral converges.

2.

[15 pts] A student estimates the solution to a differential equation from an initial time t=0 out to a final time $t=4\pi$. She acquires the following results for three runs.

h	$y(4\pi)$
$\frac{\pi}{16}$	1.07455259139819
$\frac{\pi}{32}$	1.04051345264161
$\frac{\pi}{64}$	1.02121990730698

A) Use Richardson extrapolation once to improve the best run.

No problem.

2(1.02121990730698)-1.04051345264161 = 1.001926362

B) She needs four digits of accuracy. What step size should she use?

To make the relative error smaller than 10^-3, so that the first four digits are credible, we expect

|I_h-I_h/2| = D h + O(h²).

We can measure these quantities.

$$\begin{eqnarray*}\vert I_{\frac{\pi}{16}}-I_{\frac{\pi}{32}}\vert & = & 0.034039... ...t I_{\frac{\pi}{32}}-I_{\frac{\pi}{64}}\vert & = & 0.0192935453. \end{eqnarray*}$$

Using these measurements, we find that D is about 0.1734 or 0.1965 depending upon which measurement you use. Taking the larger as the more conservative estimate, we need

0.1965 h < 10^-3.

So, we need h<0.005. If you like $\pi$ divided by powers of 2, $\frac{\pi}{1024}$ will do it.

3.

[20 pts] Find a quadrature to estimate $\int_0^1 f(x) dx$ knowing f(0), $f'\left(\frac{1}{4}\right)$ and f(1). In other words, determine nonzero constants A₁, A₂ and A₃ where

$$\int_0^1 f(x) dx \approx A_1 f(0) + A_2 f'\left(\frac{1}{4}\r... ...\frac{1}{4}}\right)+ A_3 f(x_{i+1})\right]. \end{displaymath}}$$

What is the order of accuracy for this quadrature?

This is taken from your homework (Set #1, question 5). We expect the quadrature to have third order accuracy because it involves three pieces of information about the function. We found in the homework and subsequent discussions in class that the method is a third order method if $\xi = 1/2$, but this is not the case in this exam question.

Many students made the mistake of trying to use divided differences or Lagrange polynomials for this problem. Remember, that these are interpolation techniques, but this problem involves mixed data about f and f' so these methods are not applicable. You must proceed with undetermined coefficients: Try

p(x) = c₀ + c₁ x + c₂ x².

If you fit the polynomial through the given data, you find by solving the resulting linear system that

$$\begin{eqnarray*}c_0 & = & f(0), \\ c_1 & = & f(0)-f(1) + 2 f'\left(\frac{1}{4}... ...), \\ c_2 & = & 2 f(1) - 2 f(0) - 2 f'\left(\frac{1}{4}\right). \end{eqnarray*}$$

Now, if we plug this into

$$\int_0^1 p(x) dx = A_1 f(0) + A_2 f'\left(\frac{1}{4}\right) + A_3 f(1),$$

we find that

$$\begin{eqnarray*}A_1 & = & \frac{5}{6}, \\ A_2 & = & \frac{1}{3}, \\ A_3 & = & \frac{1}{6}. \end{eqnarray*}$$

Thus, we find that

$$\int_0^1 f(x) dx = \frac{5}{6} f(0) + \frac{1}{3} f'\left(\frac{1}{4}\right)+ \frac{1}{6} f(1) + O(h^3).$$

4.

[15 pts] Consider the differential equation

$$y' = f(y,t), \ \ \ \ y(0) = y_0.$$

Suppose computing f(y,t) is expensive, so one chooses to modify Euler's Method

$$u_{k+1} = u_k + h f(u_k,t_k), \ \ \ \ u_0 = y_0.$$

Instead, one chooses to approximate f(u_k,t_k) with the finite difference quotient,

$$y'(t_k) = f(u_k,t_k) \approx \frac{u_k - u_{k-1}}{h}$$

from the previous timestep to create the new method:

$$u_{k+1} = u_k + h \left( \frac{u_k - u_{k-1}}{h}\right), \ \ \ \ u_0 = y_0.$$

What is the local truncation error for this new method?

Recall: local truncation error is the error in one time step assuming that you have exact initial data. That is, we assume that u_k = y_kand u_k-1 = y_k-1. We need to try to estimate the difference between y_k+1 and u_k+1. After some simplication, the numerical method is

u_k+1 = 2 u_k - u_k-1.

Inspecting the method, we see that we need to Taylor expand near u_kand u_k-1.

$$\begin{eqnarray*}y_{k+1} & = & y_k + y'_k h + \frac{1}{2} y''_k h^2 + \frac{1}{3... ...+ \frac{1}{2} y''_{k-1} (2h)^2 + \frac{1}{3!} y'''(\xi_2) (2h)^3 \end{eqnarray*}$$

So

$$\begin{eqnarray*}y_{k+1}-u_{k+1} & = & y_{k+1}-2 u_k + u_{k-1} \\ & = & y_{k+1}... ...+ y'''(\xi_3) h^3 + \frac{1}{2} y'''(\xi_4) h^3 \\ & = & O(h^2) \end{eqnarray*}$$

5.

[15 pts] Prove that the backward Euler method converges to the exact solution of

$$y' = f(y,t), \ \ \ \ y(0) = y_0.$$

What is the order of convergence?

This problem follows the same steps as forward Euler which was done in class.

First, we find the local truncation error. Using a Taylor expansion, we seen that

$$y_{k+1} - u_{k+1} = \frac{1}{2} f''(\xi) h^2 = h\tau_k = O(h^2)$$

Second, we have to bound the error globally. If we let e_k = y_k - u_k, we can find out how the errors grow globally. It works similarly to the same way it did in class except for the first line.

$$e_{n+1} = e_n + h[f(y_{n+1},t_{n+1})-f(u_{n+1},t_{n+1})] \\$$

Notice that the derivative is evaluated at step n+1 rather than nfor forward Euler. But, you Taylor expand exactly the same way we did in class about this step. Applying Lipschitz continuity, we get

$$\begin{eqnarray*}(1-hK) \vert e_{N+1}\vert & \leq & \vert e_N\vert + h \tau \\ ... ...rt & \leq & \frac{1}{1-hK}\vert e_N\vert + \frac{1}{1-hK} h \tau \end{eqnarray*}$$

Now, if we work backward in time toward e₀, we get

$$\begin{eqnarray*}\vert e_{N+1}\vert & \leq & \frac{1}{(1-hK)^{T/h}} \vert e_0\ve... ...au \\ & \leq & e^{KT} \vert e_0\vert + \frac{e^{KT}-1}{K} \tau. \end{eqnarray*}$$

This is the same result as for forward Euler... more or less.

6.

[20 pts] Consider the numerical method

$$u_{k+1} = u_k + h\left[\frac{2}{3} f(u_k,t) + \frac{1}{3} f(u_{k+1},t)\right].$$

Using the scalar test equation,

$$y' = -\lambda y,$$

determine the stability of the method when $\lambda>0$. For what values of $h\lambda$ does the numerical solution oscillate? For what values of $h\lambda$ is the method unstable?

This is direct application of the scalar test equation. If we substitute $y' = -\lambda y$ into our numerical method, we obtain the recurrence relation

$$u_{k+1} = u_k - \frac{2}{3} h \lambda u_k - \frac{1}{3} h \lambda u_{k+1}$$

Simplifying, we get

$$u_{k+1} = \left(\frac{1-\frac{2}{3}h\lambda}{1+\frac{1}{3}h\lambda}\right)u_k.$$

The solution to the recurrence relation (when u₀=1) is

$$u_k = \left(\frac{1-\frac{2}{3}h\lambda}{1+\frac{1}{3}h\lambda}\right)^k,$$

so it is just a matter of figuring out when the quantity in parentheses is less than zero or has modulus greater than one. Routine analysis shows that it is negative when $h\lambda>3/2$, so we expect these solutions to oscillate. Similarly, the modulus is greater than 1 when $h \lambda > 6$. At this point, solutions will be unstable.