Math428: Problem Set #2 solution guide

1.

Write the differential equation

$$y'' + (y')^2/y - \sin(y) = 0, \ \ \ \ y(0) = 1, \ \ \ \ y'(0) = 1$$

as a first order system of equations. Note: it will not necessarily be linear. Suggest a general class of $n^{\rm th}$-order nonlinear differential equations that can be converted to a first order system. Try to make your class as broad as possible.

Following the steps laid out in class, we define

$$\begin{eqnarray*}u_0 & = & y \\ u_1 & = & y' \end{eqnarray*}$$

Therefore, the new system becomes

$$\left[ \begin{array}{c} u_0 \\ u_1 \end{array}\right]' = \lef... ...y}\right] = \left[ \begin{array}{c} 1 \\ 1 \end{array}\right].$$

In generality, we see that it would be easy to perform this procedure on any nonlinear system of the form

$$y^{(n)} + f(y^{(n-1)},y^{(n-2)},\ldots,y,t) = 0,$$

by the same procedure. However, we can get into trouble if we cannot solve explicitly for y⁽ⁿ⁾. For instance, consider

$$(y^{(n)})^5 + (y^{(n)})^5 +(y^{(n)})^4 +(y^{(n)})^3 +(y^{(n)})^2 (y^{(n)}) + f(y^{(n-1)},y^{(n-2)},\ldots,y,t) = 0.$$

Since it is impossible to find all the roots of a quintic analytically, our little trick will not work nicely here. Speaking in generalities, if we had a system of the form,

$$g(y^{(n)}) + f(y^{(n-1)},y^{(n-2)},\ldots,y,t) = 0,$$

we need to be able to invert g, or we might be in trouble. Of course, I look forward to reading about your thoughts.

2.

Using Euler's method, solve

$$y' = -y, \ \ \ \ y(0) = 1$$

for $0 \leq t \leq 10$ with step sizes of h=0.5,0.75,1,1.5,2,2.5. Plot the errors in the numerical computation as a function of time for the different step sizes all on the same graph. (Make sure your plot has a legend.) How does the computed solution compare to the exact solution for different values of h? Would a Taylor series method improve this computation at all?

I generated these plots using matlab. See the course web site for the details of this calculation. Here, e_k = y_k-u_k. I plotted the abosolute value of the error, but this was my choice. On the left, I used linear axes which can make it difficult to see the errors as the solutions decay. On the right, I used log-normal axes. For exponential data, logarithmic plots reduce the function to lines which can make it easy for the eye to follow. Of course, when the error is zero and the logarithm is undefined, the data point is excluded. Notice that the unstable solutions (h > 1) stand out as particularly poor.

3.

Repeat problem #2 using the backward Euler's method. Explain the difference in behavior by examining the numerical method as a difference equation.

Here I used the same procedure as in Problem 2, and you can see similar performance. To do backward Euler, you will need to be able to implicitly u_k+1 - h f(u_k+1,t_k+1). This is easy for the scalar test equation.

4.

Repeat problem #2 for the differential equation

$$y' = -\left[1+\frac{9}{10}\cos(t)\right] y, \ \ \ \ y(0) = 1$$

for $0 \leq t \leq 20$. If you cannot find an exact solution, compute a reference solution by setting h to be very, very small. Explain the difference in behavior for different values of h by examining the varying decay rate of the exact system as a function of time.

Here, the solution can be found by solving it for small h, or exactly via separation.

$$y = \exp \left[-\left(t+\frac{9}{10}\sin t\right)\right].$$

Here, we can think of

$$\lambda \approx 1+\frac{9}{10}\cos t$$

as a variable decay constant. Thus, when $(1+\frac{9}{10}\cos t)h > 1$, we can expect some unstable oscillations. Notice that, at some times, a larger step size works better. For example, we see that the h=2 works better than h=1.5 most of the time. This is the curse of instabilities.