\documentclass[12pt]{article}
% These packages help with equations, math symbols, etc.
\usepackage{amsmath,amssymb}
% This EXCELLENT resource is your guide to writing equations in LaTeX!
% ftp://ftp.ams.org/pub/tex/doc/amsmath/short-math-guide.pdf
% The following gives us better margins and spacing
\usepackage[margin=1in]{geometry}
\setlength{\topskip}{0in} % No blank space on first page
\setlength{\parindent}{0in} % 0 = no paragraph indentation
\setlength{\parskip}{0.5em} % space between paragraphs
%% remove ligatures (e.g. ``ff") to allow text copying from PDF
\usepackage{microtype}
\DisableLigatures{encoding = *, family = * }
% In case you want to hyperlink anything...
\usepackage[colorlinks=true, urlcolor=blue, anchorcolor=black, linkcolor=black, citecolor=black]{hyperref}
% For more compact numbered lists...
\usepackage{enumitem}
\setlist{noitemsep, topsep=-0.5em}
%% DRAFT watermark across each page:
%\usepackage[printwatermark]{xwatermark}
%\usepackage{xcolor}
%\newwatermark[allpages,color=gray!30,angle=55,scale=6,xpos=0,ypos=0]{DRAFT}
% If you want to include an image
\usepackage{graphicx}
\begin{document}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% Configure some details for knitr...
%<>=
%# smaller font size for chunks
%knitr::opts_chunk$set(size = 'footnotesize', concordance=TRUE)
%@
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\textbf{\large MATH 420 Homework \#5 \hfill Due Mon, 4 Dec 2017} % Replace Due date with your name.
\vspace{1pt}
\hrule
\vspace{0.5em} % this just forces a blank line.
%\textbf{Instructions:} A printed copy of your homework should be handed in at \textbf{the start of class} the day it is due. Supplementary electronic files (e.g. R scripts or wxMaxima files;) should be emailed to the instructor prior to class with file name format: \texttt{Lastname-hwX.ext} \\[0.5em]
Recall that for 1D systems $\dot{x}=f(x)$, the local asymptotic stability of an equilibrium point $x_*$ ($f(x_*)=0$) is determined by the slope of $f(x)$ at that point: if $f'(x_*)$ is negative then $x_*$ is LAS, but if $f'(x_*)$ is positive then $x_*$ is unstable.
For systems of multiple equations, where $\dot{\mathbf{x}}=f(\mathbf{x})$ for $\mathbf{x}\in\mathbb{R}^n$ and $f:\mathbb{R}^n \to \mathbb{R}^n$ is differentiable, we determine the stability of an equilibrium point $\mathbf{x}_*=(x_{1*}, x_{2*},\ldots,x_{n*})$ by determining whether or not all of the eigenvalues $\lambda_i$ of the Jacobian Matrix $\mathbf{J}=(\partial f_i(\mathbf{x})/\partial x_j)$ have negative real part. In practice, we do this using computer algebra systems like Maxima, Mathematica or Maple to check the Routh-Hurwitz Criteria which allows us to assess whether or not all $\lambda_i$ have negative real part (or not) without explicitly calculating the eigenvalues.
The Routh-Hurwitz Criteria are a set of necessary and sufficient conditions for whether or not the roots of a polynomial have negative real part. You may recall that the \textbf{characteristic equation} of an $n\,$x$\,n$ matrix $\mathbf{A}$ is the $n^\text{th}$-order polynomial $p(x)=\det(A-x\,I)$, and its roots are by definition the eigenvalues of $\mathbf{A}$. \\
\fbox{\begin{minipage}{0.91\textwidth}
\centerline{\textbf{Routh-Hurwitz Criteria}}
\vspace{0.5em}
All roots of the polynomial (with real coefficients $c_i$) $$p(x)=c_n + c_{n-1}\,x + \cdots + c_{1}\,x^{n-1} + x^n$$ have negative real parts if and only if the determinants of all the corresponding Hurwitz matrices are positive. This result provides an algorithm for computing stability criteria, which gives these equivalent conditions for small values of $n$: \begin{align*}
n=2 \quad & \quad c_i>0\; \\
n=3 \quad & \quad c_i>0,\; c_1\,c_2\, > \,c_3\\
n=4 \quad & \quad c_i>0,\; c_1\,c_2\,c_3 > \, c_3^2 + c_1^2\,c_4
\end{align*}
\textbf{Further reading:} See \href{http://wwwhome.math.utwente.nl/~meinsmag/reports/routh.html}{Meinsma (1995)}, Gantmacher (1989), or Ch. 4 of \textit{Introduction to Mathematical Biology} by Allen (2007).
\end{minipage} }
~\\[2em]
\textbf{Exercise 1:} Find all equilibria of the following 1D models, and determine their stability (and/or state any parameter conditions that would yield a stable eq. pt.) by calculating $f'(x_*)$ and inferring whether or not the equilibrium point is locally asymptotically stable.\\
\begin{enumerate}
\item $\dot{x}=1-a\,x$ \\
\item $\dot{x}=(1-a\,x)(x-\epsilon)$ \\
\item $\dot{x}=\cos(x)$ \\
\end{enumerate}
\clearpage
\textbf{Exercise 2:} Find all equilibria for the following system of equations, and determine each of their stability, by following the steps below. Assume $a,r>0$ and $\theta>1$ $$\dot{x}=f(x,y)=1-r\,x$$ $$\dot{y}=g(x,y)=r\,x - \lambda\,y^\theta$$
~\\
\begin{enumerate}
\item Solve $\dot{x}=0$ and $\dot{y}=0$ for $x$ and $y$.\\
\item Next we construct the Jacobian matrix by hand. \begin{enumerate}
\item Find the following partial derivatives:\\
$$\frac{\partial f}{\partial x}=\hspace{10cm}$$
$$\frac{\partial f}{\partial y}=\hspace{10cm}$$
$$\frac{\partial g}{\partial x}=\hspace{10cm}$$
$$\frac{\partial g}{\partial y}=\hspace{10cm}$$
~\\
\item Find the Jacobian at each equilibrium point by evaluating $$\mathbf{J}=\begin{bmatrix}
\frac{\partial f}{\partial x} & \frac{\partial f}{\partial y} \\
\frac{\partial g}{\partial x} & \frac{\partial g}{\partial y}
\end{bmatrix}$$ at each equilibrium point. Simplify as much as possible.
\end{enumerate}
\item Either apply the Routh-Hurwitz criteria to determine equilibrium stability, or use the following shortcut for 2D systems:
If $\mathbf{A}=\begin{bmatrix} a & b\\ c& d \end{bmatrix}$ then it's Trace ($Tr(\mathbf{A})=a+d$) and Determinant ($Det(\mathbf{A})=a\,d-c\,b$) are related to it's eigenvalues $\lambda_i$ by $Tr(\mathbf{A})=\lambda_1+\lambda_2$ and $Det(\mathbf{A})=\lambda_1\,\lambda_2$. Therefore the eigenvalues have negative real part iff $$Tr(\mathbf{A})<0 \qquad \text{and} \qquad Det(\mathbf{A})>0.$$
\end{enumerate}
\clearpage
\textbf{Exercise 3: } The Lorenz equations are a simplification of a fluid dynamics model, that were derived to illustrate chaotic dynamics by \href{http://dx.doi.org/10.1175/1520-0469(1963)020<0130:DNF>2.0.CO;2}{Ed Lorenz in 1963}. \begin{subequations}\label{lorenz}
\begin{align}
\dot{x} = & \sigma \, (y-x)\\
\dot{y} = & r\,x-y-x\,z\\
\dot{z} = & x\,y-b\,z
\end{align} \end{subequations}
The Routh-Hurwitz criteria can be used to find stability conditions (i.e., conditions on the values of parameters $\sigma>0$, $b>0$, $r>1$) for the equilibrium point
$$(x_*,\;y_*,\;z_*)=\left(\sqrt{b\,r-b},\;\sqrt{b\,r-b},\;r-1\right )$$ as follows:\\
The Jacobian for a general point $(x,y,z)$ for equations (1a-1c) is
\[ \mathbf{J} = \begin{bmatrix} -\sigma & \sigma & 0\cr r-z & -1 & -x\cr y & x & -b\end{bmatrix} \]
By evaluating the Jacobian at the given equilibrium point ($x_*,y_*,z_*)$ then finding its \textit{characteristic polynomial equation} we can apply the Routh-Hurwitz criteria.\\
\[ \mathbf{J_*} = \begin{bmatrix}-\sigma & \sigma & 0\cr 1 & -1 & -\sqrt{b\cdot r-b}\cr \sqrt{b\cdot r-b} & \sqrt{b\cdot r-b} & -b\end{bmatrix} \]
% Feel free to indicate c_i using an \overbrace^{} or \underbrace_{}:
% $p(x)= x^3+\overbrace{(stuff)}^{c_1}\,x^2+...$\\
Finding the characteristic polynomial and dividing by the leading coefficient (so that the $x^3$ term has coefficient 1) yields $$ {{x}^{3}}+\left( \sigma+b+1\right) \cdot {{x}^{2}}+\left( b\cdot \sigma+b\cdot r\right) \cdot x+2\cdot b\cdot r\cdot \sigma-2\cdot b\cdot \sigma = 0 $$
Therefore, $c_3=2b\,r\,\sigma-2b\,\sigma$, $c_2=b\,\sigma+b\,r$ and $c_1=\sigma+b+1$ applying the criteria above, the Routh-Hurwitz criteria are that each coefficient is positive (which is always true for the given constraints $r>1$, $b,\sigma>0$) and that
$$(\sigma+b+1)(\sigma+r)b > 2b\,\sigma\,(r-1)$$
which reduces to
$$(\sigma+b+1)(\sigma+r) > 2\sigma\,(r-1).$$
Hence, the given equilibrium point is stable whenever (1) it exists in the positive real orthant (i.e., $x_*$, $y_*$, and $z_*$ are all positive real numbers) and (2) whenever $$(\sigma+b+1)(\sigma+r) > 2\sigma\,(r-1).$$
Interestingly, when this equilibrium point is \textit{unstable}, these equations can exhibit \textbf{chaos}: dynamics characterized by two trajectories near the asymptotic attractor will eventually diverge regardless of how close their initial conditions. In this exercise, you will modify some R code to illustrate this divergence of trajectories near the chaotic attractor in the Lorenz system.
\textbf{(1)} Modify the code \textbf{lorenz.R} on the website to create a 3x3 panel of figures: The first column are the $x$, $y$, and $z$ trajectories given in the code. The second column is the same three curves but for a different numerical solution of the model that begins at slightly different initial conditions: $Y0=c(x=10.01, y=10.99, z=12.01$). The third column is the two state-space plots (two versions of the scatterplot3d plot in the code, one for each trajectory) followed by a plot that shows the Euclidian distance between the two trajectories at each time point, illustrating whether or not the trajectories diverge over time.
Submit a hard copy of this figure, and an electronic copy of your code, to receive full points for this problem. Please feel free to change colors, etc. or otherwise make your code and figure unique so that it doesn't look like anyone else's in class.
\end{document}