\documentclass[12pt]{article}
\usepackage{fullpage,amsmath,amsfonts}

\newcommand{\x}{\boldsymbol{x}}
\newcommand{\y}{\boldsymbol{y}}
\newcommand{\f}{\boldsymbol{f}}
\renewcommand{\v}{\boldsymbol{v}}
\newcommand{\w}{\boldsymbol{w}}
\renewcommand{\r}{\boldsymbol{r}}
\renewcommand{\a}{\boldsymbol{a}}
\renewcommand{\b}{\boldsymbol{b}}
\renewcommand{\c}{\boldsymbol{c}}
\renewcommand{\k}{\boldsymbol{k}}
\newcommand{\e}{\mathrm{e}}
\renewcommand{\d}{\mathrm{d}}
\renewcommand{\Re}{\operatorname{Re}}
\renewcommand{\Im}{\operatorname{Im}}
\newcommand{\Range}{\operatorname{Range}}
\newcommand{\Ker}{\operatorname{Ker}}
\newcommand{\Span}{\operatorname{Span}}
\renewcommand{\dim}{\operatorname{dim}}
\newcommand{\R}{\mathbb{R}}
\newcommand{\C}{\mathbb{C}}

\def\nl{\vphantom\vdots\\}

\begin{document}

\title{Numerical Methods I}
\author{Problem Set 10}
\date{due in class, December 3, 2003}
\maketitle


\begin{enumerate}


\item \textbf{Project:} Implement each of the following solvers as an
\textsf{Octave} function.  Each solver should take arguments of the
form
\begin{verbatim}
  ode_rk4 ('f', [t0,t1], y0, N)
\end{verbatim}
where \texttt{f} is the function $f(t, y)$ on the right hand side of
the differential equation, \texttt{t0} is the initial time,
\texttt{t1} is the final time, \texttt{y0} the initial condition, and
\texttt{N} the number of steps.

\begin{enumerate}
\item \verb+ode_midpoint+, the implicit midpoint method from Lab~8,
where
\[
  \y_{n+1} = \y_n + h \, \f \Bigl(\frac{t_n+t_{n+1}}2,
             \frac{\y_n+\y_{n+1}}2 \Bigr) \,.
\]
Note: Use simple iteration to solve for $\y_{n+1}$.

\item \verb+ode_rk4+, one of the fourth order explicit Runge--Kutta
methods given by
\begin{align*}
  \k_1 & = \f(t_n, \y_n) \,, \\
  \k_2 & = \f(t_n + \tfrac12 h, \y_n + \tfrac12 h \, \k_1) \,, \\
  \k_3 & = \f(t_n + \tfrac12 h, \y_n + \tfrac12 h \, \k_2) \,, \\
  \k_4 & = \f(t_{n+1}, \y_n + h \, \k_3) \,, \\
  \y_{n+1} & = y_n + \tfrac16 h \, (\k_1 + 2\,\k_2 + 2\,\k_3 + \k_4) \,. 
\end{align*}

\item \verb+ode_bdf4+, the fourth order backward differentiation
formula
\[
  25 \, \y_{n+4} - 48 \, \y_{n+3} + 36 \, \y_{n+2} - 16 \, \y_{n+1}
  + 3 \, \y_n = 12 h \, \f(t_{n+4}, \y_{n+4}) \,.
\]
Note: Use \verb+ode_rk4+ to generate the required three extra initial
values.
\end{enumerate}

\item \textbf{Project:} Verify, as on the previous homework, the order
of each of the integrators on a model problem.

\item Analyze the linear stability of critical points of the
Volterra--Lotka model
\begin{align*}
  \dot x = x \, (y-1) \equiv f (x,y) \,, \\
  \dot y = y \, (1-x) \equiv g (x,y) \,,
\end{align*}
as follows.
\begin{enumerate}
\item Find the critical points, i.e.\ the zeros of the right hand side
of the equation.  (There are two critical points in this case.)

\item Compute the Jacobi Matrix of the right hand side, i.e., the
matrix of partial derivatives
\[
DF =
\begin{pmatrix}
\frac{\partial f}{\partial x} & \frac{\partial f}{\partial y} \\
\frac{\partial g}{\partial x} & \frac{\partial g}{\partial y}
\vphantom{\Bigl(}
\end{pmatrix} \,.
\]

\item Find the eigenvalues of $DF$ evaluated at each of the critical
points.  If all the eigenvalues have negative real part, then the
critical point is \emph{linearly stable}, if at least one eigenvalue
has positive real part, the critical point is \emph{linearly
unstable}, if all the eigenvalues have zero real part, the critical
point is \emph{linearly neutrally stable}.

(You should fine one linearly neutrally stable and one linearly
unstable critical point.)

\end{enumerate}

\item \textbf{Project:} It can be shown that the critical point at
$(1,1)$ of the Volterra--Lotka model is actually stable in the sense
that if a solution starts near this point, it will always remain in a
small neighborhood of this point.

Test the three numerical methods you have implemented to see if they
retain this property.

\item \textbf{Project:} (See QS, Section 7.9) Solve the Van der Pol
equation
\begin{align*}
  \dot x & = y \,, \\
  \dot y & = \mu \, (1-x^2) \, y - x \,,
\end{align*}
with initial data $x (0) = 2$ and $y(0) = 0$.  Which of the methods
works best when $\mu$ is large?

\end{enumerate}

\end{document}
