\documentclass[12pt,a4paper,titlepage]{article} %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% \usepackage{amsmath} \usepackage{amsmath} \usepackage{amsmath} \usepackage{amsmath} \usepackage{amsmath} \usepackage{amsmath} \usepackage{amsmath} \usepackage{amsmath} \usepackage{amsmath} \usepackage{amsmath} \usepackage{amsmath} \usepackage{amsmath} \usepackage{amsmath} \usepackage{amsmath} \usepackage{amsmath} \usepackage{amsmath} \usepackage{amsmath} \usepackage{amsmath} \usepackage{amsmath} \usepackage{amsmath} \usepackage{amsmath} \usepackage{amsmath} \usepackage{amsmath} \usepackage{amsmath} \usepackage{amsmath} \usepackage{amsmath} \usepackage{amsmath} \usepackage{amsmath} \usepackage{amsmath} \usepackage{amsmath} \usepackage{graphicx} \usepackage{amsmath} \usepackage{harvard} \setcounter{MaxMatrixCols}{10} %TCIDATA{OutputFilter=LATEX.DLL} %TCIDATA{Version=5.50.0.2960} %TCIDATA{} %TCIDATA{BibliographyScheme=Manual} %TCIDATA{Created=Mon Nov 05 11:50:44 2001} %TCIDATA{LastRevised=Tuesday, December 14, 2010 09:51:23} %TCIDATA{} %TCIDATA{} %TCIDATA{Language=American English} %TCIDATA{CSTFile=LaTeX article (bright).cst} %TCIDATA{ComputeDefs= %$\delta =.2$ %$h=2$ %$q=1/5$ %1$s_{e}=1$ %1$s_{u}=-2$ %$\alpha =1/2$ %} \newtheorem{theorem}{Theorem} \newtheorem{acknowledgement}[theorem]{Acknowledgement} \newtheorem{algorithm}[theorem]{Algorithm} \newtheorem{axiom}[theorem]{Axiom} \newtheorem{case}[theorem]{Case} \newtheorem{claim}[theorem]{Claim} \newtheorem{conclusion}[theorem]{Conclusion} \newtheorem{condition}[theorem]{Condition} \newtheorem{conjecture}[theorem]{Conjecture} \newtheorem{corollary}[theorem]{Corollary} \newtheorem{criterion}[theorem]{Criterion} \newtheorem{definition}{Definition} \newtheorem{example}[theorem]{Example} \newtheorem{exercise}[theorem]{Exercise} \newtheorem{lemma}[theorem]{Lemma} \newtheorem{notation}[theorem]{Notation} \newtheorem{problem}[theorem]{Problem} \newtheorem{proposition}[theorem]{Proposition} \newtheorem{remark}[theorem]{Remark} \newtheorem{result}[theorem]{Result} \newtheorem{solution}[theorem]{Solution} \newtheorem{summary}[theorem]{Summary} \newenvironment{proof}[Proof]{\textbf{#1.} }{\ \rule{0.5em}{0.5em}} \input{tcilatex} \begin{document} \title{MathII at SU and SSE} \author{John Hassler \\ %EndAName IIES, Stockholm University} \maketitle \section{Introduction\protect\footnote{% I am grateful for comments and corrections provided by Dirk Niepelt.}} This course is about dynamic systems, i.e., systems that evolve over time. The analysis of the dynamic evolvement of economic systems is of core importance in many areas of economics. Growth, business cycles, asset pricing, and dynamic game theory are just a few examples. \subsection{Solving a simple dynamic system} Very often, our economic models provide a difference or differential equation for the endogenous variables. Take a very simple example; an arbitrage asset pricing model. Suppose there is a safe asset, a bond, that provides a return $r$ every period. Suppose a share, giving rights to dividend flow of $d$ per period$,$ is introduced. Now, arbitrage theory says that the share should also yield a return $r$ in equilibrium. Defining the price on the share as $p$, arbitrage theory thus implies, \begin{equation} \frac{p_{t+1}+d}{p_{t}}=1+r. \end{equation} This is a simple difference equation, \begin{equation} p_{t+1}=\left( 1+r\right) p_{t}-d. \end{equation} Note that we $p_{t}$ is an \emph{endogeneous} variable in this system while $% r$ and $d$ are \emph{exogeneous}. This would be the case also if $d$ and $r$ varies over time, cases we will study below. One straightforward way of solving it is to substitute forward or backward, e.g., noting that \begin{eqnarray} p_{t} &=&\left( 1+r\right) p_{t-1}-d \\ &=&\left( 1+r\right) \left( \left( 1+r\right) p_{t-2}-d\right) -d \\ &=&\left( 1+r\right) ^{2}p_{t-2}-d\left( 1+\left( 1+r\right) \right) \end{eqnarray}% and so on. A more general approach is to first characterize all possible paths consistent with the law-of-motion. Here, this is quite simple. You will learn that set of possible paths is \begin{equation} p_{t}=c\left( 1+r\right) ^{t}+\frac{d}{r} \end{equation}% for any constant $c.$ As we see, there is an infinite number of solutions, i.e., we need more information. If, for example, we know that the value of $% p_{0}$, we can solve for the constant \begin{eqnarray} c &=&p_{0}-\frac{d}{r} \\ &\rightarrow &p_{t}=\left( p_{0}-\frac{d}{r}\right) \left( 1+r\right) ^{t}+% \frac{d}{r}. \end{eqnarray} In finance, the solution $p=\frac{d}{r}$ is called the fundamental solution, and we see that if $r>0,$ the solution explodes as $t$ goes to infinity if $% p_{0}\neq \frac{d}{r}.$ This gives us another way of solving the difference equation. Suppose, we have reason to believe that the solution should remain bounded. Then, if $r>0$, the only solution left is when $c=0.$ Note that $% \left( 1+r\right)$ is called the root of the system. Note, the importance of whether the root is bigger or smaller than unity (in absolute values). We will also work in continuous time. Then,\ I usually use put the time variable in parenthesis and use the dot symbol to indicate time derivatives. In continuous time, non-existence of arbitrage means that capital gains, i.e., the change in the price per unit of time plus dividends per unit of time should equal the opportunity cost, i.e., interest rate on the price of the assets. Thus, non-existence of arbitrage implies \begin{equation} \dot{p}\left( t\right) +d=rp\left( t\right) . \end{equation}% This is a simple linear first-order differential equation. The set if solution is:% \begin{equation} p\left( t\right) =ce^{rt}+\frac{d}{r}. \end{equation} Also here, we have a term that is explosive if $r>0.$ Later in the course, we will learn how to solve more complicated dynamic systems, involving, e.g., several endogenous variables and varying parameters. \subsection{Two approaches to Dynamic Optimization} The second part of the course, is to solve maximization problems in dynamic systems. Suppose there is a potentially infinite set of paths $x\left( t\right) _{0}^{T},$ each denoting a particular continuous function $x\left( t\right)$ for $t\in \left[ 0,T\right] .$ Suppose also that we can evaluate them, i.e., they give different payoffs. Then, we will learn how to derive difference or differential equations, that are necessarily satisfied for optimal paths. If we can solve these equations, we can find the optimal path. We will use two approaches in this course. \subsubsection{Dynamic Programming (Bellman).} Suppose we want to find an optimal investment plan in discrete time and let $% x_{t}$ denote the stock of capital at time $t.$ Also, let $u_{t}$ denote investments and assume \begin{equation} x_{t+1}=g(x_{t},u_{t}), \end{equation}% which we call the law-of-motion for $x_{t}.$ Each period, the payoff is given by $F\left( t,x_{t},u_{t}\right)$ and the problem is to solve \begin{eqnarray} &&\max_{x_{0}^{T+1},u_{0}^{T}}\sum_{t=0}^{T}F\left( t,x_{t},u_{t}\right) , \\ \text{ s.t. }x_{t+1} &=&g(x_{t},u_{t})\forall t, \\ x_{0} &=&x. \end{eqnarray} Note that this is a dynamic problem; The choice of investment at time $t,$ $% u_{t}$ may affect payoffs in many future period. First, the payoff at $t,$ $% F\left( t,x_{t},u_{t}\right)$ is affected directly and also the next periods payoff since \begin{equation*} F\left( t+1,x_{t+1},u_{t+1}\right) =F\left( t+1,g(x_{t},u_{t}),u_{t+1}\right) . \end{equation*} Furthermore, also payoffs further away can be affected since, for example, $% x_{t+2}=$ $g(x_{t+1},u_{t+1})=g(g(x_{t},u_{t}),u_{t+1})$, affecting the payoff in period $t+2.$ The choice of $u_{t}$ must thus take into account all future payoff relevant effects. Sometimes the dynamic problem degenerates in the sense that this dynamic link breaks. To illustrate this, let us simplify and take the example where $% F\left( t,x_{t},u_{t}\right) =\left( \frac{1}{1+r}\right) ^{t}\left( f_{t}\left( x_{t}\right) -u_{t}\right)$ and $g\left( x_{t},u_{t}\right) =\left( 1-\delta \right) x_{t}-u_{t}.$ We can interpret $\delta$ as the rate of capital depreciation. If $\delta =1$, we have $x_{t+1}=u_{t}$ and by substituting from the law-of-motion, we can write the problem as% \begin{eqnarray} &&\max_{x_{1}^{T}}\sum_{t=0}^{T}\left( \frac{1}{1+r}\right) ^{t}\left( f_{t}\left( x_{t}\right) -x_{t+1}\right) , \\ x_{0} &=&x. \end{eqnarray} The first-order condition for choosing $x_{t}$ for any $t>0$ is $% f_{t}^{\prime }\left( x_{s}\right) =1+r$, so to know how much to invest in period $t,$ we only need to know the marginal productivity of capital next period, i.e., we maximize period-by-period. With $\delta <1,$ we cannot do this. Then dynamic programing is a handy way to attack the problem. We will use the Bellman's principle of optimality, saying that there exists a sequence of functions $V_{t}\left( x_{t}\right)$ such that \begin{equation} V_{t}\left( x_{t}\right) =\max_{u_{t}}F\left( t,x_{t},u_{t}\right) +V_{t+1}\left( g(x_{t},u_{t})\right) \label{Bellman} \end{equation}% and where% \begin{eqnarray} V_{t}\left( x_{t}\right) &\equiv &\max_{x_{t}^{T+1},u_{t}^{T}}\sum_{s=t}^{T}F\left( s,x_{s},u_{s}\right) \\ \text{ s.t. }x_{t+1} &=&g(x_{t},u_{t}) \\ &&x_{t}\text{ given.} \end{eqnarray} $V_{t}\left( x_{t}\right)$ is called the value function and $x_{t}$ a state variable. We find that through the Bellman principle, we have split up the problem into a sequence of one period problems. We can then solve the problem more easily since the first-order condition for maximizing the RHS of (\ref{Bellman}) \begin{equation} F_{u}\left( t,x_{t},u_{t}\right) +V_{t+1}^{\prime }\left( g(x_{t},u_{t})\right) g_{u}(x_{t},u_{t})=0. \end{equation} implicitly yields a system of difference equations in $x_{t}$ and $u_{t}$ that we may be able to solve. \subsubsection{Optimal control (Pontryagin)} The other way of solving dynamic optimization problems that we will use is called Optimal Control. We will use it for continuous time problems. Suppose we want to solve \begin{eqnarray} &&\max_{x\left( t\right) _{0}^{T},u\left( t\right) _{0}^{T}}\int_{0}^{T}F\left( t,x\left( t\right) ,u\left( t\right) \right) dt \\ s.t.\,\dot{x}\left( t\right) &=&g\left( x\left( t\right) ,u\left( t\right) \right) \\ x\left( 0\right) &=&x_{0} \\ x\left( T\right) &=&x_{T}. \end{eqnarray} Then, Pontryagin's maximum principle says that for each point in time, the optimal control, call it $u^{\ast }\left( t\right) ,$ satisfies \begin{equation} u^{\ast }\left( t\right) =\arg \max_{u\left( t\right) }F\left( t,x\left( t\right) ,u\left( t\right) \right) +\lambda \left( t\right) g\left( x\left( t\right) ,u\left( t\right) \right) , \label{Pontryagin} \end{equation} where $\lambda \left( t\right)$ can be interpreted as the shadow value of the state variable (capital). The sum $F\left( t,x\left( t\right) ,u\left( t\right) \right) +\lambda \left( t\right) g\left( x\left( t\right) ,u\left( t\right) \right)$, is called the Hamiltonian. Again, we have turned the dynamic problem into a sequence of static problems. The first order condition for (\ref{Pontryagin}) will implicitly define a system of differential equations. Note the similarity between (\ref{Pontryagin}) and (% \ref{Bellman}).\newpage \section{Some basics} \subsection{Taylor series} We will often need to approximate a continuous and differentiable function around some point. Specifically, suppose we know $f\left( x_{0}\right)$ and some of its derivatives and want to approximate the value in some other point $x\neq x_{0}$. An efficient way of doing such an approximation is to use Taylors formula. This can be seen as an attempt to fit a polynomial (we will talk about such below) to a curve. \begin{result} The Taylor approximation is% \begin{eqnarray} f\left( x\right) &\approx &f\left( x_{0}\right) +\frac{1}{1!}f^{\prime }\left( x_{0}\right) \left( x-x_{0}\right) +\frac{1}{2!}f^{\prime \prime }\left( x_{0}\right) \left( x-x_{0}\right) ^{2} \\ &&+\frac{1}{3!}f^{\prime \prime \prime }\left( x_{0}\right) \left( x-x_{0}\right) ^{3}+... \notag \end{eqnarray} \end{result} Usually we will just use the first order approximation $f\left( x\right) \approx f\left( x_{0}\right) +f^{\prime }\left( x_{0}\right) \left( x-x_{0}\right)$ but sometimes a higher order approximation can be useful. It can furthermore be shown, that the approximation error of an n'th order Taylor approximation is given by \begin{equation} \frac{1}{\left( n+1\right) !}\frac{d^{\left( n+1\right) }f\left( c\right) }{% dx^{\left( n+1\right) }}\left( x-x_{0}\right) ^{n+1} \end{equation}% where $c\in \lbrack x_{0},x].$ As we see, when $x-x_{0}$ is small, two forces imply that the approximation error tends to be small when $n$ is large. In the denominator, $\left( n+1\right) !$ is large and $\left( x-x_{0}\right) ^{n+1}$ is small. Although higher-order approximations often are better, we cannot say that for any $\left( x-x_{0}\right)$ a higher order approximation is better. \subsection{Integration\label{integration}} %TCIMACRO{\TeXButton{TeX field}{\setcounter{equation}{0}}}% %BeginExpansion \setcounter{equation}{0}% %EndExpansion If $b>a$, the expression \begin{equation} \int\limits_{a}^{b}f\left( t\right) dt, \label{integral1} \end{equation}% can be interpreted as the area under the graph $y=f(x)$ for $x\in \lbrack a,b]$. How should one compute such an area? The most natural way would be to divide the interval $[a,b]$ into (many) sub-intervals by choosing numbers $% a=t_{1}0)$ . Here this gives an end condition, since it requires that \begin{equation} \lim_{t\rightarrow \infty }e^{-rt}y\left( t\right) =0. \end{equation}% Using this, we get \begin{equation} c=\lim_{t\rightarrow \infty }\int_{t_{0}}^{t}e^{-rs}q\left( s\right) ds\equiv \int_{t_{0}}^{\infty }e^{-rs}q\left( s\right) ds, \end{equation}% giving \begin{eqnarray} e^{-rt}y\left( t\right) &=&\int_{t_{0}}^{\infty }e^{-rs}q\left( s\right) d-\int_{t_{0}}^{t}e^{-rs}q\left( s\right) ds \\ &=&\int_{t}^{\infty }e^{-rs}q\left( s\right) ds \\ y\left( t\right) &=&e^{rt}\int_{t}^{\infty }e^{-rs}q\left( s\right) ds=\int_{t}^{\infty }e^{-r\left( s-t\right) }q\left( s\right) ds, \end{eqnarray}% i.e., that no arbitrage and \textquotedblright No-Ponzi\textquotedblright\ implies that the price must equal the discounted sum of future dividends. \subsubsection{Variable coefficients} If also the coefficient on the endogenous variable ($y\left( t\right)$) is varying over time, we need a more general integrating factor to make the LHS od the equation being the time differential of a known function. For example, consider the no-arbitrage equation for an asset price $y\left( t\right)$ with dividends $q\left( t\right)$ and interest rate $r\left( t\right)$% \begin{equation} \dot{y}\left( t\right) +q\left( t\right) =r\left( t\right) y\left( t\right) \Rightarrow \dot{y}\left( t\right) -r\left( t\right) y\left( t\right) =-q\left( t\right) . \end{equation} Here the integrating factor is \begin{equation} e^{-\int_{t_{0}}^{t}r(s)ds}. \end{equation}% with \begin{equation} \frac{de^{-\int_{t_{0}}^{t}r(s)ds}}{dt}=-r\left( t\right) e^{-\int_{t_{0}}^{t}r(s)ds}, \end{equation}% alternatively expressed as% \begin{equation*} \frac{de^{^{-\int r(t)dt}}}{dt}=-r\left( t\right) e^{-\int r(t)dt}. \end{equation*} Approximating the integral as a sum of rectangles with base $\Delta t$ as in section \ref{integration}, (which is exact if $r(t)$ were piecewise constant), and defining $r\left( t_{0}+s\Delta t\right) \equiv r_{t_{s}}$, for the integers $s\in \left\{ 0,S\right\} ,$ $S=(t-t_{0})/\Delta t$, the integrating factor, can be written \begin{equation} e^{-\int_{t_{0}}^{t}r(s)ds}\approx e^{-r_{0}\Delta t}e^{-r_{1}\Delta t}...e^{-r_{S}\Delta t}\approx \left( \frac{1}{1+r_{1}\Delta t}\right) \left( \frac{1}{1+r_{1}\Delta t}\right) ...\left( \frac{1}{1+r_{S}\Delta t}% \right) , \end{equation}% i.e., it is product of all short run discount factors between $t_{0}$ and $% t.$ To save on notation, this product is denoted as \begin{equation} e^{-\int_{t_{0}}^{t}r(s)ds}\equiv R\left( t;t_{0}\right) , \end{equation}% or if the starting point is suppressed as $R\left( t\right) .$ This has a straightforward interpretation. Suppose the variable discount rate is given by $r\left( s\right)$, then, the the discounted value of a payment $y\left( t\right)$ at $t$ seen from $t_{0}$ is $R\left( t;t_{0}\right) y\left( t\right) .$Using, the integrating factor, we get \begin{eqnarray} R\left( t\right) \left( \dot{y}\left( t\right) -r\left( t\right) y\left( t\right) \right) &=&-R\left( t\right) q\left( t\right) \\ \frac{dR\left( t\right) y\left( t\right) }{dt} &=&-R\left( t\right) q\left( t\right) ,\rightarrow \\ R\left( t\right) y\left( t\right) &=&-\int_{t_{0}}^{t}R\left( s\right) q\left( s\right) ds+c \end{eqnarray} Suppose again, $\lim {}_{t\rightarrow \infty }R\left( t\right) y\left( t\right) =0,$ implying $c=$ $\int_{t_{0}}^{\infty }R\left( s\right) q\left( s\right) ds$. Then, noting that $R\left( t\right) ^{-1}R\left( s\right) =e^{\int_{t_{0}}^{t}r(v)dv-\int_{t_{0}}^{s}r(v)dv}=e^{-\int_{t}^{s}r(v)dv}$ $% =R\left( s;t\right)$ we have \begin{eqnarray} R\left( t\right) y\left( t\right) &=&-\int_{t_{0}}^{t}R\left( s\right) q\left( s\right) ds+\int_{t_{0}}^{\infty }R\left( s\right) q\left( s\right) ds=\int_{t}^{\infty }R\left( s\right) q\left( s\right) ds \\ y\left( t\right) &=&\int_{t}^{\infty }R\left( t\right) ^{-1}R\left( s\right) q\left( s\right) ds=\int_{t}^{\infty }e^{-\int_{t}^{s}r(v)dv}q\left( s\right) ds, \end{eqnarray}% i.e., that the asset price equals the discounted value of future dividends. As another example, consider money on a bank account with variable interest rate and deposits $q\left( t\right) ,$then \begin{eqnarray} \dot{y}\left( t\right) &=&r\left( t\right) y\left( t\right) +q\left( t\right) \\ \dot{y}\left( t\right) -r\left( t\right) y\left( t\right) &=&q\left( t\right) \\ \frac{dR\left( t,t_{0}\right) y\left( t\right) }{dt} &=&R\left( t,t_{0}\right) q\left( t\right) \\ &\rightarrow &R\left( t,t_{0}\right) y\left( t\right) =\int_{t_{0}}^{t}R\left( s,t_{0}\right) q\left( s\right) ds+c \\ y\left( t\right) &=&\int_{t_{0}}^{t}R\left( t,t_{0}\right) ^{-1}R\left( s,t_{0}\right) q\left( s\right) +R\left( t,t_{0}\right) ^{-1}c \end{eqnarray} Since $s\leq t$ here, the term $R\left( t,t_{0}\right) ^{-1}R\left( s,t_{0}\right)$ is more conveniently\footnote{% But remember $-\int_{t}^{s}r(v)dv=\int_{s}^{t}r(v)dv.$} written $% e^{\int_{s}^{t}r(v)dv}.$ Using also an initial condition, $y\left( t_{0}\right) ,$we have \begin{equation} y\left( t\right) =\int_{t_{0}}^{t}e^{\int_{s}^{t}r(v)dv}q\left( s\right) ds+e^{\int_{t_{0}}^{t}r(v)dv}y\left( t_{0}\right) , \end{equation}% where the first term is the period $t$ value of all deposits up until $t$ and the second term is the period $t$ value of the initial amount on the bank. \subsubsection{Separating variables} Sometimes, we can write a differential equation such that the LHS only contains functions of $x$ and $\dot{x}$ and the RHS only a function of $t$. For example \begin{eqnarray} \dot{x}\left( t\right) &=&\frac{h\left( t\right) }{g\left( x\right) } \\ \dot{x}\left( t\right) g\left( x\right) &=&h\left( t\right) . \end{eqnarray} In this case, we can use the following trick. Let $G\left( x\right)$ be any primitive of $g\left( x\right)$, i.e., $\frac{dG\left( x\right) }{dx}% =g\left( x\right) .$ Then, \begin{eqnarray} \dot{x}\left( t\right) g\left( x\left( t\right) \right) &=&\frac{dG(x\left( t\right) )}{dt}=h\left( t\right) \\ G(x) &=&\int h\left( s\right) ds+c \end{eqnarray} We can then recover $x$ by inverting $G$. An example; \begin{eqnarray} \dot{x}\left( t\right) &=&\left( x\left( t\right) t\right) ^{2} \\ \dot{x}\left( t\right) x\left( t\right) ^{-2} &=&t^{2}. \end{eqnarray}% Now let $g\left( x\right) =x^{-2}$ implying $G\left( x\right) =-x^{-1}$. Then since $\frac{dG\left( x\left( t\right) \right) }{dt}=\dot{x}\left( t\right) g\left( x\left( t\right) \right) =\dot{x}\left( t\right) x\left( t\right) ^{-2}$, we have \begin{eqnarray} \frac{dG\left( x\left( t\right) \right) }{dt} &=&t^{2} \\ G\left( x\left( t\right) \right) &=&-x\left( t\right) ^{-1}=\int t^{2}dt+c=% \frac{t^{3}}{3}+c. \end{eqnarray}% Using the fact that $z=G\left( x\right) =-x^{-1}\Rightarrow x=-z^{-1}$, so the inverse function is given by $G\left( z\right) =-z^{-1}$, we get \begin{equation*} x\left( t\right) =G^{-1}\left( \frac{t^{3}}{3}+c\right) =-\left( \frac{t^{3}% }{3}+c\right) ^{-1}. \end{equation*} \subsection{Linear differential equations of higher order} \subsubsection{Linear second order differential equations} A linear second order differential equation has the form \begin{equation} \ddot{y}\left( t\right) +p\left( t\right) \dot{y}\left( t\right) +q\left( t\right) y\left( t\right) =R\left( t\right) \end{equation} This cannot be solved directly by transformations in the simple way we did with first order equations. Instead we use a more general method (that would have worked above also). First some definitions \begin{definition} The homogeneous part of differential equation is achieved by setting all exogenous variables, including constants, to zero. \end{definition} \begin{definition} Two functions $y_{1}(t)$and $y_{2}(t)$ are linearly independent in a region $% \Omega$ iff there is no $c_{1},c_{2}\neq \{0,0\}$ s.t. \begin{equation} c_{1}y_{1}(t)=c_{2}y_{2}(t)\text{ }\forall t\in \Omega . \end{equation} \end{definition} \begin{result} $\label{solutiondiff1}$The general solution (the complete set of solutions) to a differential equation is the general solution to the homogeneous part plus any particular solution to the complete equation. The general solution of the homogeneous part of a second order linear differential equation can be expressed as $c_{1}y_{1}(t)+c_{2}y_{2}(t)$, where $y_{1}(t)$ and $y_{2}(t)$ are two linearly independent particular solutions to the homogeneous equations and $c_{1}$ and $c_{2}$ are some arbitrary integration constants. \end{result} \subsubsection{Homogeneous Equations with Constant Coefficients} Consider the homogeneous part of a differential equation, given by \begin{equation} \ddot{y}\left( t\right) +p\dot{y}\left( t\right) +qy\left( t\right) =0. \end{equation} To solve this equation, we first define the \emph{characteristic equation, }% for this equation; \begin{equation} r^{2}+pr+q=0. \end{equation} Since this is a second-order polynomial, it has two roots \begin{equation} r_{1,2}=-\frac{1}{2}p\pm \frac{1}{2}\sqrt{\left( p^{2}-4q\right) }. \end{equation} Now, it is straightforward to see that $e^{r_{1}t}$ and $e^{r_{2}t}$ both are solutions to the homogeneous equation, by noting that for $i\in \left\{ 1,2\right\}$% \begin{eqnarray} \frac{de^{r_{i}t}}{dt} &=&r_{i}e^{r_{i}t} \\ \frac{d^{2}e^{r_{i}t}}{dt^{2}} &=&r_{i}^{2}e^{r_{i}t}\rightarrow \\ \frac{d^{2}e^{r_{i}t}}{dt^{2}}+p\frac{de^{r_{i}t}}{dt}+qe^{r_{i}t} &=&\left( r_{i}^{2}+pr_{i}+q\right) e^{r_{i}t}=0. \end{eqnarray} So then, result \ref{solutiondiff1} tells us that the general solution to the homogenous equation is \begin{equation} y_{h}\left( t\right) =c_{1}e^{r_{1}t}+c_{2}e^{r_{2}t}, \end{equation}% \emph{provided }$e^{r_{1}t}$ and $e^{r_{2}t}$ are linearly independent. It is easy to verify that this is the case, if and only if $r_{1}\neq r_{2}.$ \subsubsection{Complex roots} Complex roots pose no particular difficulty, we simply have to recall that for any real number $b$% \begin{eqnarray} e^{bi} &=&\cos \left( b\right) +i\sin \left( b\right) , \\ e^{-bi} &=&\cos \left( b\right) -i\sin \left( b\right) , \end{eqnarray} yielding for the complex roots $r_{1,2}=a\pm bi,$% \begin{gather} c_{1}e^{\left( a+bi\right) t}+c_{2}e^{\left( a-bi\right) t}= \\ c_{1}e^{at}\left( \cos \left( bt\right) +i\sin \left( bt\right) \right) +c_{2}e^{at}\left( \cos \left( bt\right) -i\sin \left( bt\right) \right) \\ =e^{at}\left( \left( c_{1}+c_{2}\right) \cos \left( bt\right) +\left( c_{1}-c_{2}\right) i\sin \left( bt\right) \right) \end{gather} Note here, that $c_{1}$ and $c_{2}$ are any constants in the space of complex numbers are. Defining \begin{eqnarray} c_{1}+c_{2} &\equiv &\bar{c}_{1}, \\ \left( c_{1}-c_{2}\right) i &\equiv &\bar{c}_{2}, \end{eqnarray} It turns out that for any $\bar{c}_{1}$ and $\bar{c}_{2}$ on the real line, we can find $c_{1}$ and $c_{2}$ satisfying this definition. Since we, at least in economics, are only interested in solutions in the real space, we can use the restricted set of constants satisfying $\bar{c}_{1}$ and $\bar{c}% _{2}$ being on the real line. We can then write the general solution \emph{% in the real space} as \begin{equation} y_{h}\left( t\right) =e^{at}\left( \bar{c}_{1}\cos \left( bt\right) +\bar{c}% _{2}\sin \left( bt\right) \right) . \end{equation} \subsubsection{Repeated roots} The general solution to the homogenous equation in the case when the roots are repeated, i.e., $r_{1}=r_{2}\equiv r$ is \begin{equation} y_{h}\left( t\right) =c_{1}e^{rt}+c_{2}te^{rt}. \end{equation} Convince yourself that they are linearly independent and check that they are both solutions if the roots are repeated but not otherwise! \subsubsection{Non-Homogeneous Equations with Constant Coefficients} Relying on result \ref{solutiondiff1}, the only added problem when we have a non-homogeneous equation is that we have to find \emph{one} particular solution to the complete equation. Consider \begin{equation} \ddot{y}\left( t\right) +p\dot{y}\left( t\right) +qy\left( t\right) =R\left( t\right) . \end{equation} Typically we guess a form of this solution and then use the \emph{method of undetermined coefficients}. Often a good guess is a solution of a form similar to $R(t)$, e.g., if it is a polynomial of degree $n$ we guess on a general $n^{\prime }th$ degree polynomial with unknown coefficients. The simplest example is if $R\left( t\right)$ equals a constant $R$, we then guess on a constant, $y_{p}\left( t\right) =y^{ss}$, i.e., a \emph{steady state}, in which case $\ddot{y}\left( t\right)$ and $\dot{y}\left( t\right)$ both are zero. To satisfy the differential equation, \begin{eqnarray} qy^{ss} &=&R \\ y^{ss} &=&\frac{R}{q}. \end{eqnarray} Another example is \begin{equation} \ddot{y}\left( t\right) -2\dot{y}\left( t\right) +y\left( t\right) =3t^{2}+t, \end{equation} in which case we guess \begin{equation} y_{p}\left( t\right) =At^{2}+Bt+C, \end{equation} for some, yet \emph{undetermined coefficients }$A,B$ and $C$\emph{. }We then solve for these constants by substituting into the differential equation \begin{equation} 2A-2\left( 2At+B\right) +At^{2}+Bt+C=3t^{2}+t. \end{equation} For this to hold for each $t$, we need \begin{eqnarray} A &=&3 \\ -4A+B &=&1 \\ 2A-2B+C &=&0 \end{eqnarray} yielding $A=3,B=13,$ and $C=20.$ So a particular solution is \begin{equation} y_{p}\left( t\right) =3t^{2}+13t+20. \end{equation} The characteristic equation is \begin{equation} r^{2}-2r+1=\left( r-1\right) ^{2}\rightarrow r_{1,2}=1. \end{equation} So the general solution is \begin{equation} y\left( t\right) =c_{1}e^{t}+c_{2}te^{t}+3t^{2}+13t+20. \end{equation} \subsubsection{Linear nth Order Differential Equations} Consider the nth order differential equation \begin{equation} y^{n}\left( t\right) +P_{1}y^{n-1}\left( t\right) +...+P_{n}y\left( t\right) =R\left( t\right) , \label{norderdiff} \end{equation} where \begin{equation} y^{n}\left( t\right) \equiv \frac{d^{n}y\left( t\right) }{dt^{n}}. \end{equation} The solution technique is here exactly analogous to the second order case. First, some definitions \begin{definition} For any differentiable function $y\left( t\right)$, the differential operator $D$ is defined as \begin{equation} Dy\left( t\right) \equiv \frac{dy\left( t\right) }{dt}, \end{equation} satisfying $D^{n}y\left( t\right) \equiv \frac{d^{n}y\left( t\right) }{dt^{n}% }.$ \end{definition} Using this definition, and defining the following polynomial, \begin{equation} P\left( r\right) \equiv r^{n}+P_{1}r^{n-1}+...+P_{n}, \end{equation}% we have \begin{eqnarray} P\left( D\right) y\left( t\right) &=&\left( D^{n}+P_{1}D^{n-1}+...+P_{n}\right) y\left( t\right) \\ &=&D^{n}y\left( t\right) +P_{1}D^{n-1}y\left( t\right) +...+P_{n}y\left( t\right) \end{eqnarray}% so we can write (\ref{norderdiff}) in the condensed form, \begin{equation} P\left( D\right) y\left( t\right) =R\left( t\right) \end{equation}% and the characteristic equation is \begin{equation} P\left( r\right) =0, \end{equation}% with roots $r_{1,...,n}$ The general solution to the homogenous part is then the sum of the $n$ solutions corresponding to each of $n$ roots. The only thing to note is that if one root, $r$, is repeated $k\geq 2$ times, the solutions corresponding to this root are given by \begin{equation} c_{1}e^{rt}+c_{2}te^{rt}+...+c_{k}t^{k-1}e^{rt}. \end{equation} Repeated complex roots are handled the same way. Say the root $a\pm bi$ is repeated in $k$ pairs. Their contribution to the general solution is given by \begin{eqnarray} &&e^{at}\left( c_{1}\cos \left( bt\right) +c_{2}\sin \left( bt\right) +t\left( c_{3}\cos \left( bt\right) +c_{4}\sin \left( bt\right) \right) \right) +... \\ &&+e^{at}t^{k-1}\left( c_{2k-1}\cos \left( bt\right) +c_{2k}\sin \left( bt\right) \right) . \end{eqnarray} A particular solution to the complete equations can often be solved by the method of undetermined coefficients. \subsection{Stability} From the solutions to the differential equations we have seen we find that the terms corresponding to roots that have positive real parts (unstable, or explosive roots) tend to explode as $t$ goes to infinity. This means that also the solution explodes unless the corresponding integration constants are zero. Terms with roots that have negative real parts (stable roots), on the other hand, always converge to zero. If all roots are \emph{strictly} negative for a differential equation, $% P\left( D\right) y\left( t\right) =r$, the system converges to a unique point as time goes to infinity, wherever it starts. This point is often called a \emph{sink} and the system is defined as globally stable. If all roots are strictly positive the system is globally unstable, but there is still a steady-state, this is sometimes called the \emph{origin}. The reason for this is that if the system is unstable in all dimensions when time goes forward, it will be \emph{stable} in all dimensions if time goes in reverse. Starting from \emph{any} point and going \emph{backward, }one reaches, in the limit the steady state, i.e., it is the origin of all paths. If a system has both stable and unstable roots, it is called saddle-path stable. Then, in some sub-dimensions it is stable. \subsubsection{Non-linear equations and local stability} Look at the nonlinear differential equation \begin{equation} \dot{x}\left( t\right) =x\left( t\right) ^{2}-4. \end{equation} Although we have not learned how to solve such an equation, we can say something about it. Let us plot the relation $x\left( t\right) \rightarrow \dot{x}\left( t\right) ;$ \FRAME{itbpF}{3.6478in}{3.0026in}{0in}{}{}{slide5.emf}{\special{language "Scientific Word";type "GRAPHIC";maintain-aspect-ratio TRUE;display "USEDEF";valid_file "F";width 3.6478in;height 3.0026in;depth 0in;original-width 10.0016in;original-height 7.5135in;cropleft "0";croptop "1";cropright "0.5458";cropbottom "0.4029";filename 'Figures/Slide5.EMF';file-properties "XNPEU";}} We see that $x(t)=2$ and $x(t)=-2$ are stationary points. We also see that $% x(t)=-2$ is locally stable. In the region $[-\infty ,2)$ $x$ converge to $% x=-2$. 2 is an unstable stationary point and in the region $(2,\infty ]$, $x$ explodes over time. In a plot $x\left( t\right) \rightarrow \dot{x}\left( t\right) ,$ local stability is equivalent to a negative slope at the stationary point. \subsection{Systems of Linear First Order Differential Equations} Consider the following system of two first order differential equations \begin{eqnarray} \left[ \begin{array}{c} \dot{y}_{1}\left( t\right) \\ \dot{y}_{2}\left( t\right)% \end{array}% \right] &=&\left[ \begin{array}{cc} a_{11} & a_{12} \\ a_{21} & a_{22}% \end{array}% \right] \left[ \begin{array}{c} y_{1}\left( t\right) \\ y_{2}\left( t\right)% \end{array}% \right] +\left[ \begin{array}{c} p_{1}\left( t\right) \\ p_{2}\left( t\right)% \end{array}% \right] \label{diffsystem1} \\ \mathbf{\dot{y}}\left( t\right) &=&\mathbf{Ay}\left( t\right) +\mathbf{p}% \left( t\right) \end{eqnarray} As in the one equations case we start by finding the general solutions to the homogeneous part. This plus some particular solution is the general solution to the complete system. If the off diagonal terms are zero the solution to the homogeneous part is trivial, since there is no interdependency between the equations \begin{equation} \left[ \begin{array}{c} \dot{y}_{1}\left( t\right) \\ \dot{y}_{2}\left( t\right)% \end{array} \right] =\left[ \begin{array}{cc} a_{11} & 0 \\ 0 & a_{22}% \end{array} \right] \left[ \begin{array}{c} y_{1}\left( t\right) \\ y_{2}\left( t\right)% \end{array} \right] \label{autonomuos} \end{equation} \begin{equation} \rightarrow \left[ \begin{array}{c} y_{1}\left( t\right) \\ y_{2}\left( t\right)% \end{array} \right] =\left[ \begin{array}{cc} e^{a_{11}t} & 0 \\ 0 & e^{a_{22}t}% \end{array} \right] \left[ \begin{array}{c} c_{1} \\ c_{2}% \end{array} \right] . \end{equation} The system in (\ref{autonomuos}) has an important property, time has no direct effect on the law-of-motion. Given knowledge of $y_{1}\left( t\right)$ and $y_{2}\left( t\right)$, $\dot{y}_{1}\left( t\right)$ and $\dot{y}% _{2}\left( t\right)$ are fully determined, \emph{regardless of t. }A system that has this property is called \emph{autonomous.} The system (\ref% {diffsystem1}) is not an autonomous system, unless $\mathbf{p}\left( t\right)$ is constant. Note that the homogeneous solution is always autonomous, given, of course, that the parameters are constant. The behavior of an autonomous system can be depicted in a graph, a phase diagram. We see in the phase diagram that if the roots are stable, i.e., negative, the homogeneous part always goes to zero as $t$ goes to infinity. With only one root stable, there is just one stable path. \FRAME{itbpF}{6.1575in}{2.5002in}{0in}{}{}{slide6.emf}{\special{language "Scientific Word";type "GRAPHIC";maintain-aspect-ratio TRUE;display "USEDEF";valid_file "F";width 6.1575in;height 2.5002in;depth 0in;original-width 10.0016in;original-height 7.5135in;cropleft "0.0377";croptop "1";cropright "0.9620";cropbottom "0.5037";filename 'Figures/Slide6.EMF';file-properties "XNPEU";}} The fact that it is trivial to solve a diagonal system, suggests a way of finding the solution to the homogeneous part of (\ref{diffsystem1}). Suppose we can make a transformation of the variables so that the transformed system is diagonal. Start by defining the new set of variables \begin{eqnarray} \left[ \begin{array}{c} x_{1}\left( t\right) \\ x_{2}\left( t\right)% \end{array} \right] &\equiv &\mathbf{B}\left[ \begin{array}{c} y_{1}\left( t\right) \\ y_{2}\left( t\right)% \end{array} \right] \rightarrow \\ \left[ \begin{array}{c} \dot{x}_{1}\left( t\right) \\ \dot{x}_{2}\left( t\right)% \end{array} \right] &\equiv &\mathbf{B}\left[ \begin{array}{c} \dot{y}_{1}\left( t\right) \\ \dot{y}_{2}\left( t\right)% \end{array} \right] =\mathbf{BA}\left[ \begin{array}{c} y_{1}\left( t\right) \\ y_{2}\left( t\right)% \end{array} \right] =\mathbf{BAB}^{-1}\left[ \begin{array}{c} x_{1}\left( t\right) \\ x_{2}\left( t\right)% \end{array} \right] \end{eqnarray} If we can find a $\mathbf{B}$ such that $\mathbf{BAB}^{-1}$ is diagonal we are half way. The solutions for $\mathbf{x}$ is then \begin{eqnarray} && \\ \left[ \begin{array}{c} x_{1}\left( t\right) \\ x_{2}\left( t\right)% \end{array}% \right] &=&\left[ \begin{array}{cc} e^{r_{1}t} & 0 \\ 0 & e^{r_{2}t}% \end{array}% \right] \left[ \begin{array}{c} c_{1} \\ c_{2}% \end{array}% \right] \end{eqnarray} where $r_{1,2}$ are the diagonal terms of the matrix $\mathbf{BAB}^{-1}$.The solution for $\mathbf{y}$ then follows from the definition of $\mathbf{x}$% \begin{equation} \left[ \begin{array}{c} y_{1}\left( t\right) \\ y_{2}\left( t\right)% \end{array} \right] =\mathbf{B}^{-1}\left[ \begin{array}{cc} e^{r_{1}t} & 0 \\ 0 & e^{r_{2}t}% \end{array} \right] \left[ \begin{array}{c} c_{1} \\ c_{2}% \end{array} \right] . \end{equation} From linear algebra we know that $\mathbf{B}^{-1}$ is the matrix of eigenvectors of $\mathbf{A}$ and that the diagonal terms of $\mathbf{BAB}% ^{-1}$ are the corresponding eigenvalues. The eigenvalues are given by the characteristic equation of $\mathbf{A}$% \begin{eqnarray} \left\vert \left[ \begin{array}{cc} a_{11}-r & a_{12} \\ a_{21} & a_{22}-r% \end{array}% \right] \right\vert &=&0 \\ &\rightarrow &\left( a_{11}-r\right) \left( a_{22}-r\right) -a_{12}a_{21}=0 \\ r^{2}-rTr\left( \mathbf{A}\right) +\left\vert \mathbf{A}\right\vert &=&0, \end{eqnarray}% where $Tr\left( \mathbf{A}\right)$ is the trace of $\mathbf{A.}$ The only crux is that we need the roots to be distinct, otherwise $\mathbf{B}^{-1}$ is not always invertible. Distinct roots imply that $\mathbf{B}^{-1}$ is invertible. (If $\mathbf{A}$ is symmetric $\mathbf{B}^{-1}$ is also invertible.) Let us draw a phase diagram with the eigenvectors of $\mathbf{A}$. The dynamic system behaves as the diagonal one but the eigenvectors have replaced the standard, orthogonal axes. \FRAME{itbpF}{3.8043in}{2.5002in}{0in}{}{}{slide7.emf}{\special{language "Scientific Word";type "GRAPHIC";maintain-aspect-ratio TRUE;display "USEDEF";valid_file "F";width 3.8043in;height 2.5002in;depth 0in;original-width 10.2385in;original-height 7.5135in;cropleft "0";croptop "1";cropright "0.5563";cropbottom "0.5037";filename 'Figures/Slide7.EMF';file-properties "XNPEU";}} What is remaining is to find a particular solution of the complete system. One way is here to use the method of undetermined coefficients. Sometimes, we can find a steady state of the system, i.e., a point where the time derivatives are all zero. This is easy if the exogenous part is constant. We then set the differential equal to zero so \begin{eqnarray} 0 &=&\mathbf{Ay}\left( t\right) +\mathbf{p} \\ &\rightarrow &\mathbf{y}_{p}\left( t\right) =\mathbf{y}^{ss}=-\mathbf{A}^{-1}% \mathbf{p.} \end{eqnarray} Given that we know the value of $\mathbf{y}\left( 0\right)$ we can now give the general solution. The formula is given in matrix form and is valid for any dimension of the system. First define \begin{equation} \mathbf{r}\left( t\right) \equiv \left[ \begin{array}{cccc} e^{r_{1}t} & 0 & .. & 0 \\ 0 & e^{r_{2}t} & .. & 0 \\ . & . & . & . \\ 0 & 0 & . & e^{r_{n}t}% \end{array} \right] , \end{equation} then we have \begin{eqnarray} \mathbf{y}\left( t\right) &=&\mathbf{B}^{-1}\mathbf{x}\left( t\right) +% \mathbf{y}^{ss}=\mathbf{B}^{-1}\mathbf{r}\left( t\right) \mathbf{c+y}^{ss} \\ \mathbf{y}\left( 0\right) &=&\mathbf{B}^{-1}\mathbf{x}\left( 0\right) +% \mathbf{y}^{ss}=\mathbf{B}^{-1}\mathbf{c+y}^{ss} \\ \mathbf{c} &=&\mathbf{B}\left( \mathbf{\mathbf{y}\left( 0\right) -y}% ^{ss}\right) \\ &\rightarrow &\mathbf{y}\left( t\right) =\mathbf{B}^{-1}\mathbf{r}\left( t\right) \mathbf{B}\left( \mathbf{\mathbf{y}\left( 0\right) -y}^{ss}\right) \mathbf{+y}^{ss}. \end{eqnarray} The method outlined above works also in the case of complex roots of the characteristic equation. If the roots are $a\pm bi$ we have \begin{equation} \mathbf{y}\left( t\right) =\mathbf{B}^{-1}\left[ \begin{array}{cc} e^{\left( a+bi\right) t} & 0 \\ 0 & e^{\left( a-bi\right) t}% \end{array} \right] \mathbf{c+y}^{ss} \end{equation} Example; \begin{eqnarray} \left[ \begin{array}{c} \dot{y}_{1}\left( t\right) \\ \dot{y}_{2}\left( t\right)% \end{array} \right] &=&\left[ \begin{array}{cc} -1 & -1 \\ 1 & -1% \end{array} \right] \left[ \begin{array}{c} y_{1}\left( t\right) \\ y_{2}\left( t\right)% \end{array} \right] +\left[ \begin{array}{c} 1 \\ 1% \end{array} \right] \\ r_{1,2} &=&-1\pm i \\ \mathbf{B}^{-1} &=&\left[ \begin{array}{cc} i & -i \\ 1 & 1% \end{array} \right] . \end{eqnarray} So, \begin{gather} \rightarrow \mathbf{y}\left( t\right) = \\ \left[ \begin{array}{cc} i & -i \\ 1 & 1% \end{array} \right] \left[ \begin{array}{cc} e^{-t}\left( \cos t+i\sin t\right) & 0 \\ 0 & e^{-t}\left( \cos t-i\sin t\right)% \end{array} \right] \left[ \begin{array}{c} c_{1} \\ c_{2}% \end{array} \right] \mathbf{+y}^{ss} \\ =e^{-t}\left[ \begin{array}{c} i\left( c_{1}-c_{2}\right) \cos t-\left( c_{1}+c_{2}\right) \sin t \\ \left( c_{1}+c_{2}\right) \cos t+i\left( c_{1}-c_{2}\right) \sin t% \end{array} \right] \mathbf{+y}^{ss} \\ =e^{-t}\left[ \begin{array}{c} \tilde{c}_{1}\cos t+\tilde{c}_{2}\sin t \\ -\tilde{c}_{2}\cos t+\tilde{c}_{1}\sin t% \end{array} \right] \mathbf{+y}^{ss} \end{gather} The steady-state is found from, \begin{eqnarray} \mathbf{0} &=&\left[ \begin{array}{cc} -1 & -1 \\ 1 & -1% \end{array} \right] \left[ \begin{array}{c} y_{1}\left( t\right) \\ y_{2}\left( t\right)% \end{array} \right] +\left[ \begin{array}{c} 1 \\ 1% \end{array} \right] \\ \left[ \begin{array}{c} y_{1}^{ss} \\ y_{2}^{ss}% \end{array} \right] &=&-\left[ \begin{array}{cc} -1 & -1 \\ 1 & -1% \end{array} \right] ^{-1}\left[ \begin{array}{c} 1 \\ 1% \end{array} \right] =\left[ \begin{array}{c} 0 \\ 1% \end{array} \right] . \end{eqnarray} If we know $\mathbf{y}\left( 0\right)$, and using $\cos 0=1$ and $\sin 0=0,$% we get \begin{gather} \left[ \begin{array}{c} y_{1}\left( 0\right) \\ y_{2}\left( 0\right)% \end{array} \right] =\left[ \begin{array}{c} \tilde{c}_{1} \\ -\tilde{c}_{2}% \end{array} \right] \mathbf{+}\left[ \begin{array}{c} 0 \\ 1% \end{array} \right] \\ \rightarrow \mathbf{y}\left( t\right) =e^{-t}\left[ \begin{array}{c} y_{1}\left( 0\right) \cos t-\left( y_{2}\left( 0\right) -1\right) \sin t \\ \left( y_{2}\left( 0\right) -1\right) \cos t+y_{1}\left( 0\right) \sin t% \end{array} \right] \mathbf{+}\left[ \begin{array}{c} 0 \\ 1% \end{array} \right] \end{gather} A phase-diagram of this system is an inward spiral. \FRAME{itbpF}{3.2534in}{2.4984in}{0in}{}{}{slide8.emf}{\special{language "Scientific Word";type "GRAPHIC";maintain-aspect-ratio TRUE;display "USEDEF";valid_file "F";width 3.2534in;height 2.4984in;depth 0in;original-width 10.1866in;original-height 7.5299in;cropleft "0";croptop "1";cropright "0.4786";cropbottom "0.5038";filename 'Figures/Slide8.EMF';file-properties "XNPEU";}} How would it look like if the real part of the root was $-1$ or $+1$? \subsubsection{Equivalent Systems} In the repeated root case the matrix of the eigenvectors may be singular, so that we cannot find $\mathbf{B}^{-1}$. Then we use the method of equivalent systems, described in this section. A linear \emph{nth} order differential equation is equivalent to a system of $n$ first order differential equations. Consider,\qquad \qquad \begin{equation} \dddot{y}\left( t\right) +a_{1}\ddot{y}\left( t\right) +a_{2}\dot{y}\left( t\right) +a_{3}y\left( t\right) =R\left( t\right) \end{equation} We can transform this into a system of first order differential equation by defining \begin{eqnarray} \dot{y}\left( t\right) &\equiv &x_{1}\left( t\right) , \\ \ddot{y}\left( t\right) &\equiv &x_{2}\left( t\right) =\dot{x}_{1}\left( t\right) , \\ \dddot{y}\left( t\right) &=&\dot{x}_{2}\left( t\right) , \end{eqnarray} giving \begin{equation} \left[ \begin{array}{c} \dot{y}\left( t\right) \\ \dot{x}_{1}\left( t\right) \\ \dot{x}_{2}\left( t\right)% \end{array} \right] =\left[ \begin{array}{ccc} 0 & 1 & 0 \\ 0 & 0 & 1 \\ -a_{3} & -a_{2} & -a_{1}% \end{array} \right] \left[ \begin{array}{c} y\left( t\right) \\ x_{1}\left( t\right) \\ x_{2}\left( t\right)% \end{array} \right] +\left[ \begin{array}{c} 0 \\ 0 \\ R\left( t\right)% \end{array} \right] \end{equation} Since the equations are equivalent they consequently have the same solutions. Sometimes one of the transformations is more convenient to solve. Let us also transform a two dimensional system into a second order equation \begin{equation} \left[ \begin{array}{c} \dot{x}_{1}\left( t\right) \\ \dot{x}_{2}\left( t\right)% \end{array} \right] =\left[ \begin{array}{cc} a_{11} & a_{12} \\ a_{21} & a_{22}% \end{array} \right] \left[ \begin{array}{c} x_{1}\left( t\right) \\ x_{2}\left( t\right)% \end{array} \right] +\left[ \begin{array}{c} k_{1} \\ k_{2}% \end{array} \right] . \label{simplesys} \end{equation} First use the first equation to express $x_{2}$ and then take the time derivative of the first. Then we can eliminate $x_{2}$ and its time derivative \begin{gather} \dot{x}_{1}\left( t\right) =a_{11}x_{1}\left( t\right) +a_{12}x_{2}\left( t\right) +k_{1} \\ x_{2}\left( t\right) =\frac{\dot{x}_{1}\left( t\right) -a_{11}x_{1}\left( t\right) -k_{1}}{a_{12}} \\ \dot{x}_{2}\left( t\right) =\frac{\ddot{x}_{1}\left( t\right) -a_{11}\dot{x}% _{1}\left( t\right) }{a_{12}} \\ \frac{\ddot{x}_{1}\left( t\right) -a_{11}\dot{x}_{1}\left( t\right) }{a_{12}}% =a_{21}x_{1}\left( t\right) +a_{22}\frac{\dot{x}_{1}\left( t\right) -a_{11}x_{1}\left( t\right) -k_{1}}{a_{12}}+k_{2} \\ \ddot{x}_{1}\left( t\right) -\left( a_{11}+a_{22}\right) \dot{x}_{1}\left( t\right) +\left( a_{22}a_{11}-a_{12}a_{21}\right) x_{1}\left( t\right) =-a_{22}k_{1}+a_{12}k_{2} \end{gather} Note that the characteristic equation of this second order equation is the same as the one for the system in (\ref{simplesys}). Consequently, the roots and thus the dynamics are identical. \subsection{Non-linear systems and Linearization\label{phasediagram}} Phase diagrams are convenient to analyze the behavior of a 2 dimensional system qualitatively that we cannot or prefer not to solve explicitly. E.g., if\qquad \begin{eqnarray} \dot{c}\left( t\right) &=&g_{1}\left( c\left( t\right) ,k\left( t\right) \right) \\ \dot{k}\left( t\right) &=&g_{2}\left( c\left( t\right) ,k\left( t\right) \right) \end{eqnarray} The first step here is to find the two curves in the $c,k$-space where $c$ and $k$, respectively are constant. Setting the time derivatives equal to zero defines two relations between $c$ and $k$, which we denote by $G_{1}$ and $G_{2}.$% \begin{eqnarray} g_{1}\left( c\left( t\right) ,k\left( t\right) \right) &=&0\rightarrow c=G_{1}\left( k\right) \\ g_{2}\left( c\left( t\right) ,k\left( t\right) \right) &=&0\rightarrow c=G_{2}\left( k\right) \end{eqnarray}% We then draw these curves in the $c,k$-space. For example, you will in the macro course analyze the Ramsey optimal consumption problem, where output is $f\left( k\right) ,$ interest rate is $f^{\prime }\left( k\right) ,$ the subjective discount rate is $\theta$ and utility is $u\left( c\right)$ where $u$ and $f$ are assumed to be concave functions. The model will deliver the following system of differential equations \begin{eqnarray} \dot{c}\left( t\right) &=&-\frac{u^{\prime }\left( c\left( t\right) \right) }{u^{\prime \prime }\left( c\left( t\right) \right) }\left( f^{\prime }\left( k\left( t\right) \right) -\theta \right) \label{Ramsey} \\ \dot{k}\left( t\right) &=&f\left( k\left( t\right) \right) -c\left( t\right) . \end{eqnarray}% Setting the time derivatives to zero we get \ \qquad \begin{eqnarray} f^{\prime }\left( k\left( t\right) \right) &=&\theta , \\ c &=&f\left( k\left( t\right) \right) . \end{eqnarray} Draw these curves in the $c,k$ space \FRAME{ftbpF}{4.651in}{3.2041in}{0pt}{}{}{slide10.emf}{\special{language "Scientific Word";type "GRAPHIC";maintain-aspect-ratio TRUE;display "USEDEF";valid_file "F";width 4.651in;height 3.2041in;depth 0pt;original-width 10.0016in;original-height 7.5135in;cropleft "0";croptop "1";cropright "0.6971";cropbottom "0.3626";filename 'Figures/Slide10.EMF';file-properties "XNPEU";}} We then have to find the signs of time derivatives, and above and below their respective zero motion curves. From (\ref{Ramsey}), we see that \begin{eqnarray} \frac{\partial \dot{c}\left( t\right) }{\partial k} &=&-\frac{u^{\prime }\left( c\right) }{u^{\prime \prime }\left( c\right) }f^{\prime \prime }\left( k\right) <0 \\ \frac{\partial \dot{k}\left( t\right) }{\partial c} &=&-1. \end{eqnarray} This means that $\dot{c}$ is positive to the left of and negative to the right of the curve $\dot{c}=0$. For $\dot{k}$, we find that it is positive below and negative above $\dot{k}=0$. Then draw these motions as arrows in the phase diagram. Note that no paths ever can cross. \FRAME{itbpF}{4.651in}{3.0026in}{0in}{}{}{slide11.emf}{\special{language "Scientific Word";type "GRAPHIC";maintain-aspect-ratio TRUE;display "USEDEF";valid_file "F";width 4.651in;height 3.0026in;depth 0in;original-width 10.0016in;original-height 7.5135in;cropleft "0";croptop "1";cropright "0.6971";cropbottom "0.4029";filename 'Figures/Slide11.EMF';file-properties "XNPEU";}} We conclude that this system has saddle point characteristics and thus has only one stable trajectory towards the steady state. The behavior close to the steady state should also be evaluated by means of linearization around the steady state. We do that by approximating in the following way \begin{equation} \left[ \begin{array}{c} \dot{c}\left( t\right) \\ \dot{k}\left( t\right)% \end{array}% \right] \approx \left[ \begin{array}{cc} \frac{\partial g_{1}\left( c^{ss},k^{ss}\right) }{\partial c} & \frac{% \partial g_{1}\left( c^{ss},k^{ss}\right) }{\partial k} \\ \frac{\partial g_{2}\left( c^{ss},k^{ss}\right) }{\partial c} & \frac{% \partial g_{2}\left( c^{ss},k^{ss}\right) }{\partial c}% \end{array}% \right] \left[ \begin{array}{c} c\left( t\right) -c^{ss} \\ k\left( t\right) -k^{ss}% \end{array}% \right] . \end{equation}% with an obvious generalization to higher dimensions. We now evaluate the roots of the matrix of derivatives. In the example we find that the coefficient matrix is \begin{eqnarray} &&\left[ \begin{array}{cc} \frac{\partial g_{1}\left( c^{ss},k^{ss}\right) }{\partial c} & \frac{% \partial g_{1}\left( c^{ss},k^{ss}\right) }{\partial k} \\ \frac{\partial g_{2}\left( c^{ss},k^{ss}\right) }{\partial c} & \frac{% \partial g_{2}\left( c^{ss},k^{ss}\right) }{\partial c}% \end{array} \right] \\ &=&\left[ \begin{array}{cc} -\left( f^{\prime }-\theta \right) \frac{\partial }{\partial c}\left( \frac{% u^{\prime }}{u^{\prime \prime }}\right) & -\frac{u^{\prime }}{u^{\prime \prime }}f^{\prime \prime } \\ -1 & f^{\prime }% \end{array} \right] , \end{eqnarray} where all functions are evaluated at the steady state. There, $f^{\prime }=\theta$, implying that the matrix simplifies to \begin{equation} \left[ \begin{array}{cc} 0 & -\frac{u^{\prime }}{u^{\prime \prime }}f^{\prime \prime } \\ -1 & f^{\prime }% \end{array} \right] \end{equation} with eigenvalues \begin{equation} \frac{1}{2}\left( f^{\prime }\pm \sqrt{(f^{\prime })^{2}+4\frac{u^{\prime }}{% u^{\prime \prime }}f^{\prime \prime }}\right) \end{equation} which clearly are of opposite signs. \subsection{Example: Steady-state asset distributions} In this example, we use our derived skills to find the steady-state wealth distribution in a simple model. As we will see, here the model gives us a system of differential equations for the wealth distribution that we can solve easily. We will see that the methods we have learned work as well when the differential equations are in wealth, rather than time. In (Hassler and Mora, JPub, 1999), we analyze preferences over unemployment insurance in a very simple continuous time economy where agents can borrow and save at rate $r$. In the model, employed individuals earn a wage $w$ per unit of time and become unemployed with an instantaneous probability $q.$ This means that over a small (infinitesimal) interval of time $dt,$ the probability of becoming unemployed is $qdt.$ Unemployed individuals receive a flow of unemployment benefits $b$ and become rehired with instantaneous probability $hdt.$ In addition, we assume that there is a instantaneous death probability of $\delta$ and that there is an inflow of newborn unemployed with zero assets of $\delta$ so that the total population is constant. In the paper, we show that if individuals have CARA utility ($U=-e^{\gamma c}$) and wages and benefits are constant, individuals choose constant savings amounts for each of the two employment states. Denoting these $s_{e}$ and $s_{u}$, where $s_{e}>0$ and $s_{u}<0$, we have that individual asset accumulation for the two types, conditional on surviving is \begin{equation} \dot{A}_{t}=s_{e}, \end{equation}% for employed and \begin{equation} \dot{A}_{t}=s_{u} \end{equation}% for unemployed. In the paper, we don't calculate the steady state wealth distribution of assets. That's the purpose of this exercise. First, we calculate the steady state share of unemployed, $\mu _{u}$. For this purpose, we note that in steady state, the inflow and the outflow to the stock of unemployed must be constant. Thus, \begin{eqnarray} \left( 1-\mu _{u}\right) q+\delta &=&\mu _{u}\left( h+\delta \right) \\ &\rightarrow &\mu _{u}=\frac{q+\delta }{q+h+\delta } \end{eqnarray} Now, we want to calculate the steady state distribution of assets among employed and unemployed in this economy. Let us denote these densities, by $% f_{e}(A)$ and $f_{u}(A).$ We will derive these by solving a system of differential equations. Consider first the number (density) of employed individuals with assets $A\neq 0$, given by $f_{e}\left( A\right)$. Over a small period $dt$, a number $f_{e}(A)\left( 1-\left( \delta +q\right) dt\right)$ of them remain employed and alive. Over the same time, a number $% f_{u}(A)hdt$ of unemployed with assets $A$ become hired. Finally, over the period $dt$ these individuals add to their assets an amount $s_{e}dt.$ Writing this down yields, \begin{eqnarray} f_{e}(A+s_{e}dt) &=&f_{e}(A)\left( 1-\left( \delta +q\right) dt\right) +f_{u}(A)hdt, \\ f_{e}(A)+f_{e}^{\prime }(A)s_{e}dt &=&f_{e}(A)\left( 1-\left( \delta +q\right) dt\right) +f_{u}(A)hdt \\ f_{e}^{\prime }(A) &=&-f_{e}(A)\frac{\delta +q}{s_{e}}+f_{u}(A)\frac{h}{s_{e}% }. \end{eqnarray}% By the same reasoning, \begin{eqnarray} f_{u}(A+s_{e}dt) &=&f_{u}(A)\left( 1-\left( \delta +h\right) dt\right) +f_{e}(A)qdt, \\ f_{u}(A)+f_{u}^{\prime }(A)s_{e}dt &=&f_{u}(A)\left( 1-\left( \delta +h\right) dt\right) +f_{e}(A)qdt \\ f_{u}^{\prime }(A) &=&-f_{u}(A)\frac{\delta +h}{s_{u}}+f_{e}(A)\frac{q}{s_{u}% }, \end{eqnarray}% yielding the system \begin{eqnarray} \left[ \begin{array}{c} f_{e}^{\prime }(A) \\ f_{u}^{\prime }(A)% \end{array}% \right] &=&\left[ \begin{array}{cc} -\frac{\delta +q}{s_{e}} & \frac{h}{s_{e}} \\ \frac{q}{s_{u}} & -\frac{\delta +h}{s_{u}}% \end{array}% \right] \left[ \begin{array}{c} f_{e}(A) \\ f_{u}(A)% \end{array}% \right] . \\ &\equiv &\mathbf{A}\left[ \begin{array}{c} f_{e}(A) \\ f_{u}(A)% \end{array}% \right] \end{eqnarray} Let us now, assign some numbers to the parameters. Say $\delta =1/40,h=2,q=1/5,s_{e}=1$ and $s_{u}=-2$, all measured as probabilities per year. Then, we can calculate the roots and the eigenvalues \begin{eqnarray} r_{1} &=&\frac{63}{160}+\frac{1}{160}\sqrt{4681}\approx .821 \\ r_{2} &=&\frac{63}{160}-\frac{1}{160}\sqrt{4681}\approx -0.0339 \\ \mathbf{B}^{-1} &=&\left[ \begin{array}{cc} 1 & 1 \\ \frac{99}{320}+\frac{1}{320}\sqrt{4681} & \frac{99}{320}-\frac{1}{320}\sqrt{% 4681}% \end{array}% \right] \approx \left[ \begin{array}{cc} 1.0 & 1.0 \\ .523 & 0.0956% \end{array}% \right] . \end{eqnarray} Thus, except at $A=0$, where there is an inflow of newborn, that we have not considered, we have \begin{gather} \left[ \begin{array}{c} f_{e}(A) \\ f_{u}(A)% \end{array} \right] =\mathbf{B}^{-1}\mathbf{r}\left( t\right) \mathbf{c} \\ =\left[ \begin{array}{cc} 1.0 & 1.0 \\ .523 & 0.0956% \end{array} \right] \left[ \begin{array}{cc} e^{.821A} & 0 \\ 0 & e^{-0.0339A}% \end{array} \right] \left[ \begin{array}{c} c_{1} \\ c_{2}% \end{array} \right] \\ =\left[ \begin{array}{c} c_{1}e^{.821A}+c_{2}e^{-.0339A} \\ .523c_{1}e^{.821A}+.0956c_{2}e^{-.0339A}% \end{array} \right] \end{gather} Now, we first note that $f_{e}(A)$ and $f_{u}(A)$ cannot be explosive in any of the directions. That would violate that these functions are densities, i.e., the sum of their respective integrals over the real line must be unity. Thus, for $A<0$, $c_{2}$ must be zero and for $A>0,$ $c_{1}=0.$ Furthermore, we know \begin{eqnarray} \int_{-\infty }^{\infty }f_{e}(A)dA &=&1-\mu _{u}=\frac{h}{q+h+\delta }% \approx 0.899 \\ \int_{-\infty }^{\infty }f_{u}(A)dA &=&\mu _{u}=\frac{q+\delta }{q+h+\delta }% \approx 0.101. \end{eqnarray} Using this, we can calculate the integrations constants. \begin{eqnarray} \int_{-\infty }^{\infty }f_{e}(A)dA &=&\int_{-\infty }^{0}c_{1}e^{.821A}dA+\int_{0}^{\infty }c_{2}e^{-.0339A}dA \\ &=&\frac{c_{1}}{0.821}+\frac{c_{2}}{.0339}=0.899 \\ \int_{-\infty }^{\infty }f_{u}(A)dA &=&\int_{-\infty }^{0}.523c_{1}e^{.821A}dA+\int_{0}^{\infty }.0956e^{-.0\,\allowbreak 339A}dA \\ &=&\frac{.523c_{1}}{0.821}+\frac{0.0956c_{2}}{.0339}=0.101 \end{eqnarray} Yielding,% \begin{eqnarray} c_{1} &=&0.0289, \\ c_{2} &=&0.0293. \end{eqnarray} This concludes our calculations, \begin{eqnarray} f_{e}(A) &=&\left\{ \begin{array}{c} 0.0289e^{.\,821A}\text{ for }A<0 \\ 0.0293e^{-.0339A}\text{ for }A>0% \end{array}% \right. , \\ f_{u}(A) &=&\left\{ \begin{array}{c} 0.0151e^{.821A}\text{ for }A<0 \\ 0.00280e^{-.0\,339A}\text{ for }A>0% \end{array}% \right. . \end{eqnarray}% \FRAME{itbpFX}{2.9996in}{2.3896in}{0in}{}{}{Plot}{\special{language "Scientific Word";type "MAPLEPLOT";width 2.9996in;height 2.3896in;depth 0in;display "USEDEF";plot_snapshots TRUE;mustRecompute FALSE;lastEngine "Maple";xmin "-5.000000";xmax "0";xviewmin "-5.000000";xviewmax "0";yviewmin "0";yviewmax "0.04";viewset"XY";plottype 4;numpoints 49;plotstyle "patch";axesstyle "normal";xis \TEXUX{A};var1name \TEXUX{$A$};function \TEXUX{$0.0289e^{.\,\allowbreak 821A}$};linecolor "black";linestyle 1;pointstyle "point";linethickness 1;lineAttributes "Solid";var1range "-5.000000,0";num-x-gridlines 49;curveColor "[flat::RGB:0000000000]";curveStyle "Line";function \TEXUX{$0.0151e^{.\,\allowbreak 821A}$};linecolor "black";linestyle 3;pointstyle "point";linethickness 1;lineAttributes "Dots";var1range "-5,0";num-x-gridlines 49;curveColor "[flat::RGB:0000000000]";curveStyle "Line";rangeset"X";valid_file "T";tempfilename 'KBKCSE02.wmf';tempfile-properties "XPR";}}\FRAME{itbpFX}{2.9996in}{2.3896in% }{0in}{}{}{Plot}{\special{language "Scientific Word";type "MAPLEPLOT";width 2.9996in;height 2.3896in;depth 0in;display "USEDEF";plot_snapshots TRUE;mustRecompute FALSE;lastEngine "Maple";xmin "0.000000";xmax "100.000000";xviewmin "0.000000";xviewmax "100.000000";yviewmin "0";yviewmax "0.04";viewset"XY";plottype 4;numpoints 49;plotstyle "patch";axesstyle "normal";xis \TEXUX{A};var1name \TEXUX{$A$};function \TEXUX{$0.02\allowbreak 93e^{-.0\,\allowbreak 339A}$};linecolor "black";linestyle 1;pointstyle "point";linethickness 1;lineAttributes "Solid";var1range "0.000000,100.000000";num-x-gridlines 49;curveColor "[flat::RGB:0000000000]";curveStyle "Line";function \TEXUX{$0.002\allowbreak 80e^{-.0\,\allowbreak 339A}$};linecolor "cyan";linestyle 3;pointstyle "point";linethickness 1;lineAttributes "Dots";var1range "0,100";num-x-gridlines 49;curveColor "[flat::RGB:0x00008080]";curveStyle "Line";rangeset"X";valid_file "T";tempfilename 'KBKCSF03.wmf';tempfile-properties "XPR";}} Plotting the densities, using a solid (dotted) line to denote wealth of employed (unemployed), we see that unemployed are more concentrated to the left. We can also calculate average assets among the two types from \begin{gather} A_{e}=\int_{-\infty }^{\infty }Af_{e}(A)dA= \\ \int_{-\infty }^{0}A\cdot 0.0289e^{.821A}dA+\int_{0}^{\infty }A\cdot 0.0293e^{-.0\,\allowbreak 339A}dA\approx 25.45 \\ A_{u}=\int_{-\infty }^{\infty }Af_{u}(A)dA= \\ \int_{-\infty }^{0}A\cdot 0.0151e^{.821A}dA+\int_{0}^{\infty }A\cdot 0.00280e^{-.0\,339A}dA\approx 2.41 \end{gather} \newpage \section{Difference equations} \subsection{Sums, forward and backward solutions} \subsubsection{Sums vs. Integrals} %TCIMACRO{\TeXButton{TeX field}{\setcounter{equation}{0}}}% %BeginExpansion \setcounter{equation}{0}% %EndExpansion Difference equations can be solved in ways very similar to how we solve differential equations. First, we will look at the analogy to integrals. Consider the difference equation \begin{equation} A_{t}-A_{t-1}\equiv \Delta A_{t}=q_{t} \label{differenceeq1} \end{equation}% where $\Delta A_{t}$ is the change in $A$ per unit (interval) of time. We sum both sides from some date $t_{0}$ until $t$ to get \begin{equation} \sum_{s=t_{0}}^{t}\Delta A_{s}=A_{t}-A_{t_{0}-1}=\sum_{s=t_{0}}^{t}q_{s}. \end{equation} We can then write the solution as \begin{eqnarray} \Delta A_{t} &=&q_{t} \\ &\rightarrow &A_{t}=\sum_{s=t_{0}}^{t}q_{s}+A_{t_{0}-1}. \label{difference1} \end{eqnarray} As we see, this is very much like integrals \begin{eqnarray} \frac{dA\left( t\right) }{dt} &=&q\left( t\right) \\ &\rightarrow &A\left( t\right) =\int_{t_{0}}^{t}q\left( s\right) ds+A_{t_{0}}, \end{eqnarray} and the relation between $q_{t}$ and $\sum_{s=t_{0}}^{t}q_{s}$ is the same as between $q\left( t\right)$ and its primitive, since \begin{equation} \Delta \sum_{s=t_{0}}^{t}q_{s}=q_{t}. \end{equation} \subsubsection{Forward and backward solutions} The part $\sum_{s=t_{0}}^{t}q_{t}$ of the RHS of (\ref{difference1}) is exogenous, i.e., independent of $A_{t}$. Sometimes, it converges in one or both directions. Then we can write the solutions in another way. Suppose the following limit exists \begin{equation} \lim_{T\rightarrow \infty }\sum_{s=-T}^{t}q_{s}\equiv \sum_{s=-\infty }^{t}q_{s}. \end{equation} Then, the other part of RHS of (\ref{difference1}) should also have a well-defined limit \begin{equation} \lim_{T\rightarrow \infty }A_{-T}\equiv \underline{A}, \end{equation} so that the solution is \begin{equation} A_{t}=\sum_{s=-\infty }^{t}q_{s}+\underline{A}. \label{backsol} \end{equation} Clearly, (\ref{backsol}) solves (\ref{differenceeq1}), \begin{equation} A_{t}-A_{t-1}=\sum_{s=-\infty }^{t}q_{s}+\underline{A}-\sum_{s=-\infty }^{t-1}q_{s}-\underline{A}=q_{t}. \end{equation}% If the solution in (\ref{backsol}) exists, it is called the \emph{backward solution. } Analogously, the limit \begin{equation} \lim_{T\rightarrow \infty }\sum_{s=t}^{T}q_{s}\equiv \sum_{s=t}^{\infty }q_{s}, \end{equation}% might exist, in which case \begin{equation} \lim_{T\rightarrow \infty }A_{T}\equiv \bar{A} \end{equation} also exists. Then, we have \begin{eqnarray} \lim_{T\rightarrow \infty }\left( A_{T}-A_{t}\right) &=&\sum_{s=t+1}^{\infty }q_{s}, \\ &\rightarrow &A_{t}=\bar{A}-\sum_{s=t+1}^{\infty }q_{s}, \end{eqnarray} which is called the \emph{forward solution}. Example. Suppose $q_{t}$ follows a simple $AR(1)$ process \begin{eqnarray} \Delta A_{t} &=&q_{t}, \\ q_{t} &=&rq_{t-1.} \end{eqnarray} If $r<1,$ we can use the forward solution and if $r>1$, the backward solution works. In the former case, \begin{eqnarray} \sum_{s=t+1}^{\infty }q_{s} &=&\sum_{s=0}^{\infty }q_{t+1}r^{s}=\frac{q_{t+1}% }{1-r}, \\ &\rightarrow &A_{t}=\bar{A}-\frac{q_{t+1}}{1-r}, \end{eqnarray}% where $\bar{A}$ is determined from, for example, an initial condition. In the latter case, \begin{eqnarray} \sum_{s=-\infty }^{t}q_{s} &=&\sum_{s=0}^{\infty }q_{t}r^{-s}=q_{t}\frac{r}{% r-1} \\ A_{t} &=&q_{t}\frac{r}{r-1}+\underline{A}. \end{eqnarray} \subsubsection{First order difference equations with constant coefficients} A first order difference equation with constant coefficients has the following form. \begin{equation} x_{t}-ax_{t-1}=c \label{differenceeq2} \end{equation} As we see, the LHS is not a pure difference, as in (\ref{differenceeq1}), so we cannot simply sum over $t$. Instead we rely on the following result. \begin{result} The general solution (the complete set of solutions) to a difference equation is the general solution to the homogeneous part plus any particular solution to the complete equation. The general solution to the homogeneous first order difference equation with coefficient $a$ can be written \begin{equation} x_{t}^{h}=Aa^{t} \end{equation} where $A$ is an arbitrary constant. \end{result} A particular solution is sometimes the steady state which exists for (\ref% {differenceeq2}) if $a\neq 1.$% \begin{eqnarray} x^{ss}-ax^{ss} &=&c \\ &\rightarrow &x^{ss}=\frac{c}{1-a}. \end{eqnarray} \subsubsection{Non-constant RHS} Now look at \begin{equation} x_{t}-ax_{t-1}=q_{t}. \label{differenceeq3} \end{equation} A way to solve (\ref{differenceeq3}) is to use $a^{-t}$ as the analogue to the integrating factor. Multiplying both sides by $a^{-t}$ yields \begin{equation} a^{-t}\left( x_{t}-ax_{t-1}\right) =a^{-t}q_{t}. \end{equation}% Now, we see that the LHS\ can be written $\Delta \left( a^{-t}x_{t}\right)$% , implying \begin{eqnarray} \Delta \left( a^{-t}x_{t}\right) &=&a^{-t}q_{t}, \\ &\Rightarrow &a^{-t}x_{t}=\sum_{s=t_{0}}^{t}a^{-s}q_{s}+A, \\ a^{-t}x_{t} &=&\sum_{s=t_{0}}^{t}a^{-s}q_{s}+A, \\ x_{t} &=&Aa^{t}+\sum_{s=t_{0}}^{t}a^{t-s}q_{s}. \end{eqnarray}% Again, we should verify that this satisfies our original difference equation. \subsubsection{Stable growth -- the Solow growth model} Often in macro, the variables in the model grow in a way that precludes the existence of a steady state. However, some transformation of the variables might possess a steady state. The simplest example of this is the Solow growth model \ We will first solve the model in trending variables and then do it in transformed (detrended) variables. As we know, the savings rate is exogenous and denoted $S$ and the labor supply $N_{t}$ follows \begin{equation} N_{t}=e^{g}N_{t-1}\approx \left( 1+g\right) N_{t-1}. \end{equation} There is one good, used for consumption and as capital, which follows the law-of-motion \begin{equation} K_{t+1}=Sf\left( K_{t},N_{t}\right) , \end{equation} where $f\left( K_{t},N_{t}\right)$ is a concave production function. Let us specify production as the CRS Cobb-Douglas function, \begin{eqnarray} f\left( K_{t},N_{t}\right) &=&N_{t}^{1-\alpha }K_{t}^{\alpha },\text{ } \label{Solow} \\ 1 &>&\alpha >0. \end{eqnarray} Letting lower case letters denote natural logarithms, we have \begin{eqnarray} k_{t+1} &=&s+\left( 1-\alpha \right) n_{t}+\alpha k_{t}, \\ k_{t+1}-\alpha k_{t} &=&s+\left( 1-\alpha \right) n_{t}, \end{eqnarray}% and \begin{eqnarray} \Delta n_{t} &=&g \\ &\rightarrow &n_{t}=\sum_{0}^{t}g+n=tg+n, \end{eqnarray}% for some constant $n.$ Here, we can guess on a particular solution of the same form as the RHS, i.e., \begin{eqnarray} k_{t} &=&k+tg_{k} \\ k+\left( t+1\right) g_{k}-\alpha \left( k+tg_{k}\right) &=&s+\left( 1-\alpha \right) n_{t} \\ k\left( 1-\alpha \right) +t\left( 1-\alpha \right) g_{k}+g_{k} &=&s+\left( 1-\alpha \right) n+t\left( 1-\alpha \right) g \end{eqnarray} For this to hold for all $t$, we need \begin{eqnarray} g_{k} &=&g, \\ k &=&\frac{s-g}{1-\alpha }+n. \end{eqnarray} giving a particular solution \begin{equation} k_{p,t}=\frac{s-g}{1-\alpha }+n+tg. \end{equation} This solution is in economics often called a \emph{balanced growth path}, in which all variables grow at the same rate. The complete solution is \begin{eqnarray} k_{t} &=&Aa^{t}+\frac{s-g}{1-\alpha }+n+tg, \label{balancedgrwth} \\ k_{t} &=&\left( k_{0}-\left( \frac{s-g}{1-\alpha }+n\right) \right) a^{t}+% \frac{s-g}{1-\alpha }+n+tg. \end{eqnarray} In this case, a convenient alternative is to define a new variable, capital per capita, $C_{t}\equiv K_{t}/N_{t}$ $\rightarrow$ $c_{t}=k_{t}-n_{t}$. Using this, and dividing (\ref{Solow}) by $N_{t}$, we get \begin{eqnarray} \frac{K_{t+1}}{N_{t}} &=&\frac{Sf\left( K_{t},N_{t}\right) }{N_{t}} \\ \frac{K_{t+1}N_{t+1}}{N_{t+1}N_{t}} &=&\frac{SN_{t}^{1-\alpha }K_{t}^{\alpha }}{N_{t}}= \\ C_{t+1}e^{g} &=&S\left( \frac{K_{t}}{N_{t}}\right) ^{\alpha } \\ c_{t+1}+g &=&s+\alpha c_{t} \\ c_{t+1}-\alpha c_{t} &=&s-g \\ &\rightarrow &c^{ss}=\frac{s-g}{1-\alpha } \\ c_{t} &=&\left( c_{0}-\frac{s}{1-\alpha }\right) \alpha ^{t}+\frac{s}{% 1-\alpha }, \end{eqnarray}% coinciding with the solution in (\ref{balancedgrwth}) but now expressed as steady state in capital per capita, rather than a balanced growth path for capital. We can also look at the log of output per capita, \begin{eqnarray} y_{t} &=&\alpha c_{t}, \\ y_{t} &=&\left( y_{0}-\alpha \frac{s}{1-\alpha }\right) \alpha ^{t}+\alpha \frac{s}{1-\alpha }. \end{eqnarray} This is the basis for the so-called growth regressions, pioneered by Barro, \begin{equation} y_{t}-y_{0}=-y_{0}\left( 1-\alpha ^{t}\right) +\frac{s}{1-\alpha }\alpha \left( 1-\alpha ^{t}\right) , \end{equation} where growth over a sample period $0$ through $t$ is seen to depend negatively on initial output and positively on savings. Also, \begin{eqnarray} y_{t}-y_{t-1} &=&\left( y_{0}-\alpha \frac{s}{1-\alpha }\right) \alpha ^{t}+\alpha \frac{s}{1-\alpha } \\ &&-\left( y_{0}-\alpha \frac{s}{1-\alpha }\right) \alpha ^{t-1}-\alpha \frac{% s}{1-\alpha } \\ &=&-\left( y_{0}-\alpha \frac{s}{1-\alpha }\right) \alpha ^{t-1}\left( 1-\alpha \right) \\ &=&\left( \alpha \frac{s}{1-\alpha }-y_{t-1}\right) \left( 1-\alpha \right) . \end{eqnarray} As we see, the growth rate (the log difference) of output per capita is a fraction $\left( 1-\alpha \right)$ of the difference between the steady state and current output per capita. Note that convergence is \emph{slower }% the larger is the capital's share of output. In the limit when $\alpha \rightarrow 1$, the model shows no convergence and has become an \emph{% endogenous growth} model, with \begin{eqnarray} \Delta k_{t} &=&s, \\ \Delta c_{t} &=&\Delta y_{t}=s-g. \end{eqnarray} \subsection{Linear difference equations of higher order} \subsubsection{Higher order homogeneous difference equations with constant coefficients} Consider the homogeneous difference equation \begin{equation} x_{t+n}+a_{1}x_{t+n-1}+...a_{n}x_{t}=0, \label{diffpoly} \end{equation} \begin{definition} The forward operator E is defined by \begin{equation} E^{s}x_{t}\equiv Ex_{t+s} \end{equation} where $s$ is any integer, positive or negative. \end{definition} We can then write (\ref{diffpoly}) in a condensed polynomial form \begin{equation} P\left( E\right) x_{t}=0. \end{equation} We then have to find the roots of the equation \begin{equation} P\left( r\right) =r^{n}+a_{1}r^{n-1}+...a_{n}=0 \end{equation} Each root contributes to the general solution with one term that is independent of the others, exactly as with differential equations. \begin{result} Let $r_{s}$ denote the roots to the polynomial $P\left( r\right)$, i.e., all solutions to $P(r)=0$. Let the first $k\geq 0$ roots be distinct and the remaining $l=n-k$ roots repeated. Then, the general solution to (\ref% {diffpoly}) is \begin{equation} x_{t}=c_{1}r_{1}^{t}+...+c_{k}r_{k}^{t}+r_{k+1}^{t}\left( c_{k+1}+tc_{k+2}+...+t^{l-1}c_{k+l}\right) \label{solutiondiffpoly} \end{equation} If there is more than one set of repeated roots, each set of size $m\geq 2$ contributes with $m$ linearly independent terms \begin{equation} r_{k+1}^{t}\left( c_{k+1}+tc_{k+2}+...+t^{m-1}c_{k+m}\right) . \end{equation} \end{result} In the case of complex roots, express the complex number in polar form% \footnote{% Note that the formula $\theta =\tan ^{-1}\frac{b}{a}$ is valid only for $% \theta \in \left( -\pi /2,\pi /2\right] .$ For example, if $r_{1,2}=-\frac{1% }{2}\pm \frac{1}{2}\sqrt{3}i$ we note that \begin{eqnarray*} -\frac{1}{2}+\frac{1}{2}\sqrt{3}i &=&e^{-i\frac{4}{3}\pi }, \\ -\frac{1}{2}-\frac{1}{2}\sqrt{3}i &=&e^{i\frac{4}{3}\pi }. \end{eqnarray*}% } \begin{eqnarray} r &=&a+bi=\left\vert r\right\vert \left( \cos \theta +i\sin \theta \right) , \\ \left\vert r\right\vert &=&\sqrt{a^{2}+b^{2},} \\ \theta &=&\tan ^{-1}\frac{b}{a}. \end{eqnarray} We then use the fact that \begin{eqnarray} e^{i\theta } &=&\left( \cos \theta +i\sin \theta \right) \\ r &=&\left| r\right| e^{i\theta }\rightarrow r^{t}=\left| r\right| ^{t}e^{i\theta t}=\left| r\right| ^{t}\left( \cos t\theta +i\sin t\theta \right) \end{eqnarray} to get for the complex conjugates $r_{1,2}=a\pm bi,$ \begin{equation} c_{1}r_{1}^{t}+c_{2}r_{2}^{t}=\left| r\right| ^{t}\left( \tilde{c}_{1}\cos t\theta +\tilde{c}_{2}\sin t\theta \right) \label{solutioncomplex} \end{equation} which we can see is a generalization of (\ref{solutiondiffpoly}). Complex roots thus give us oscillating solutions. \subsubsection{Stability} From the solutions (\ref{solutiondiffpoly}) and (\ref{solutioncomplex}), it is clear that roots such that $\left| r\right| <1,$give converging terms. \subsubsection{Higher order non-homogeneous difference equations with constant coefficients} \begin{result} The general solution to a non-homogeneous difference equation with constant coefficients is given by the general solution to the homogeneous part plus any solution to the full equation. \end{result} To solve non-homogeneous equations we thus have to find particular solution to the complete equation to add to the general solution of the homogeneous part. The simplest non-homogeneous difference equation is \begin{equation} P\left( x_{t}\right) =c. \end{equation} Here, we try a steady state \begin{eqnarray} P\left( x^{ss}\right) &=&P\left( 1\right) x^{ss}=c \\ x^{ss} &=&\frac{c}{P\left( 1\right) } \end{eqnarray} provided $P\left( 1\right) \neq 0.$ In the more general case we often have to guess a particular solution. \subsection{Systems of linear first order difference equations} Systems of first order difference equations are solved with the diagonalization method that we also used for the differential equation. \begin{eqnarray} \left[ \begin{array}{c} x_{1,t+1} \\ x_{2,t+1} \\ . \\ x_{n,t+1}% \end{array} \right] &=&\mathbf{A}\left[ \begin{array}{c} x_{1,t} \\ x_{2,t} \\ . \\ x_{n,t}% \end{array} \right] +\mathbf{P} \\ \mathbf{x}_{t+1} &=&\mathbf{Ax}_{t}+\mathbf{P.} \end{eqnarray} First we find the diagonalizing matrix of eigenvectors $\mathbf{B}^{-1}$. Then we solve the homogeneous equation by defining \begin{eqnarray} \mathbf{y}_{t+1} &=&\mathbf{Bx}_{t+1} \\ &=&\mathbf{BAx}_{t} \\ &=&\mathbf{BA\mathbf{B}^{-1}y}_{t} \\ &=&\left[ \begin{array}{cccc} r_{1}^{t} & 0 & . & 0 \\ 0 & r_{2}^{t} & . & 0 \\ . & . & . & . \\ 0 & 0 & . & r_{n}^{t}% \end{array}% \right] \left[ \begin{array}{c} c_{1} \\ c_{2} \\ . \\ c_{n}% \end{array}% \right] \\ &\equiv &\mathbf{r}^{t}\mathbf{c.} \end{eqnarray} Then we transform back\footnote{% Since the vector of constants, $\mathbf{c}$, is arbitrary, we may equally well write $\mathbf{x}_{t}=\mathbf{B}^{-1}\mathbf{r}^{t}\mathbf{c}$.} \begin{equation} \mathbf{x}_{t+1}=\mathbf{B}^{-1}\mathbf{y}_{t+1}=\mathbf{B}^{-1}\mathbf{r}% ^{t}\mathbf{c.} \label{Diff1} \end{equation} A particular solution to the complete equation then has to be added. Here we might try to find a steady state as a particular solution \begin{eqnarray} \mathbf{x}^{ss} &=&\mathbf{Ax}^{ss}+\mathbf{P,} \\ \mathbf{x}^{ss} &=&\left( \mathbf{I-A}\right) ^{-1}\mathbf{P.} \end{eqnarray} If we have the initial conditions we find that \begin{eqnarray} \mathbf{x}_{t+1} &=&\mathbf{B}^{-1}\mathbf{r}^{t}\mathbf{c+x}^{ss} \\ \mathbf{x}_{1} &=&\mathbf{B}^{-1}\mathbf{r}^{0}\mathbf{c+x}^{ss} \\ &=&\mathbf{B}^{-1}\mathbf{c+x}^{ss} \\ \mathbf{c} &=&\mathbf{B}\left( \mathbf{x}_{1}-\mathbf{x}^{ss}\right) \\ \mathbf{x}_{t+1} &=&\mathbf{B}^{-1}\mathbf{r}^{t}\mathbf{B}\left( \mathbf{x}% _{1}-\mathbf{x}^{ss}\right) \mathbf{+x}^{ss} \end{eqnarray} An example; \begin{eqnarray} \mathbf{x}_{t+1} &=&\left[ \begin{array}{cc} 0 & 1 \\ -1 & 0% \end{array}% \right] \mathbf{x}_{t}+\left[ \begin{array}{c} -1 \\ -1% \end{array}% \right] . \\ \mathbf{x}^{ss} &=&\left( \left[ \begin{array}{cc} 1 & 0 \\ 0 & 1% \end{array}% \right] -\left[ \begin{array}{cc} 0 & 1 \\ -1 & 0% \end{array}% \right] \right) ^{-1}\left[ \begin{array}{c} -1 \\ -1% \end{array}% \right] =\left[ \begin{array}{c} -1 \\ 0% \end{array}% \right] \end{eqnarray} Eigenvalues are \begin{eqnarray} r_{1,2} &=&i,-i\rightarrow r_{1,2}^{t}=1^{t}\left( \cos \left( \frac{\pi }{2}% t\right) \pm i\sin \left( \frac{\pi }{2}t\right) \right) \\ &=&\left( \cos \left( \frac{\pi }{2}t\right) \pm i\sin \left( \frac{\pi }{2}% t\right) \right) \end{eqnarray} and \begin{equation} \mathbf{B}^{-1}=\left[ \begin{array}{cc} 1 & 1 \\ i & -i% \end{array} \right] . \end{equation} The solution is \begin{gather} \mathbf{x}_{t}=\left[ \begin{array}{cc} 1 & 1 \\ i & -i% \end{array}% \right] \left[ \begin{array}{cc} i^{t} & 0 \\ 0 & \left( -i\right) ^{t}% \end{array}% \right] \left[ \begin{array}{c} c_{1} \\ c_{2}% \end{array}% \right] +\left[ \begin{array}{c} -1 \\ 0% \end{array}% \right] \\ =\left[ \begin{array}{cc} i^{t} & \left( -i\right) ^{t} \\ i^{t+1} & \left( -i\right) ^{t+1}% \end{array}% \right] \left[ \begin{array}{c} c_{1} \\ c_{2}% \end{array}% \right] +\left[ \begin{array}{c} -1 \\ 0% \end{array}% \right] \\ =\left[ \begin{array}{c} c_{1}i^{t}+c_{2}\left( -i\right) ^{t} \\ c_{1}i^{t+1}+c_{2}\left( -i\right) ^{t+1}% \end{array}% \right] +\left[ \begin{array}{c} -1 \\ 0% \end{array}% \right] \\ =\left[ \begin{array}{c} c_{1}\left( \left( \cos \left( \frac{\pi }{2}t\right) +i\sin \left( \frac{% \pi }{2}t\right) \right) \right) \\ +c_{2}\left( \cos \left( \frac{\pi }{2}t\right) -i\sin \left( \frac{\pi }{2}% t\right) \right) \\ c_{1}\left( \left( \cos \left( \frac{\pi }{2}\left( t+1\right) \right) +i\sin \left( \frac{\pi }{2}\left( t+1\right) \right) \right) \right) \\ +c_{2}\left( \cos \left( \frac{\pi }{2}\left( t+1\right) \right) -i\sin \left( \frac{\pi }{2}\left( t+1\right) \right) \right)% \end{array}% \right] \\ +\left[ \begin{array}{c} -1 \\ 0% \end{array}% \right] \\ =\left[ \begin{array}{c} \left( c_{1}+c_{2}\right) \left( \cos \left( \frac{\pi }{2}t\right) +\left( c_{1}-c_{2}\right) i\sin \left( \frac{\pi }{2}t\right) \right) \\ \left( c_{1}+c_{2}\right) \left( \cos \left( \frac{\pi }{2}\left( t+1\right) \right) +\left( c_{1}-c_{2}\right) i\sin \left( \frac{\pi }{2}\left( t+1\right) \right) \right)% \end{array}% \right] +\left[ \begin{array}{c} -1 \\ 0% \end{array}% \right] \\ \left[ \begin{array}{c} \tilde{c}_{1}\left( \cos \left( \frac{\pi }{2}t\right) +\tilde{c}_{2}\sin \left( \frac{\pi }{2}t\right) \right) \\ \tilde{c}_{1}\left( \cos \left( \frac{\pi }{2}\left( t+1\right) \right) +% \tilde{c}_{2}\sin \left( \frac{\pi }{2}\left( t+1\right) \right) \right)% \end{array}% \right] +\left[ \begin{array}{c} -1 \\ 0% \end{array}% \right] . \end{gather} \subsubsection{Non-invertible eigenvectors} If some roots are repeated $\mathbf{B}^{-1}$ may be non-invertible. In this case we cannot use the diagonalization method. Instead we can use the existence of a higher order single difference equation that is equivalent to the system of first order difference equations we want to solve. This is exactly analogous to the case of differential equations. Look at the following example: \begin{equation} \mathbf{x}_{t+1}=\left[ \begin{array}{cc} -1 & 1 \\ -4 & 3% \end{array}% \right] \mathbf{x}_{t} \label{differencesystem1} \end{equation} The eigenvalues of the coefficient matrix are both 1 and the matrix of eigenvectors are non-invertible. Now we transform the system into a second order single difference equation. From the first row we have that \begin{eqnarray} x_{1,t+1} &=&-x_{1,t}+x_{2,t} \label{subst} \\ &\rightarrow &x_{2,t}=x_{1,t+1}+x_{1,t}, \\ x_{2,t+1} &=&x_{1,t+2}+x_{1,t+1} \end{eqnarray} Using the second row, \begin{eqnarray} x_{2,t+1} &=&-4x_{1,t}+3x_{2,t} \label{equivsystem} \\ &\rightarrow &x_{1,t+2}+x_{1,t+1}=-4x_{1,t}+3\left( x_{1,t+1}+x_{1,t}\right) \\ 0 &=&x_{1,t+2}-2x_{1,t+1}+x_{1,t} \end{eqnarray} Note that the polynomial associated with this second order difference equation is identical to the characteristic equation of the coefficient matrix in (\ref{differencesystem1}). Consequently, they have the same (repeated) root 1. The solution to (\ref{equivsystem}) is \begin{equation} x_{1,t}=\left( c_{1}+tc_{2}\right) 1^{t}=\left( c_{1}+tc_{2}\right) \end{equation} With knowledge of $x_{1,0}$ and $x_{1,1}$, we have \begin{equation} x_{1,t}=\left( x_{1,0}+t\left( x_{1,1}-x_{1,0}\right) \right) \end{equation} By using (\ref{subst}), we can express the solution in terms of $x_{1,t}$ and $x_{2,t}$% \begin{eqnarray} x_{1,t} &=&x_{1,0}+t\left( x_{1,1}-x_{1,0}\right) \\ x_{2,t} &=&x_{1,t+1}+x_{1,t} \\ &=&2x_{1,0}+\left( 2t+1\right) \left( x_{1,1}-x_{1,0}\right) \\ &\rightarrow &\left[ \begin{array}{c} x_{1,t} \\ x_{2,t}% \end{array}% \right] =\left[ \begin{array}{cc} 1 & t \\ 2 & 2t+1% \end{array}% \right] \left[ \begin{array}{c} x_{1,0} \\ x_{1,1}-x_{1,0}% \end{array}% \right] . \end{eqnarray}% \newpage \section{Dynamic Optimization in Discrete Time} \subsection{Non-Stochastic Dynamic Programming} %TCIMACRO{\TeXButton{TeX field}{\setcounter{equation}{0}}}% %BeginExpansion \setcounter{equation}{0}% %EndExpansion Consider the dynamic problem \begin{eqnarray} &&\max_{\left\{ k_{t},c_{t}\right\} _{t=1}^{T}}\sum_{t=1}^{T}u\left( k_{t},c_{t},t\right) \label{Dynoptproblem1} \\ \text{s.t.\thinspace }k_{1} &=&\underline{k} \label{constraint1} \\ k_{t+1} &=&f\left( k_{t},c_{t},t\right) ,t=1...,T, \\ k_{T+1} &=&\bar{k} \label{constraint3} \\ c_{t} &\in &\Omega _{t} \end{eqnarray} Before trying to solve this, notice \begin{enumerate} \item Per period payoff is additive over time. \item $k_{t}$ cannot be changed in period $t$, but its future values, its law-of-motion can be changed by $c_{t}.$ We will call $k$ a \emph{state variable} (to be more properly defined later) and $c$ a \emph{control variable }that must be chosen in the control set $\Omega _{t}$. A sequence $% k_{1},k_{2},...,k_{T+1}$ is said to be \emph{admissible} if and only if it satisfies the constraints (\ref{constraint1}-\ref{constraint3}) for some sequence $\left\{ c_{t}\in \Omega _{t}\right\} _{t=1}^{T}$. \end{enumerate} A direct way to solve this would be to form the Lagrangian \begin{eqnarray} L &=&\sum_{t=1}^{T}u\left( k_{t},c_{t},t\right) +\sum_{t=1}^{T}\lambda _{t}\left( f\left( k_{t},c_{t},t\right) -k_{t+1}\right) \\ &&+\mu _{1}\left( \underline{k}-k_{1}\right) +\mu _{T+1}\left( k_{T+1}-\bar{k% }\right) \end{eqnarray} with first order conditions \begin{eqnarray} u_{k}\left( k_{t},c_{t},t\right) +\lambda _{t}f_{k}\left( k_{t},c_{t},t\right) -\lambda _{t-1} &=&0,\forall t=2,..,T, \\ u_{c}\left( k_{t},c_{t},t\right) +\lambda _{t}f_{c}\left( k_{t},c_{t},t\right) &=&0,\forall t=1,..,T, \\ u_{k}\left( k_{1},c_{1},1\right) +\lambda _{1}f_{k}\left( k_{1},c_{1},1\right) -\mu _{1} &=&0, \\ -\lambda _{T}+\mu _{T+1} &=&0, \end{eqnarray} and (\ref{constraint1}-\ref{constraint3}). This works, at least in principle, if $T$ is finite. An alternative way is to recognize that in a problem like this, each sub-section of the path must be optimal in itself. This means that the problem has a recursive formulation i.e., it can be set up sequentially. We can thus solve the problem backwards starting from the last period. In any period, the remaining problem only depends on earlier actions through the ''inherited'' value of $k$. For example, if the problem is over three periods ($T=3$) we can rewrite (% \ref{Dynoptproblem1}) \begin{gather} \max_{c_{1},k_{2}|_{k_{1}}}\left( u\left( k_{1},c_{1},1\right) +\max_{c_{2},k_{3}|_{k_{2}}}\left( u\left( k_{2},c_{2},2\right) +\max_{c_{3},k_{4}|_{k_{3}}}u\left( k_{3},c_{3},3\right) \right) \right) \\ \text{s.t. }k_{t+1}=f\left( k_{t},c_{t},t\right) ,t=1...,3 \\ k_{4}=\bar{k} \end{gather} In the final period ($T=3$), the problem is trivial; simply set $c_{3}$ so that $k_{4}=\bar{k}$. The value of $c_{3}$ that solves $\bar{k}=f\left( k_{3},c_{3},3\right)$ is a function of $k_{3}$.\footnote{% For now, we just assume it is unique.} Denote that function $c_{3}\left( k_{3}\right) .$ We can then define \begin{equation} u\left( k_{3},c_{3}\left( k_{3}\right) ,3\right) \equiv W\left( k_{3},3\right) . \end{equation} The interpretation of $W\left( k_{3},3\right)$, is the \emph{maximum remaining pay-off in period 3,} being a function of the state variable $% k_{3}$. In period 2, we then want to solve \begin{eqnarray} &&\max_{c_{2},k_{3}}\left( u\left( k_{2},c_{2},2\right) +u\left( k_{3},c_{3}\left( k_{3}\right) ,3\right) \right) \\ \text{s.t.\thinspace }k_{3} &=&f\left( k_{2},c_{2},2\right) . \end{eqnarray} Using $W\left( k_{3},3\right) ,$we can write this \begin{equation} \max_{c_{2}}\left( u\left( k_{2},c_{2},2\right) +W\left( f\left( k_{2},c_{2},2\right) ,3\right) \right) \end{equation}% The solution and the maximized value depends on $k_{2}$ only and we called the latter the value function and $k_{2}$ the state variable. We can then define \begin{equation} W\left( k_{2},2\right) \equiv \max_{c_{2}}\left( u\left( k_{2},c_{2},2\right) +W\left( f\left( k_{2},c_{2},2\right) ,3\right) \right) . \end{equation} The interpretation of this Bellman equation is straightforward. It says that the maximum remaining pay-off in period 2, being a function of $k_{2}$, is identically (i.e., for all $k_{2}$) equal to the maximum over the control in period 2, $c_{2},$ over period 2 pay-off and the maximum remaining pay-off in period 3 with period 3 state variable given by $f\left( k_{2},c_{2},2\right) .$ The trade-off between generation of current and future pay-off is optimized by \emph{one} simple FOC \begin{equation} u_{c}\left( k_{2},c_{2},2\right) +W_{k}\left( f\left( k_{2},c_{2},2\right) ,3\right) f_{c}\left( k_{2},c_{2},2\right) =0. \end{equation} Finally, in the first period, \begin{equation} W\left( k_{1},1\right) \equiv \max_{c_{1}}\left( u\left( k_{1},c_{1},1\right) +W\left( f\left( k_{1},c_{1},1\right) ,2\right) \right) . \end{equation} If we know the value functions, the multidimensional problem has become much simpler. The Bellman equation provides a way of verifying that the value function we use is correct. It is of course straightforward to extend the analysis to any finite horizon problem, yielding \begin{eqnarray} W\left( k_{t},t\right) &\equiv &\max_{c_{t}}\left( u\left( k_{t},c_{t},t\right) +W\left( k_{t+1},t+1\right) \right) \\ \text{s.t. }k_{t+1} &=&f\left( k_{t},c_{t},t\right) \end{eqnarray} \subsubsection{Discounting and the Current Value Bellman equation} Very often in macroeconomics, the objective function is a discounted sum of pay-offs, i.e., (\ref{Dynoptproblem1}) can be written% \begin{equation} \max_{\left\{ c\right\} _{t=1}^{T}}\sum_{t=1}^{T}\beta ^{t-1}u\left( k_{t},c_{t},t\right) . \end{equation}% In this case, it is convenient to work with \emph{current value functions, }$% V\left( k,t\right)$ which should be interpreted as the maximum remaining value that can be achived from time $t$ and onward, given $k_{t}$ and \emph{% seen from period }$t.$ In other words, given a problem \begin{eqnarray} &&\max_{\left\{ c\right\} _{t=1}^{T}}\sum_{t=1}^{T}\beta ^{t-1}u\left( k_{t},c_{t},t\right) \\ \text{s.t.~}k_{t+1} &=&f\left( k_{t},c_{t},t\right) ,t=1...,T \\ k_{1} &=&\underline{k} \\ k_{T+1} &=&\bar{k} \notag \end{eqnarray}% we define for any $t\in \left\{ 1,...,T\right\}$% \begin{eqnarray} V\left( k_{t},t\right) &\equiv &\max_{\left\{ c,k\right\} _{s=t}^{T}}\sum_{s=t}^{T}\beta ^{s-t}u\left( k_{s},c_{s},s\right) \\ \text{s.t.~}k_{t+1} &=&f\left( k_{t},c_{t},t\right) ,t=s...,T. \\ k_{T+1} &=&\bar{k}. \notag \end{eqnarray} From this follows that \begin{equation*} V\left( k,t\right) \equiv \beta ^{-\left( t-1\right) }W\left( k,t\right) . \end{equation*} Using this in the Bellman equation, (where I substitute for $k_{t+1}$ from the law-of-motion) we get the \emph{current value Bellman equation}% \begin{eqnarray} W\left( k_{t},t\right) &\equiv &\max_{c_{t}}\left( \beta ^{t-1}u\left( k_{t},c_{t},t\right) +W\left( f\left( k_{t},c_{t},t\right) ,t+1\right) \right) \\ \beta ^{t-1}V\left( k,t\right) &\equiv &\max_{c_{t}}\left( \beta ^{t-1}u\left( k_{t},c_{t},t\right) +\beta ^{t}V\left( f\left( k_{t},c_{t},t\right) ,t+1\right) \right) \\ V\left( k,t\right) &\equiv &\max_{c_{t}}\left( u\left( k_{t},c_{t},t\right) +\beta V\left( f\left( k_{t},c_{t},t\right) ,t+1\right) \right) \label{CurrentBellman} \end{eqnarray}% In practice, the current value Bellman equation is the most used variant in macroeconomics and, therefore, you will often see the word \emph{current }% dropped \ and (\ref{CurrentBellman}) is simple referred to as the Bellman equation and $V\left( k,t\right)$ is referred to as the value function. \subsubsection{Infinite Horizon and Autonomous Problems} In an infinite horizon problem we cannot use the method of starting from the last period. Still, if the problem has a well-defined value function, \emph{% it satisfies the Bellman equation}. Furthermore, under some conditions (which we will return to later), there is\emph{\ only one} function that solves the Bellman equation, so if we find one function that solves the Bellman equation, we have a solution to the dynamic optimization problem. Since geometric discounting will prove to be important for showing uniqueness, we will use that from now on. An important case when there is necessarily exists a unique solution to Bellman equation is if there is strict discounting ($\beta \in \left( 0,1\right)$) and $u\left( k_{t},c_{t},t\right)$ is bounded for any admissible sequence of $\left\{ k_{t},c_{t}\right\} .$ To find a solution, we will use two different approaches. \begin{enumerate} \item Guess on a value function and verify that it satisfies the Bellman equation. \item Iterate on the Bellman equation until it converges. \end{enumerate} Guessing is often feasible when the problem is autonomous (stationary). Then, the problem is independent of time in the sense that given an initial condition on the state variable(s), the solution and the maximized objective is independent of the starting date. A problem is autonomous if \begin{enumerate} \item Time is infinite, \item the law of motion for the state (including constraints on the control) is independent of time, and \item the per-period return function is the same over time, except for possibly a geometric discount factor, i.e., $u\left( k,c,t\right) =\beta ^{t}u\left( k,c\right)$. \end{enumerate} (Think about what would happen if any of these conditions is not satisfied). In this case, the \emph{current }value function turns out to be independent of time.\footnote{% The present value function $W(k,t)$ is not independent ot time, but separable so that we can write $W\left( k,t\right) =W\left( k\right) \beta ^{t-1}.$} We can then write the current Bellman equation in terms of the current value function \begin{eqnarray} V\left( k_{t}\right) &\equiv &\max_{c_{t}}\left( u\left( k_{t},c_{t}\right) +\beta V\left( k_{t+1}\right) \right) \label{Bellman1} \\ \text{s.t. }k_{t+1} &=&f\left( k_{t},c_{t}\right) . \end{eqnarray} Suppose we find a solution to the maximization problem in the RHS\ of (\ref% {Bellman1}). This will be a time invariant function $c(k)$, since $u,$ $V$ and $f$ are time-invariant. Plugging $c\left( k\right)$ into (\ref{Bellman1}% ), we get rid of the max-operator: \begin{equation} V\left( k\right) =u\left( k,c\left( k\right) \right) +\beta V\left( f\left( k,c\left( k\right) \right) \right) \label{Bellmanidentity1} \end{equation} If (\ref{Bellmanidentity1}) is satisfied for \emph{all values of admissible }% $\emph{k,}$ we have a solution to the value function, otherwise our guess was incorrect. Note that (\ref{Bellmanidentity1}) is a functional equation, i.e., the LHS\ and RHS have to be the same functions. It is convenient to define the RHS\ as \begin{equation} T(V\left( k\right) )\equiv \max_{c}u\left( k,c\right) +\beta V\left( f\left( k,c\right) \right) \end{equation}% where $T$ operates on \emph{functions} rather than on values.\footnote{% See section \ref{Tmap} for more on the $T$-mapping.} In the autonomous case, when the value function is unchanged over time, the Bellman equation defines a \emph{fixed point} for $T$ in the space of functions $V$; \begin{equation} V\left( k\right) =T\left( V\left( k\right) \right) . \label{Fixedpoint} \end{equation} (\ref{Fixedpoint}) means that if we plug in some function of $k$ in the RHS of (\ref{Fixedpoint}) we must get out the same function on the LHS. The Bellman equation is a necessary condition for $V\left( k\right)$ being a correctly specified value function, we will later discuss conditions under which it is also sufficient. Typically the value function and the objective function are of similar form. This is intuitive in the light of (\ref{Bellmanidentity1}). For example, if the $u$ function is logarithmic, we guess that the value function is of the logarithmic form too. For HARA utility functions (e.g., CRRA, CARA and quadratic) the value functions are generally of the same type as the utility function (Merton, 1971). \subsubsection{An example of guessing} In the problem (\ref{Dynoptproblem1}) let time be infinite and \begin{eqnarray} u(k,c,t) &=&\beta ^{t}\ln (c), \\ f(k,c,t) &=&k^{\alpha }-c,0<\alpha <1, \end{eqnarray} and let the end-condition $k_{T+1}=\bar{k}$ be replaced by $k_{t}>0$ $% \forall t.$ This is an autonomous problem, so we have \begin{eqnarray} V\left( k_{t}\right) &=&\max_{c,k}\left( u\left( c_{t}\right) +\beta V\left( k_{t+1}\right) \right) \\ \text{s.t. }k_{t+1} &=&f\left( k_{t},c_{t}\right) , \\ &\Rightarrow &V\left( k_{t}\right) =\max_{c}\left( \ln c_{t}+\beta V\left( k_{t}^{\alpha }-c_{t}\right) \right) \end{eqnarray} Now, guess that $V$ is of the same form as $u$, here $X\ln k_{t}+Y$, for some unknown constants $X$ and $Y,$ giving first order conditions \begin{eqnarray} u^{\prime }\left( c_{t}\right) &=&\beta V^{\prime }\left( k_{t}^{\alpha }-c_{t}\right) \\ \frac{1}{c_{t}} &=&\beta \frac{X}{k_{t}^{\alpha }-c_{t}} \\ &\Rightarrow &c_{t}=\frac{k_{t}^{\alpha }}{1+\beta X}\equiv c\left( k_{t}\right) . \\ k_{t+1} &=&k_{t}^{\alpha }-c_{t}=k_{t}^{\alpha }-\frac{k_{t}^{\alpha }}{% 1+\beta X}=\frac{\beta X}{1+\beta X}k_{t}^{\alpha } \end{eqnarray} Plugging $c\left( k_{t}\right)$ into the Bellman equation yields \begin{eqnarray} X\ln k_{t}+Y &=&\ln \frac{k_{t}^{\alpha }}{1+\beta X}+\beta \left( X\ln \frac{\beta X}{1+\beta X}k_{t}^{\alpha }+Y\right) \\ &=&\ln \frac{k_{t}^{\alpha }}{1+\beta X}+\beta \left( X\ln \frac{% k_{t}^{\alpha }}{1+\beta X}+X\ln \beta X+Y\right) \\ &=&\left( 1+\beta X\right) \ln \frac{k_{t}^{\alpha }}{1+\beta X}+\beta \left( X\ln \beta X+Y\right) \\ &=&\alpha \left( 1+\beta X\right) \ln k_{t}-\left( 1+\beta X\right) \ln \left( 1+\beta X\right) \\ &&+\beta X\ln \beta X+\beta Y. \notag \end{eqnarray} This is true for all values of $k$, if and only if \begin{eqnarray} X &=&\alpha \left( 1+\beta X\right) \\ Y &=&-\left( 1+\beta X\right) \ln \left( 1+\beta X\right) +\beta X\ln \beta X+\beta Y. \end{eqnarray} giving \begin{eqnarray} X &=&\frac{\alpha }{1-\alpha \beta } \\ Y &=&\frac{\beta \alpha \ln \left( \beta \alpha \right) +\left( 1-\beta \alpha \right) \ln \left( 1-\alpha \beta \right) }{\left( 1-\beta \right) \left( 1-\alpha \beta \right) } \\ &\Rightarrow &V\left( k\right) =\frac{\alpha }{1-\alpha \beta }\ln k+\frac{% \beta \alpha \ln \left( \beta \alpha \right) +\left( 1-\beta \alpha \right) \ln \left( 1-\alpha \beta \right) }{\left( 1-\beta \right) \left( 1-\alpha \beta \right) }. \end{eqnarray} Having $V\left( k\right) ,$ it is easy to find the \emph{optimal control}, or the \emph{policy rule}, \begin{eqnarray} c\left( k_{t}\right) &=&\frac{k_{t}^{\alpha }}{1+\beta \frac{\alpha }{% 1-\alpha \beta }}=\left( 1-\alpha \beta \right) k_{t}^{\alpha }, \\ \ln k_{t+1} &=&\ln \alpha \beta +\alpha \ln k_{t.} \end{eqnarray} \subsubsection{Iteration} As mentioned above, an alternative way to find the infinite horizon value function is to find the limit of finite horizon Bellman equation as the horizon goes to infinity. Under for our purposes quite general conditions this limit exists and is equal to the value function for the infinite horizon problem. Let us change notation slightly, measuring time as the number of remaining periods $s$ until the final period. We then denote the (current) value function with $s$ periods left by \begin{equation} V(k,s). \end{equation}% \qquad\ We also assume geometric discounting and that both pay-offs and the law-of-motion for the state variable are time-independent ($u(k,c,s)=u\left( k,c\right)$ and $f\left( k,c,s\right) =f\left( k,c\right)$) so that the infinite horizon problem is autonomous. If the following limit is well-defined, we denote \begin{equation} \lim {}_{s\rightarrow \infty }V(k,s)\equiv V\left( k\right) . \end{equation} The iteration method is usually done numerically, but it can (at some cost of messiness) be done also analytically. Using the $T$ operator and the Bellman equation, we find the limit in the following way \begin{eqnarray} V\left( k,s\right) &=&\max_{c_{s}}\left( u\left( k_{s},c_{s}\right) +\beta V\left( f\left( k_{s},c_{s}\right) ,s-1\right) \right) \\ &\equiv &T\left( V\left( k,s-1\right) \right) \\ V\left( k,s\right) &=&T^{s}V\left( k,0\right) \\ V\left( k\right) &\equiv &\lim_{s\rightarrow \infty }V\left( k,s\right) =\lim_{s\rightarrow \infty }T^{s}V\left( k,0\right) . \end{eqnarray} If the limit exists, it clearly satisfies the Bellman equation \begin{eqnarray} V\left( k\right) &=&T\left( V\left( k\right) \right) \\ \lim_{s\rightarrow \infty }T^{s}V\left( k,0\right) &=&T\left( \lim_{s\rightarrow \infty }T^{s}V\left( k,0\right) \right) =\lim_{s\rightarrow \infty }T^{s+1}V\left( k,0\right) . \end{eqnarray} The remaining issue is what function $V\left( k,0\right)$ to start the iteration with. However, suppose that we can show that the limit $V\left( k\right)$ satisfies \begin{equation} \lim_{s\rightarrow \infty }\beta ^{s}V\left( k_{s}\right) =0 \label{iterationlimit} \end{equation}% for ALL\ admissable\ values of $k_{t+s}$ that can be reached, given relevant initial conditions and other constraints. Then, it can easily be shown that the Bellman equation is \emph{sufficient. }Then, $\lim_{s\rightarrow \infty }T^{s}V\left( k,0\right)$ provides a unique solution to the Bellman equation. This means that the limit is independent of the choice of $V\left( k,0\right) .$ As we see from (\ref{iterationlimit}), $\beta <1$ and $V\left( k\right)$ bounded are sufficient for uniqueness. Intuitively, if $\beta <1$ and pay-offs are bounded, the pay-off in the infinite horizon has no impact on the value function. Let us revert to measure time in the usual way. Then, if the Bellman equation is satisfied, we have \begin{gather} V\left( k_{t}\right) =\max_{c_{t}}\left( u\left( k_{t},c_{t}\right) +\beta V\left( k_{t+1}\right) \right) \\ \text{s.t. }k_{t+1}=f\left( k_{t},c_{t}\right) \notag \\ \max_{c_{t}}\left( u\left( k_{t},c_{t}\right) +\beta \left( \max_{c_{t+1}}\left( u\left( k_{t+1},c_{t+1}\right) +\beta V\left( k_{t+2}\right) \right) \right) \right) \\ \text{s.t.}k_{t+1}=f\left( k_{t},c_{t}\right) ,k_{t+2}=f\left( k_{t+1},c_{t+1}\right) \notag \\ =\max_{\{c_{t+n}\}_{0}^{1}}\sum_{n=0}^{1}\beta ^{n}u\left( k_{t+n},c_{t+n}\right) +\beta ^{2}V\left( f\left( k_{t+1},c_{t+1}\right) \right) , \\ \text{s.t.}k_{t+1}=f\left( k_{t},c_{t}\right) ,k_{t+2}=f\left( k_{t+1},c_{t+1}\right) \end{gather} Repeating this, and taking the limit yields \begin{eqnarray} V\left( k_{t}\right) &=&\max_{\{c_{t+n}\}_{0}^{s}}\sum_{n=0}^{s}\beta ^{n}u\left( k_{t+n},c_{t+n}\right) +\beta ^{s+1}V\left( k_{t+n+1}\right) , \\ \text{s.t.}k_{t+n+1} &=&f\left( k_{t+n},c_{t+n}\right) \forall n\in \left\{ 0,s\right\} \\ V\left( k_{t}\right) &=&\max_{\{c_{t+n}\}_{0}^{\infty }}\sum_{n=0}^{\infty }\beta ^{n}u\left( k_{t+n},c_{t+n}\right) +\lim_{s\rightarrow \infty }\beta ^{s+1}V\left( k_{t+s+1}\right) , \\ \text{s.t.}k_{t+n+1} &=&f\left( k_{t+n},c_{t+n}\right) \forall n\geq 0, \end{eqnarray}% Now, if $\lim_{s\rightarrow \infty }\beta ^{s+1}V\left( k_{t+s+1}\right) =0,$ for all permissible paths of $k$ we have showed that \begin{eqnarray} V\left( k_{t}\right) &=&\max_{\{c_{t+s}\}_{0}^{\infty }}\sum_{s=0}^{\infty }\beta ^{s}u\left( k_{t+s},c_{t+s}\right) \\ \text{s.t. }k_{t+s+1} &=&f\left( k_{t+s},c_{t+s}\right) \forall s\geq 0\text{% , given}k_{t}, \end{eqnarray}% i.e., that the Bellman equation implies optimality. The iteration can easily be done numerically, either by specifying a functional form if we know that, or by just choosing a grid. In the latter case we assume that the state variable must belong to a finite set of values, say for example that in every period $k$ must be chosen from the set $\mathbf{K\equiv }\left\{ k_{1},k_{2},...,k_{n}\right\}$. Then, we can compute the corresponding set of possible controls, $c_{m,n}\in \mathbf{C}$ from the equation \begin{equation} k_{m}=f\left( k_{n},c_{m,n}\right) . \end{equation} Then, in each iteration, we solve the Bellman equation for each $k\in \mathbf{K}$, giving for iteration $s$% \begin{equation} V\left( k_{n},s\right) =\max_{c_{m,n}\in \mathbf{C}}\left( u\left( k_{n},c_{m,n}\right) +\beta V\left( k_{m},s-1\right) \right) . \end{equation} This goes quickly on a computer and the iteration is repeated until $V\left( k,s\right)$ is sufficiently close to $V\left( k,s-1\right)$ over the set of $k\in \mathbf{K.}$ \subsubsection{An envelope result} We will later have use for the following envelope result, which implies that we evaluate $dV\left( k\right) /dk\,$as the partial derivative holding $c$ constant. To see this, note that \begin{eqnarray} \frac{dV\left( k_{t}\right) }{dk_{t}} &\equiv &V^{\prime }\left( k_{t}\right) =\frac{\partial u\left( k_{t},c_{t}\right) }{\partial k_{t}}% +\beta V^{\prime }\left( k_{t+1}\right) \frac{\partial f\left( k_{t},c_{t}\right) }{\partial k_{t}} \\ &&+\frac{dc_{t}}{dk_{t}}\left( \frac{\partial u\left( k_{t},c_{t}\right) }{% \partial c_{t}}+\beta V^{\prime }\left( k_{t+1}\right) \frac{\partial f\left( k_{t},c_{t}\right) }{\partial c_{t}}\right) . \end{eqnarray}% However, the Bellman equation implies that the last term above is zero if $c$ is chosen optimally, either because the first-order condition for an interior maximum is satisfied $\frac{\partial u\left( k,c\right) }{\partial c% }+\beta V^{\prime }\left( k_{t+1}\right) \frac{\partial f\left( k_{t},c_{t}\right) }{\partial c_{t}}=0$, or in the case of a corner because $% \frac{dc_{t}}{dk_{t}}=0.$ This envelope result is often very useful. Consider the example \begin{eqnarray} u\left( k_{t},c_{t}\right) &=&v\left( c_{t}\right) , \\ f\left( k_{t},c_{t}\right) &=&y\left( k_{t}\right) -c_{t}. \end{eqnarray} The interior solution to the Bellman equation satisfies \begin{eqnarray} 0 &=&v^{\prime }\left( c_{t}\right) +\beta V^{\prime }\left( k_{t+1}\right) \frac{\partial f\left( k_{t},c_{t}\right) }{\partial c_{t}}, \label{consopt} \\ &\rightarrow &v^{\prime }\left( c_{t}\right) =\beta V^{\prime }\left( k_{t+1}\right) . \end{eqnarray} The envelope condition yields \begin{eqnarray} V^{\prime }\left( k_{t}\right) &=&\frac{\partial u\left( k_{t},c_{t}\right) }{\partial k_{t}}+\beta V^{\prime }\left( k_{t+1}\right) \frac{\partial f\left( k_{t},c_{t}\right) }{\partial k_{t}} \\ &=&v^{\prime }\left( c_{t}\right) y^{\prime }\left( k_{t}\right) . \\ V^{\prime }\left( k_{t+1}\right) &\rightarrow &v^{\prime }\left( c_{t+1}\right) y^{\prime }\left( k_{t+1}\right) . \end{eqnarray} Using this in (\ref{consopt}) yields the \emph{Euler equation} \begin{equation} v^{\prime }\left( c_{t}\right) =\beta v^{\prime }\left( c_{t+1}\right) y^{\prime }\left( k_{t+1}\right) . \end{equation} \subsubsection{State Variables} We often solve the dynamic programming problem by guessing a form of the value function. The first thing to determine is then which variables should enter, i.e., which variables are the state variables. The state variables must satisfy both following conditions: 1. To enter the value function at time they must be realized at $t$. Note, however, that it sometimes may be convenient to use a conditional expectation $E_{t}(z_{t+s})$ as a state variable. The expectation as of $t$ is certainly realized at $t$ even if the stochastic variable $z_{t+s}$ is not realized. 2. The set of variables chosen as state variables must together give sufficient information so that the value of the program from $t$ and onwards when the optimal control is chosen can be calculated.\footnote{% Can you figure out what do we need if the per period utility function in (% \ref{Dynoptproblem1}) were $u\left( c_{t},c_{t-1}\right) ?$} Note, we should try to find the smallest such set. If, for example we have an investment problem with several assets to invest in and without any costs of adjusting the portfolio, total wealth may be a sufficient as a state variable. \subsection{Stochastic Dynamic Programming} As long as the recursive structure of the problem is intact adding a stochastic element to the transition equation does not change the Bellman equation. Consider the problem \begin{eqnarray} &&\max_{\{c_{t}\}_{0}^{\infty }}E_{0}\sum_{t=0}^{\infty }\beta ^{t}u\left( k_{t},c_{t}\right) \\ \text{s.t. }k_{t+1} &=&f\left( k_{t},c_{t},\varepsilon _{t+1}\right) \forall t\geq 0,k_{0}\,\text{given, and} \\ k_{t} &>&0\text{ }\forall t\text{.} \end{eqnarray}% where $E_{0}$ is the expectations operator, conditional on time $0$ information and we assume that $c_{t}$ can be chosen conditional on information about $\varepsilon _{s}$ for all $s\leq t.$ Furthermore, let us assume that the distribution of $\varepsilon _{t}$ is $i.i.d.$ over time. Then, the Bellman equation becomes \begin{equation} V\left( k_{t}\right) =\max_{c_{t}}\left( u\left( k_{t},c_{t}\right) +\beta E_{t}V\left( f\left( k_{t},c_{t},\varepsilon _{t+1}\right) \right) \right) , \end{equation}% with a first-order condition \begin{equation} 0=u_{c}\left( k_{t},c_{t}\right) +\beta E_{t}\left( V^{\prime }\left( k_{t+1}\right) f_{c}\left( k_{t},c_{t},\varepsilon _{t+1}\right) \right) . \end{equation} Note that, in general $E_{t}\left( V^{\prime }\left( k_{t+1}\right) f_{c}\left( k_{t},c_{t},\varepsilon _{t+1}\right) \right) \neq E_{t}V^{\prime }\left( k_{t+1}\right) E_{t}f_{c}\left( k_{t},c_{t},\varepsilon _{t+1}\right)$. \subsubsection{A Stochastic Consumption Example} Consider the following problem \begin{eqnarray} &&\max_{\{c_{t}\}_{0}^{\infty }}E_{0}\sum_{t=0}^{\infty }\beta ^{t}\ln c_{t} \\ \text{s.t. }A_{t+1} &=&\left( A_{t}-c_{t}\right) \left( 1+\tilde{m}% _{t+1}\right) \forall t\geq 0, \\ A_{t} &\geq &0\text{, }\forall t\geq 0\text{, }A_{0}\,\text{given.} \end{eqnarray} The consumer decides how much to consume each period. The savings is placed in a risky asset with gross return $\left( 1+\tilde{m}_{t+1}\right)$, that is drawn from an i.i.d. distribution with $E\left( \ln \left( 1+\tilde{m}% _{t+1}\right) \right) =$ $m>-\infty$. If, for example the gross return is log-normal, with mean $\bar{m}$ and variance $\sigma ^{2}$ then, $E\left( 1+% \tilde{m}_{t+1}\right) =e^{\bar{m}+\frac{\sigma ^{2}}{2}}.$ The problem is autonomous so we write the current value Bellman equation with time independent value function $V$% \begin{equation} V\left( A_{t}\right) =\max_{c_{t}}\left\{ \ln c_{t}+\beta E_{t}V\left( \left( A_{t}-c_{t}\right) \left( 1+\tilde{m}_{t+1}\right) \right) \right\} . \end{equation} The necessary first order condition for $c_{t}$ yield \begin{equation} \frac{1}{c_{t}}=\beta E_{t}V^{^{\prime }}\left( A_{t+1}\right) \left( 1+% \tilde{m}_{t+1}\right) . \label{FOCStoch} \end{equation} Now we use Merton's result and guess that the value function is \begin{equation} V\left( A_{t}\right) =Y+X\ln A_{t}, \end{equation} for some constants $Y$ and $X$. Substituting into (\ref{FOCStoch}), we get \begin{eqnarray} \frac{1}{c_{t}} &=&\beta E_{t}\frac{X}{A_{t+1}}\left( 1+\tilde{m}% _{t+1}\right) , \\ &=&\beta E_{t}\frac{X\left( 1+\tilde{m}_{t+1}\right) }{\left( A_{t}-c_{t}\right) \left( 1+\tilde{m}_{t+1}\right) }, \\ &=&\beta \frac{X}{A_{t}-c_{t}}, \\ &\rightarrow &c_{t}=\frac{A_{t}}{1+\beta X}, \\ A_{t}-c_{t} &=&\frac{\beta XA_{t}}{1+\beta X}. \end{eqnarray} Now we have to solve for the constant $X$. This is done by substituting the solutions to the first order conditions and the guess into the Bellman equation, \begin{gather} Y+X\ln A_{t}=\max_{c_{t}}\left\{ \ln c_{t}+\beta E_{t}V\left( \left( A_{t}-c_{t}\right) \left( 1+\tilde{m}_{t+1}\right) \right) \right\} \\ =\ln \frac{A_{t}}{1+\beta X}+\beta E_{t}V\left( \frac{\beta XA_{t}}{1+\beta X% }\left( 1+\tilde{m}_{t+1}\right) \right) \\ =\ln \frac{A_{t}}{1+\beta X}+\beta E_{t}\left( Y+X\ln \left( \frac{\beta XA_{t}}{1+\beta X}\left( 1+\tilde{m}_{t+1}\right) \right) \right) \\ =\left( 1+\beta X\right) \ln A_{t}-\left( 1+\beta X\right) \ln \left( 1+\beta X\right) +\beta Y \\ +\beta X\ln \beta X+\beta Xm. \end{gather} This is satisfied for all $A_{t}$ iff \begin{eqnarray} X &=&\left( 1+\beta X\right) \\ &\rightarrow &X=\frac{1}{1-\beta }, \\ Y &=&-\frac{1}{1-\beta }\ln \left( \frac{1}{1-\beta }\right) +\beta Y+\beta \frac{1}{1-\beta }\ln \beta \frac{1}{1-\beta }+\beta \frac{m}{1-\beta } \\ &\rightarrow &Y=-\frac{1}{1-\beta }\ln \frac{1}{1-\beta }+\frac{\beta }{% \left( 1-\beta \right) ^{2}}\left( m+\ln \beta \right) . \end{eqnarray} Thus, \begin{eqnarray} c_{t} &=&\frac{A_{t}}{1+\beta \frac{1}{1-\beta }}=\left( 1-\beta \right) A_{t}, \\ A_{t}-c_{t} &=&\beta A_{t}. \\ A_{t+1} &=&A_{t}\beta \left( 1+\tilde{m}_{t+1}\right) , \\ \ln A_{t+1} &=&\ln A_{t}+\ln \beta +\ln \left( 1+\tilde{m}_{t+1}\right) \\ E_{t}\ln A_{t+1} &=&\ln A_{t}+\ln \beta +m. \end{eqnarray}% Note that since $\ln \left( 1+\tilde{m}_{t+1}\right)$ is normally distributed, $\left( 1+\tilde{m}_{t+1}\right) >0,$ for all $t,$ implying $% A_{t}>0$ for all $t.$ If, on the other hand, $\left( 1+\tilde{m}% _{t+1}\right)$ can be negative with positive probability, $E_{t}\ln \left( 1+\tilde{m}_{t+1}\right)$ is minus infinity implying that the value function is ill-defined. \subsection{Contraction mappings} In the previous section we discussed guessing on solutions to the Bellman equation. However, we would like to know whether there exists a solution and whether it is unique. If the latter is not the case, it is not in principle sufficient to guess and verify, since we might have other value functions that also satisfy the Bellman equation. To prove existence and uniqueness we will apply a contraction mapping argument.\footnote{% An alternative is sometimes to look for the limit $\lim_{s\rightarrow \infty }T^{s}\left( V(k,0)\right)$, which typically is the solution we are interested in (at least in macroeconomics).} For this purpose, we first have to define some concepts. Before, let's look at a simple example of a contaction mapping. Below, I show a graph of Sweden. Some of you might recognize the dots as representing some cities and towns in our country. \FRAME{itbpFU}{3.1609in}{4.0049in}{0in}{\Qcb{Map of Sweden}}{}{slide13.emf}{% \special{language "Scientific Word";type "GRAPHIC";maintain-aspect-ratio TRUE;display "USEDEF";valid_file "F";width 3.1609in;height 4.0049in;depth 0in;original-width 10.0016in;original-height 7.5135in;cropleft "0.2265";croptop "0.9495";cropright "0.6979";cropbottom "0.1510";filename 'Figures/Slide13.EMF';file-properties "XNPEU";}} Next, let's put the same map in a smaller scale on top of the first. This is done in the next picture. \FRAME{itbpF}{3.1462in}{4.0058in}{0in}{}{}{slide14.emf}{\special{language "Scientific Word";type "GRAPHIC";maintain-aspect-ratio TRUE;display "USEDEF";valid_file "F";width 3.1462in;height 4.0058in;depth 0in;original-width 10.0016in;original-height 7.5135in;cropleft "0.2270";croptop "0.9496";cropright "0.6972";cropbottom "0.1511";filename 'Figures/Slide14.EMF';file-properties "XNPEU";}} If you look closely, you will see that two dots (cities) overlap. In fact, they represent the nice city Karlstad in V\"{a}rmland. Perhaps more importantly (at least for some), you cannot find any other dots overlapping. In fact, you might be able to convince yourself that no other points in Sweden than Karlstad overlap. If you would move around the map but keeping the smaller map inside the borders of the larger, you will be able to verify for yourself that there is always one and one only points in Sweden that overlap in the two maps. One can actually prove this result -- and generalizing this proof gives us the contraction mapping proof. Here is a sketch of how the result is formalized. Consider the space of points in the large map. Put the smaller scale map ontop and this of this as representing a mapping from the all points in large map back to the same set of points. Karlstad maps back to Karlstad, Malm\"{o} maps to a point a little west of J\"{o}nk\"{o}ping, Kiruna to a point northwest of Sundsvall and so on. Lets call this mapping $T$ and the set of points in the large map $X.$ Clearly, for all $x\in X,$ $T\left( x\right) \in X.$ Furthermore, \begin{equation*} Karlstad=T\left( Karlstad\right) \end{equation*}% which is the unique fix-point of $T.$ Further more, we note that the distance between any two points $x,y$ in the domain (large map), is strictly larger than the distance between $T\left( x\right)$ and $T\left( y\right) ,$ i.e., between the corresponding points in the small map. Now, the contraction mapping theorem says that if; \begin{enumerate} \item For all $x\in X,T\left( x\right) \in X,$i.e., it maps elements in $X$ back into $X$\thinspace\ (Which in this case means that the small scale map is completely inside the larger map), and \item there is a distance measure $d,$such that for all $x,y$ $\in X,d\left( x,y\right) >d\left( T\left( x\right) ,T\left( y\right) \right)$ (the small scale map has a strictly smaller scale), \end{enumerate} then there is exactly one fix-point $x=T\left( x\right) .$ $T$ is then called a \emph{contraction mapping.} Below, we generalize this result to another space $X$, namely the space of continuous and bounded functions. Clearly, the Bellman $T$-operator maps elements of the set bounded continuos (value)functions, which we call $C,$ into the same set. What remains is to define a distance measure. We use the so called sup-norm. This is defined as follows. Consider two functions (in our case value functions), $w\left( s\right)$ and $v\left( s\right) ,s\in S.$ This two functions maps elements in a state space $S\mathbf{\subseteq R}% ^{n}\rightarrow \mathbf{R.}$ $w$ and $v$ are themselves elements of $C.$ Then, the distance $d\left( w,v\right)$ between the two functions is the maximum absolute difference, \begin{equation*} d\left( w,v\right) \equiv \sup_{s\in \mathbf{S}}\left\Vert w\left( s\right) -v\left( s\right) \right\Vert =\sup_{s\in \mathbf{S}}\sqrt{\left( w\left( s\right) -v\left( s\right) \right) ^{2}} \end{equation*} Typically, our Bellman \thinspace $T$-operator is a contraction mapping if we make sure we have bounded state spaces and there is strict discounting. The next sections make these claims a little more formally. \subsubsection{Complete Metric Spaces and Cauchy Sequences} Let $\mathbf{X}$ be a vector space, i.e., a set on which addition and scalar multiplication is defined. Also define an operator $d$: which we can think of as measuring the (generalized) distance between any two elements of $% \mathbf{X}$. We call $d$ a norm assumed to satisfy \begin{enumerate} \item Positivity $\forall x,y\,\in \mathbf{X}$, $d\left( x,y\right) \geq 0$ and $d\left( x,y\right) =0\Rightarrow x=y.$ \item Symmetry $\forall x,y\,\in \mathbf{X}$, $d\left( x,y\right) =d\left( y,x\right) .$ \item Triangle inequality $\forall x,y,z\,\in \mathbf{X}$, $d\left( x,z\right) \leq d\left( x,y\right) +d\left( y,z\right)$ \end{enumerate} Now, we call $(\mathbf{X},d)$ a \emph{normed vector space} or a \emph{metric space}. An example of such a space would be $\mathbf{R}^{n}$ together with the Euclidean norm $d\left( x,y\right) \equiv \left\Vert x,y\right\Vert$. Another example is the space $\mathbf{C(S)}$ of continuous, bounded functions where each element is a function from $\mathbf{S\subseteq R}% ^{n}\rightarrow \mathbf{R}$ together with the \textquotedblright sup-norm\textquotedblright\ defined as follows. For any two elements in $% \mathbf{C}\left( \mathbf{S}\right)$, i.e., any two functions $w\left( s\right)$ and $v\left( s\right)$, the distance $d$ between them is the maximal euclidean distance, i.e., \begin{equation} d\left( w,v\right) \equiv \sup_{s\in \mathbf{S}}\left\Vert w\left( s\right) ,v\left( s\right) \right\Vert \end{equation} Now let us define a Cauchy sequence. Intuitively, this is a sequence of elements $\{x_{n}\}$ in a space $\mathbf{X}$ that come closer and closer to each other, using some particular norm. More precisely, $\{x_{n}\}$ is defined as a sequence of elements in $\mathbf{X}$ such that for all $% \varepsilon >0$, there exist a number $n$, such that for all $m,p\geq n$, $% d(x_{m},x_{p})<\varepsilon$. An example of such a sequence would be the sequence $\{1,1/2,1/3,\}$ which is a Cauchy sequence using the Euclidean norm. A Cauchy sequence converges if there is an element $y\in \mathbf{X}$ such that $\lim_{n\rightarrow \infty }d(x_{n},y)=0$. It may, of course, be the case that the Cauchy sequence does not converge to a point in $\mathbf{X}$. An example would be if we let $\mathbf{X}$ be the open interval $(0,1]$ and look at the Cauchy sequence $\{1,1/2,1/3,\}$ which is in $\mathbf{X}$ but converges to zero which is not in $\mathbf{X.}$ \subsubsection{Complete metric spaces} Now we are ready to define the \emph{complete metric space}. This is a metric space in which all Cauchy sequences converge to a point in the space. \subsubsection{ Contraction Mapping} Consider a metric space $(\mathbf{X},d)$ and look at an operator $T$ that maps $\mathbf{X}$ $\rightarrow$ $\mathbf{X}$. $T$ is a contraction mapping by definition if there exists a non-negative number $\rho \in \lbrack 0,1),$ such that for all elements $x,y\in \mathbf{X,}$ \begin{equation} d\left( T\left( x\right) ,T\left( y\right) \right) \leq \rho d\left( x,y\right) , \label{contractiondef} \end{equation} where we note that $\rho$ must be \emph{strictly} smaller than one. \subsubsection{The Contraction Mapping Theorem} Now we can state the very important \emph{contraction mapping theorem.} \begin{result} Consider a complete metric space $\left( \mathbf{X},d\right)$ and let $T:% \mathbf{X\rightarrow X}$ be a contraction mapping. Then, $T$ has one unique fixed point $x^{\ast }\in \mathbf{X}$, i.e., the solution to $x=T(x)$ always exists and is unique. Furthermore, the sequence $x_{0}$, $T(x_{0})$,$% T^{2}(x_{0})$,...,$T^{n}(x_{0})$ converges to $x^{\ast }$ for all $x_{0}\in \mathbf{X}$. \end{result} There are theorems that can be used to show that $T$ is a contraction mapping. \begin{result} \label{Blackwell}Let the state space $\mathbf{S}$ be a subset of $\mathbf{R}% ^{n}$ and $\mathbf{B}(\mathbf{S})$ the set of all continuous, \emph{bounded} functions from $\mathbf{S}$ to $\mathbf{R}$. Endowed with the sup-norm, $% \mathbf{B}(\mathbf{S})$, is a complete metric space. Let $T$ be a map that maps all elements of $\mathbf{B}(\mathbf{S)\rightarrow }$ $\mathbf{B}(% \mathbf{S)}$. Then, $T$ is a contraction mapping if \begin{enumerate} \item for any functions $w(s),v(s)\in \mathbf{B}(\mathbf{S})$, the following holds; if $w(s)\geq v(s)\forall s\in \mathbf{S}$ then $T\left( w(s)\right) \geq T\left( v(s)\right) \forall s\in \mathbf{S}$ (\emph{monotonicity}), and \item there is a $\beta \in \lbrack 0,1)$ such that for any constant $\kappa \in R$, and any function $w\left( s\right) \in \mathbf{B}\left( \mathbf{S}% \right) ,T\left( w(s)+\kappa \right) =T\left( w(s)\right) +\beta \kappa$. (% \emph{discounting}). \end{enumerate} \end{result} Usually it is straightforward to apply the previous result to show that if we have strict discounting, the Bellman equation is a contraction mapping. There is one major limitation which we have to live with, however, result % \ref{Blackwell} and variants of it require bounded value functions. Let us look at an example, where we apply result \ref{Blackwell}. Consider a simple growth model, \begin{eqnarray} &&\max_{\{c_{t}\}_{0}^{\infty }}E_{0}\sum_{t=0}^{\infty }\beta ^{t}u\left( c_{t}\right) \\ \text{s.t. }k_{t+1} &=&f\left( k_{t}\right) -c_{t},\forall t\geq 0,k_{0}\,% \text{given, and} \\ c_{t},k_{t} &>&0\text{ }\forall t\text{.} \end{eqnarray} where $u$ is a continuous and increasing (utility) function with $u\left( 0\right) \geq u_{\min }>-\infty$ and $0\leq \beta <1.$ To use the theorems we need to make some assumptions. First, we need boundedness. For this purpose, we assume \begin{eqnarray} f\left( 0\right) &=&0, \\ f^{\prime }\left( k\right) &\geq &0\forall k\geq 0, \\ \exists \bar{k} &>&0,\text{such that }f\left( k\right) \leq k\forall k\geq \bar{k}. \end{eqnarray} Now, define $\mathbf{S}=[0,\bar{k}]$ and note that if $k_{0}$ is in $\mathbf{% S}$, so is all admissible $k_{t}$. Then, $u\left( c_{t}\right)$ is bounded since $u_{\min }\leq u\left( c\right) \leq u\left( \bar{k}\right)$ implying that any value function must satisfy $\frac{u_{\min }}{1-\beta }\leq V(k)\leq \frac{u\left( \bar{k}\right) }{1-\beta }.$ By restricting the state space $\mathbf{S}=[0,\bar{k}]$, we can therefore restrict our search for value functions that are bounded on our state space. Now, let us establish that the following mapping is a contraction: \begin{equation} V\left( k\right) =\max_{c\geq 0}u\left( c\right) +\beta V\left( f\left( k)-c\right) \right) \equiv T\left( V\left( k\right) \right) \end{equation} Regarding condition 1, we need for any two bounded functions $v\left( k\right) ,w\left( k\right) ,k\in S,.$ \begin{equation} v\left( k\right) \geq w\left( k\right) \forall k\Rightarrow T\left( v\left( k\right) \right) \geq T\left( w\left( k\right) \right) \forall k. \end{equation}% which is satisfied. To see this, define \begin{equation} c^{\ast }=\arg \max_{c}u\left( c\right) +\beta w\left( f\left( k\right) -c\right) \end{equation}% then, \begin{eqnarray} T\left( v\left( k\right) \right) &=&\max_{c}u\left( c\right) +\beta v\left( f\left( k\right) -c\right) \geq u\left( c^{\ast }\right) +\beta v\left( f\left( k\right) -c^{\ast }\right) \\ &\geq &u\left( c^{\ast }\right) +\beta w\left( f\left( k\right) -c^{\ast }\right) =T\left( w\left( k\right) \right) . \end{eqnarray} Regarding condition 2, we note \begin{eqnarray} T\left( v\left( k\right) +\kappa \right) &=&\max_{c\geq 0}u\left( c\right) +\beta \left( v\left( f\left( k-c\right) \right) +\kappa \right) \\ &=&\max_{c\geq 0}u\left( c\right) +\beta v\left( f\left( k-c\right) \right) +\beta \kappa \\ &=&T\left( v\left( k\right) \right) +\beta \kappa , \end{eqnarray}% so the second condition is satisfied too. Thus, the Bellman equation is a contraction mapping and always has one and one only solution, $V\left( k\right) .$ \newpage \section{Dynamic Optimization in Continuous Time} \subsection{Dynamic programming in continuous time} %TCIMACRO{\TeXButton{TeX field}{\setcounter{equation}{0}}}% %BeginExpansion \setcounter{equation}{0}% %EndExpansion Consider the problem \begin{eqnarray} &&\max_{k\left( t\right) _{0}^{T},c\left( t\right) _{0}^{T}}\int_{0}^{T}e^{-rt}u\left( k,c,t\right) dt \label{Conttime1} \\ s.t.\,\dot{k} &=&f\left( k,c,t\right) \label{LoM1} \\ k\left( 0\right) &=&\underline{k}, \end{eqnarray} with \begin{gather} k\left( T\right) =\bar{k}\text{ (case 1), or} \label{endcases} \\ \text{ }k\left( T\right) \text{ free (case 2), or} \\ k\left( T\right) \geq \bar{k}\text{ (case 3).} \end{gather} Thinking of the integral in the maximand as a sum of rectangles with base $% dt$ and height $e^{-rt}u\left( k,c,t\right)$, we can approximate the problem with a discrete time problem. Noting that for a small time interval $% dt$, $k\left( t+dt\right) =k\left( t\right) +f\left( k,c,t\right) dt,$ we can write the current value Bellman equation \begin{equation} V\left( k,t\right) =\max_{c}\left\{ u\left( k,c,t\right) dt+e^{-rdt}V\left( k+f\left( k,c,t\right) dt,t+dt\right) \right\} . \end{equation} We can then make Taylor approximations; \begin{gather} V\left( k,t\right) =\max_{c}{\Large \{}u\left( k,c,t\right) dt \\ +\left( 1-rdt\right) \left( V\left( k,t\right) +V_{k}\left( k,t\right) f\left( k,c,t\right) dt+V_{t}\left( k,t\right) dt\right) {\Large \}} \end{gather} Subtracting $\left( 1-rdt\right) V\left( k,t\right)$ from both sides, dividing by $dt$ and then letting $dt$ go to zero yields, \begin{eqnarray} rV\left( k,t\right) &=&\max_{c}\left\{ u\left( k,c,t\right) +V_{k}\left( k,t\right) f\left( k,c,t\right) +V_{t}\left( k,t\right) \right\} \label{arbitrage} \\ &=&\max_{c}\left\{ u\left( k,c,t\right) +\frac{dV\left( k,t\right) }{dt}% \right\} , \label{arbitrageb} \end{eqnarray}% where we note that $dV\left( k,t\right) /dt$ is the total time derivative of the value function. Sometimes, this equation is referred to as an asset pricing equation -- and interpreted as follows; if an asset is correctly valued (not providing arbitrage opportunities) the opportunity cost of holding it (the LHS of (\ref{arbitrage})) equals the (optimal) sum of the immediate pay-off or dividend (the first term of (\ref{arbitrageb})) and the capital gain (the second term of (\ref{arbitrageb})). Sometimes we can use the guess and verify technique to solve for the value function if the problem is autonomous. An alternative is to use \emph{% Optimal Control} and \emph{Pontryagin's maximum principle.} \subsection{Optimal Control} Consider the problem in (\ref{Conttime1}) with the continuous time Bellman equation (\ref{arbitrage}). Noting that since the term $V_{t}\left( k,t\right)$ is independent of $c,$ we can rewrite (\ref{arbitrage}) \begin{equation} rV\left( k,t\right) -V_{t}\left( k,t\right) =\max_{c}\left\{ u\left( k,c,t\right) +V_{k}\left( k,t\right) f\left( k,c,t\right) \right\} \label{Arbitrage2} \end{equation} Defining the \emph{co-state} or current shadow value variable as the derivative of the value function w.r.t. $k,$ along its optimal path $k^{\ast }$, \begin{equation} \lambda \left( t\right) \equiv V_{k}\left( k^{\ast },t\right) . \end{equation} We see that a necessary condition for $c^{\ast }\left( t\right)$ to be optimal is that it is given by \begin{equation} c^{\ast }\left( t\right) =\arg \max_{c}\left\{ u\left( k,c,t\right) +\lambda \left( t\right) f\left( k,c,t\right) \right\} . \label{Pont1b} \end{equation} Now, we need to pin down $\lambda \left( t\right) .$ For this purpose we analyze how $\lambda$ develops over time by deriving a differential equation for $\lambda \left( t\right) .$Taking derivatives w.r.t. $k$ of the identity (\ref{Arbitrage2}), we get \begin{gather} rV_{k}\left( k,t\right) -V_{kt}\left( k,t\right) =u_{k}\left( k,c^{\ast },t\right) +V_{kk}\left( k,t\right) f\left( k,c^{\ast },t\right) \\ +V_{k}\left( k,t\right) f_{k}\left( k,c^{\ast },t\right) , \\ r\lambda \left( t\right) -\left( V_{tk}\left( k,t\right) +V_{kk}\left( k,t\right) \dot{k}\right) =u_{k}\left( k,c^{\ast },t\right) +\lambda \left( t\right) f_{k}\left( k,c^{\ast },t\right) , \\ r\lambda \left( t\right) -\dot{\lambda}\left( t\right) =u_{k}\left( k,c^{\ast },t\right) +\lambda \left( t\right) f_{k}\left( k,c^{\ast },t\right) , \label{Pont2} \end{gather}% where in second equation, we changed the order of diffentiation, using the fact that $\frac{\partial }{\partial k}\frac{\partial V\left( k,t\right) }{% \partial t}=\frac{\partial }{\partial t}\frac{\partial V\left( k,t\right) }{% \partial k}.$ Equation (\ref{Pont1b}) and (\ref{Pont2}) form the basis of Pontryagin's maximum principle. According to this, we derive necessary (and sometimes sufficient) conditions for an optimal control. We do this without explicitly solving for the value function. We first define the \emph{current value Hamiltoninan}. This is the sum of the instantaneous payoff and the product of costate(s) and the function determining the law-of-motion of the state variable; \begin{equation} \mathcal{H}\left( k,c,\lambda ,t\right) \mathcal{\equiv }u\left( k,c,t\right) +\lambda \left( t\right) f\left( k,c,t\right) . \end{equation} Note that the Hamiltonian has the same interpretation as the RHS\ of the Bellman equation without the max-operator. In words, it is the sum of the flow of current pay-off and the generation of future pay-off. In the Bellman equation, we use the value function to measure future payoffs while the Hamiltonian uses the shadow-value $\lambda .$ According to Pontryagin's maximum principle, the optimal control $c^{\ast }\left( t\right)$ maximizes the Hamiltonian at each instant, the co-state (or shadow value) satisfies the differential equation, $r\lambda \left( t\right) -\dot{\lambda}\left( t\right) =\mathcal{H}_{k}\left( k,c,\lambda ,t\right)$. These necessary conditions will provide differential equations which we need to solve. Typically we have one initial condition for each state variable. But we need more information to solve the system since we also have the control variable(s). In case (1) of (\ref{endcases}), we have the necessary additional information. In case (2), it must be that the shadow value of the state variable approach zero as $t\rightarrow T.$ This is called, \emph{the transversality condition. }In case (3), either the inequality is slack, in which case $\lambda \left( T\right) =0,$ or it binds, giving the necessary additional info in both cases. We can the summarize: \begin{result} According to the \emph{Pontryagin's maximum principle}\newline 1. An optimal control $c^{\ast }\left( t\right) ,$ satisfies \begin{equation} c^{\ast }\left( t\right) =\arg \max_{c}\mathcal{H}\left( k,c,\lambda ,t\right) , \label{pontr1} \end{equation}% \newline 2. where the current co-state or shadow value $\lambda \left( t\right)$ , is a continuous function of time, satisfying the differential equation. \begin{equation} r\lambda \left( t\right) -\dot{\lambda}\left( t\right) =\mathcal{H}% _{k}\left( k,c,\lambda ,t\right) , \label{pont2} \end{equation}% except at points in time where $c$ is discontinuous and with \newline 3. end-condition(s) provided by \begin{align} k\left( T\right) & =\bar{k}\text{ (case 1) or} \label{pont3a} \\ \lambda \left( T\right) & =0\text{ if }k\left( T\right) \text{ is free (case 2), or} \label{pont3b} \\ \lambda \left( T\right) & \geq 0,\text{ and }\lambda \left( T\right) \left( k\left( T\right) -\bar{k}\right) =0\text{ (case (3).} \label{pont3c} \end{align}% \newline In addition, the path for $k\left( t\right)$ must satisfy the initial condition, $k\left( 0\right) =\underline{k}$, and $\dot{k}\left( t\right) =% \mathcal{H}_{\lambda }\left( k,c,\lambda ,t\right) =f\left( k,c,t\right) .$ \end{result} \subsubsection{The consumption problem} As an example, consider the problem standard consumption-savings (Ramsey) problem \begin{eqnarray} &&\max_{k\left( t\right) _{0}^{T},c\left( t\right) _{0}^{T}}\int_{0}^{T}e^{-rt}u\left( c\right) dt \label{Ramseyob} \\ s.t.\,\dot{k} &=&f\left( k\right) -c \\ k\left( 0\right) &=&\underline{k}, \\ k\left( T\right) &\geq &0, \end{eqnarray}% with $u$ and $f$ increasing and concave. The current value Hamiltonian is \begin{equation} \mathcal{H}\left( k,c,\lambda ,t\right) =u\left( c\right) +\lambda \left( t\right) \left( f\left( k\right) -c\right) . \end{equation} Thus, $c^{\ast }\left( t\right) =\arg \max_{c}u\left( c\right) +\lambda \left( t\right) \left( f\left( k\right) -c\right)$ which is interior, implying \begin{gather} \mathcal{H}_{c}\left( k,c,\lambda ,t\right) =0. \\ \rightarrow u^{\prime }\left( c^{\ast }\right) =\lambda . \end{gather} Furthermore, from the second condition \begin{equation} \mathcal{H}_{k}\left( k,c,\lambda ,t\right) =\lambda f^{\prime }\left( k\right) =r\lambda -\dot{\lambda}. \label{lam} \end{equation} Taking time-derivatives of $u^{\prime }\left( c^{\ast }\right) =\lambda$ and substituting into (\ref{lam}) we get \begin{eqnarray} u^{\prime \prime }\left( c^{\ast }\right) \dot{c}^{\ast } &=&\dot{\lambda}, \\ u^{\prime }\left( c^{\ast }\right) f^{\prime }\left( k\right) &=&ru^{\prime }\left( c^{\ast }\right) -u^{\prime \prime }\left( c^{\ast }\right) \dot{c}% ^{\ast }, \\ &\rightarrow &\dot{c}^{\ast }=\frac{u^{\prime }\left( c^{\ast }\right) }{% -u^{\prime \prime }\left( c^{\ast }\right) }\left( f^{\prime }\left( k\right) -r\right) , \end{eqnarray}% which is the Euler equation we have seen before. To analyze the behavior of this system, we can use the phase-diagram as in section \ref{phasediagram}. To get a closed form solution, i.e., an expression for the endogenous variables in terms of only the exogenous ones, we must specify the utility and production functions. Considering first the utility functions, we have two important special cases. First, CARA utility, \begin{eqnarray} u\left( c\right) &=&-\frac{e^{-\gamma c}}{\gamma }, \\ &\rightarrow &u^{\prime }\left( c\right) =e^{-\gamma c},u^{\prime \prime }\left( c\right) =-\gamma e^{-\gamma c} \end{eqnarray}% in which case we get \begin{equation} \dot{c}^{\ast }=\frac{1}{\gamma }\left( f^{\prime }\left( k\right) -r\right) , \end{equation}% i.e., consumption increase is a linear function in the difference between the marginal return on savings and the subjective discount rate. The other case is CRRA, \begin{eqnarray} u\left( c\right) &=&\frac{c^{1-\sigma }}{^{1-\sigma }}, \\ &\rightarrow &u^{\prime }\left( c\right) =c^{-\sigma },u^{\prime \prime }\left( c\right) =-\sigma c^{-\sigma -1}, \end{eqnarray}% yielding \begin{eqnarray} \dot{c}^{\ast } &=&\frac{c^{-\sigma }}{\sigma c^{-\sigma -1}}\left( f^{\prime }\left( k\right) -r\right) =\frac{c}{\sigma }\left( f^{\prime }\left( k\right) -r\right) \label{Ramseyc} \\ \frac{\dot{c}^{\ast }}{c} &=&\frac{1}{\sigma }\left( f^{\prime }\left( k\right) -r\right) , \label{Ramseyk} \end{eqnarray}% i.e., the consumption growth \emph{rate }is a linear function in $f^{\prime }\left( k\right) -r.$ The sensitivity is given by $1/\sigma$, which we call the intertemporal elasticity of substitution. Here, we also have that $% \sigma$ is the constant of relative risk aversion. What about the transversality condition? In this case, we know $u^{\prime }\left( c^{\ast }\left( T\right) \right) =\lambda \left( T\right) .$ So since utility in the examples is unbounded, i.e., $u^{\prime }\left( c\right) >0$ for all finite $c,$ $\lambda \left( T\right)$ cannot be 0, instead $k\left( T\right)$ is zero. In other words, whenever consumption is valuable at $T$, the lower bound on $k$ should bind and nothing should be left. Let us complete the example by assuming, for simplicity, a linear (Romer type) production function $f\left( k\right) =Ak$. In the CRRA\ case, we get the linear system \begin{eqnarray} \dot{c}^{\ast } &=&\frac{A-r}{\sigma }c \\ \dot{k} &=&-c+Ak \\ \left[ \begin{array}{c} \dot{c}^{\ast } \\ \dot{k}% \end{array}% \right] &=&\left[ \begin{array}{cc} \frac{A-r}{\sigma } & 0 \\ -1 & A% \end{array}% \right] \left[ \begin{array}{c} c^{\ast } \\ k% \end{array}% \right] . \end{eqnarray} This system has roots $A$ and $\frac{A-r}{\sigma }$ and the matrix of eigenvectors is \begin{equation} \left[ \begin{array}{cc} 0 & \frac{r+A\left( \sigma -1\right) }{\sigma } \\ 1 & 1% \end{array} \right] \equiv \mathbf{B}^{-1}. \end{equation} Consequently, the solution is of the form \begin{eqnarray} \left[ \begin{array}{c} c^{\ast }\left( t\right) \\ k\left( t\right)% \end{array}% \right] &=&\mathbf{B}^{-1}\left[ \begin{array}{cc} e^{At} & 0 \\ 0 & e^{\frac{A-r}{\sigma }t}% \end{array}% \right] \left[ \begin{array}{c} \kappa _{1} \\ \kappa _{2}% \end{array}% \right] \\ &=&\left[ \begin{array}{cc} 0 & \frac{r+A\left( \sigma -1\right) }{\sigma }e^{\frac{A-r}{\sigma }t} \\ e^{At} & e^{\frac{A-r}{\sigma }t}% \end{array}% \right] \left[ \begin{array}{c} \kappa _{1} \\ \kappa _{2}% \end{array}% \right] , \end{eqnarray}% where $\kappa _{1}$ and $\kappa _{2}$ are two integration constants. We solve for the latter by using $k\left( 0\right) =\underline{k}$ and $k\left( T\right) =0.$% \begin{eqnarray} k\left( 0\right) &=&\underline{k}=\left[ \begin{array}{cc} 1 & 1% \end{array}% \right] \left[ \begin{array}{c} \kappa _{1} \\ \kappa _{2}% \end{array}% \right] \\ &=&\kappa _{1}+\kappa _{2} \\ k\left( T\right) &=&0=\left[ \begin{array}{cc} e^{T} & e^{\frac{A-r}{\sigma }T}% \end{array}% \right] \left[ \begin{array}{c} \kappa _{1} \\ \kappa _{2}% \end{array}% \right] \\ &\rightarrow &\kappa _{1}=\frac{e^{\frac{A-r-\sigma }{\sigma }T}}{e^{\frac{% A-r-\sigma }{\sigma }T}-1}\underline{k},\kappa _{2}=\frac{1}{1-e^{\frac{% A-r-\sigma }{\sigma }T}}\underline{k}. \end{eqnarray} We can now, for example, evaluate \begin{eqnarray} c^{\ast }\left( t\right) &=&\frac{r+A\left( \sigma -1\right) }{\sigma }e^{% \frac{A-r}{\sigma }t}\frac{1}{1-e^{\frac{A-r-\sigma }{\sigma }T}}\underline{k% }, \\ c^{\ast }\left( 0\right) &=&\frac{r+A\left( \sigma -1\right) }{\sigma }\frac{% 1}{1-e^{\frac{A-r-\sigma }{\sigma }T}}\underline{k}. \end{eqnarray}% $\allowbreak$ \subsection{Sufficiency} Assume that $f$ and $u$ are concave in $k,c$ and $\lambda \geq 0$. This implies that the Hamiltonian is concave in $k,c$. Then, Pontryagin's necessary conditions (\ref{pontr1}) and (\ref{pont2}) and (\ref{pont3a}), or (\ref{pont3b}) or (\ref{pont3c}) are sufficient. \subsection{Infinite horizon} Consider the infinite horizon problem \begin{eqnarray} &&\max_{c\left( t\right) _{0}^{\infty }}\int_{0}^{\infty }e^{-rt}u\left( k,c,t\right) dt \\ s.t.\,\dot{k} &=&f\left( k,c,t\right) \\ k\left( 0\right) &=&\underline{k}, \end{eqnarray} Pontryagin's conditions (\ref{pontr1}) and (\ref{pont2}) are necessary also in the infinite horizon case, provided, of course, that there is a well defined solution. If there is a binding restriction on the state variable of the type $\lim_{T\rightarrow \infty }k\left( T\right) =\bar{k}$, this can help us pin down the solution. The finite horizon transversality conditions can, however, not immediately be used in the infinite horizon case. Suppose the \emph{maximized} Hamiltonian is concave in $k$ for every $t$, then the conditions (\ref{pontr1}) and (\ref{pont2}) plus the infinite horizon transversality conditions \begin{eqnarray} \lim_{T\rightarrow \infty }e^{-rT}\lambda \left( T\right) k\left( T\right) &=&0,\text{ and } \label{trans1} \\ \lim_{T\rightarrow \infty }e^{-rT}\lambda \left( T\right) &\geq &0, \label{trans2} \end{eqnarray}% provide a \emph{sufficient} set of conditions for optimality (see de la Fuente, 1999, p 577). Often, the Hamiltonian is concave in $k,c$ together. This is sufficient for the maximized Hamiltonian to be concave in $k.$ Sometimes a so called \emph{No-Ponzi} condition helps us to make sure that the transversality conditions are satisfied. Suppose, for example, the pay-off $u\left( k,c,t\right) =u\left( c\right)$, that $k$ represents debt of the agent and for simplicity that $f\left( k,c,t\right) =c+\rho k-w$, so debt increases by the difference between consumption plus interest payments $% \rho k$ and the wage $w$. It is reasonable to assume that creditors demand to be repaid in a present value sense -- the discounted value of future repayment should always be at least as large as debt. This is the No-Ponzi condition. When in addition, the agent prefers to pay back no more than he ows, the implication is. \begin{equation} \lim_{T\rightarrow \infty }e^{-\rho T}k\left( T\right) =0. \end{equation}% To see this, solve \begin{equation} \dot{k}\left( t\right) -\rho k\left( t\right) =c\left( t\right) -w\left( t\right) \end{equation}% giving \begin{eqnarray} e^{-\rho t}\left( \dot{k}\left( t\right) -\rho k\left( t\right) \right) &=&e^{-\rho t}\left( c\left( t\right) -w\left( t\right) \right) \\ e^{-\rho t}k\left( t\right) &=&\int_{0}^{t}\left( e^{-\rho s}\left( c\left( t\right) -w\left( t\right) \right) \right) ds+k\left( 0\right) \\ \lim_{T\rightarrow \infty }e^{-\rho T}k\left( T\right) &=&\int_{0}^{\infty }\left( e^{-\rho s}\left( c\left( t\right) -w\left( t\right) \right) \right) ds+k\left( 0\right) \end{eqnarray}% So, the No-Ponzi requirement is that if the PDV of "mortgage" repayments is no smaller than initial debt, i.e., $-\int_{0}^{\infty }\left( e^{-\rho s}\left( c\left( t\right) -w\left( t\right) \right) \right) ds$ $k\left( 0\right) \geq k\left( 0\right) ,$ then $\lim_{T\rightarrow \infty }e^{-\rho T}k\left( T\right) \leq 0.$ Clearly, when marginal utility is strictly positive, the individual would never want to satisfy this with inequality, since he could then increase consumption. Therefore, $\lim_{T\rightarrow \infty }e^{-\rho T}k\left( T\right) =0.$ The second necessary condition in (\ref{pont2}) is now \begin{eqnarray} -\left( \rho -r\right) \lambda \left( t\right) &=&\dot{\lambda}\left( t\right) \\ \lambda \left( t\right) &=&\lambda \left( 0\right) e^{-\left( \rho -r\right) t}. \end{eqnarray}% So provided marginal utility is positive at $t=0$, (\ref{trans2}) is satisfied. Furthermore, \begin{eqnarray} \lim_{T\rightarrow \infty }e^{-rT}\lambda \left( T\right) k\left( T\right) &=&\lambda \left( 0\right) \lim_{T\rightarrow \infty }e^{-rT}e^{-\left( \rho -r\right) T}k\left( T\right) \\ &=&\lambda \left( 0\right) \lim_{T\rightarrow \infty }e^{-\rho T}k\left( T\right) \\ &=&0. \end{eqnarray}% where the last equality is the No-Ponzi condition. Sometimes, the sufficient conditions allow us to identify the optimal control as the stable manifold (saddle-path) leading to a saddle-point stable steady state. Consider again the problem \begin{eqnarray} &&\max_{c\left( t\right) _{0}^{\infty }}\int_{0}^{\infty }e^{-rt}u\left( c\right) dt \\ \text{s.t.}\,\dot{k} &=&f\left( k\right) -c,\text{ }k\left( 0\right) =% \underline{k}, \end{eqnarray}% which we graphically analyzed in section \ref{phasediagram} showing the existence of saddle-path and a steady state with \begin{eqnarray} f^{\prime }\left( k^{ss}\right) &=&r \\ c^{ss} &=&f\left( k^{ss}\right) \end{eqnarray} Restating the current value Hamiltonian \begin{equation} \mathcal{H}\left( k,c,\lambda ,t\right) =u\left( c\right) +\lambda \left( t\right) \left( f\left( k\right) -c\right) , \end{equation} we note that if both $u\left( c\right)$ and $f\left( k\right)$ are concave, and $\lambda \left( t\right) \geq 0,$ $\mathcal{H}\left( k,c,\lambda ,t\right)$ is concave in $k,c$ so the conditions for using the sufficiency result are satisfied. In addition to (\ref{Ramseyc}), (\ref% {Ramseyk}), and $k\left( 0\right) =\underline{k}$, we thus only need to verify that (\ref{trans1}) and (\ref{trans2}) are satisfied. This is straightforward, \begin{eqnarray} \lim_{T\rightarrow \infty }e^{-rT}\lambda \left( T\right) &=&\lim_{T\rightarrow \infty }e^{-rT}u^{\prime }\left( c^{ss}\right) =0, \\ \lim_{T\rightarrow \infty }e^{-rT}\lambda \left( T\right) k\left( T\right) &=&\lim_{T\rightarrow \infty }e^{-rT}u^{\prime }\left( c^{ss}\right) k^{ss}=0. \end{eqnarray} \subsection{Present value Hamiltonian} Sometimes, it is convenient to define the present value Hamiltonian, i.e., expressing everything in values as seen from time 0. In problem (\ref% {Conttime1}), the present value Hamiltonian is given by \begin{equation} \mathbf{H}\left( k,c,\mu ,t\right) =e^{-rt}u\left( k,c,t\right) +\mu f\left( k,c,t\right) , \end{equation} where $\mu \left( t\right)$ is the present shadow value of the state variable. In this case, the necessary conditions for optimality are \begin{eqnarray} c^{\ast }\left( t\right) &=&\arg \max_{c}\mathbf{H}\left( k,c,\mu ,t\right) \\ -\dot{\mu}\left( t\right) &=&\mathbf{H}_{k}\left( k,c,\mu ,t\right) \\ \dot{k} &=&\mathbf{H}_{\mu }\left( k,c,\mu ,t\right) \end{eqnarray} In the finite horizon case, the transversality conditions are the same in terms of $\mu \left( T\right)$ and $\lambda \left( T\right) .$In the infinite horizon case, we note that \begin{equation} \mu \left( T\right) =e^{-rT}\lambda \left( T\right) , \end{equation}% so the conditions (\ref{trans1}) and (\ref{trans2}) become \begin{eqnarray} \lim_{T\rightarrow \infty }\mu \left( T\right) k\left( T\right) &=&0,\text{ and } \\ \lim_{T\rightarrow \infty }\mu \left( T\right) &\geq &0. \end{eqnarray} Note also that Michel (Econometrica, July 1982) states \begin{equation*} \lim_{\lim_{T\rightarrow \infty }}\mathbf{H}\left( k,c,\mu ,t\right) =0 \end{equation*}% as a necessary general transversality condition in the infinite horizon. This condition seems intuitive in the sense that it can be interpreted as saying that all profit opportunities should be exhausted in the long run. However, the issue of general transversaility conditions is still in dispute, and caution is warranted although counterexamples typically seem quite "cooked". \subsection{Bang-Bang solutions} We have already noted that the condition $r\lambda \left( t\right) -\dot{% \lambda}\left( t\right) =\mathcal{H}_{k}\left( k,c,\lambda ,t\right)$ and $-% \dot{\mu}\left( t\right) =\mathbf{H}_{k}\left( k,c,\mu ,t\right)$ holds everywhere except at points in time where $c$ is discontinuous. Let us take an example when such discontinuities arise because the Hamiltonian is linear in the control. The control variable has to be chosen from a compact set and the linearity implies that the optimal control jumps from one end of the constraint to the other -- a \emph{bang-bang} solution. An individual allocates his daily time between working and enjoying leisure. When enjoying leisure, his utility increases in the amount of toys he has to play with. If he works, he gets money which he spends on increasing his stock of toys. On the other hand, while working, he cannot enjoy his toys. The problem can be written, \begin{eqnarray*} &&\max_{\{l\}_{0}^{T}}\int\limits_{0}^{T}\left( 1-l\left( t\right) \right) U\left( x\left( t\right) \right) dt \\ s.t.\text{ }\dot{x}(t) &=&l(t), \\ 1-l(t) &\geq &0, \\ l(t) &\geq &0, \\ x(0) &=&0, \end{eqnarray*}% where $l\in \lbrack 0,1]$ is the share of time per day spent on working, $x$ is the stock of toys and $U$ is a strictly increasing utility function with $% U(0)=0$. The Hamiltonian for the problem is% \begin{equation*} H=\left( 1-l\left( t\right) \right) U\left( x\left( t\right) \right) +\lambda \left( t\right) \left( l(t)\right) . \end{equation*} Since there is no terminal condition on $x\left( T\right) ,$ the transversality condition is $\lambda \left( T\right) =0.$The necessary conditions for optimality are% \begin{eqnarray} l &=&\arg \max_{l}\left( 1-l\left( t\right) \right) U\left( x\left( t\right) \right) +\lambda \left( t\right) \left( l(t)\right) , \\ H_{x} &=&\left( 1-l\left( t\right) \right) U_{x}\left( x\left( t\right) \right) =-\dot{\lambda}\left( t\right) . \label{eq_Bang2} \end{eqnarray} Clearly, $\lambda >U\left( x\right)$ implies $l^{\ast }=1$ and vice versa. From (\ref{eq_Bang2}) we see that $\lambda$ is weakly falling. Therefore, $% \lambda \left( 0\right) >0$ (Otherwise, the total payoff would be zero which cannot be optimal). Since $U\left( x\left( 0\right) \right) =0,$ $l\left( 0\right) =1.$ Furthermore, it cannot be optimal to work all time. We conclude that the optimal policy is a stopping rule -- work, i.e., set $l\left( t\right) =1$ until $t^{\ast }$ then set $l\left( t\right) =0$ for $t\in \left[ t^{\ast },T% \right] .$ From (\ref{eq_Bang2}) we see that $\dot{\lambda}\left( t\right) =0$ for $t\in \lbrack 0,t^{\ast })$ and we also have $\dot{x}(t)=1$ for $% t\in \lbrack 0,t^{\ast }).$ For $t\in (t^{\ast },T],-\dot{\lambda}\left( t\right) =-\frac{1}{2}x\left( t^{\ast }\right) ^{-1/2}$ with end condition $% \lambda \left( T\right) =0.$ This in combination with the fact that $\lambda$ is continuos and $\lambda \left( t^{\ast }\right) =U\left( t^{\ast }\right)$ is enough to solve the problem. If we specify for example $U(x)=\sqrt{x},$ we get a closed form solution. Summarizing \begin{eqnarray*} x\left( t\right) &=&\left\{ \begin{array}{c} t\text{ for }tt^{\ast }% \end{array}% \right. \\ x\left( 0\right) &=&0 \\ \dot{\lambda}\left( t\right) &=&\left\{ \begin{array}{c} 0\text{ for }tt^{\ast }% \end{array}% \right. \\ \lambda \left( t\right) &=&\left\{ \begin{array}{c} \lambda \left( t^{\ast }\right) \text{ for }tt^{\ast }% \end{array}% \right. \\ \lambda \left( t^{\ast }\right) &=&U\left( x\left( t^{\ast }\right) \right) =\left( t^{\ast }\right) ^{\frac{1}{2}} \\ \lambda \left( T\right) &=&0. \end{eqnarray*} Using $\lambda \left( t^{\ast }\right) =U\left( x\left( t^{\ast }\right) \right) =\left( t^{\ast }\right) ^{\frac{1}{2}}$ in $\lambda \left( t^{\ast }\right) -\frac{1}{2}x\left( t^{\ast }\right) ^{-1/2}\left( t-t^{\ast }\right)$ and solving for $t^{\ast }$ gives, \begin{equation*} t^{\ast }=\frac{1}{3}T. \end{equation*} \subsection{Many state variables and controls} Having several state variables and controls pose no principle problem. Neither does pointwise discontinuities in the control variable. To generalize, suppose we have $n$ state variables and $n$ controls. An optimal control maximizes the Hamiltonian over all available controls $\mathbf{c}$% \begin{equation} \mathbf{c}^{\ast }\left( t\right) =\arg \max_{\mathbf{c}}\mathcal{H}\left( \mathbf{k,c},\mathbf{\lambda },t\right) \equiv u\left( \mathbf{k,c},t\right) +\sum_{i=1}^{n}\lambda _{i}f_{i}\left( \mathbf{k,c},t\right) . \end{equation}% where $\lambda _{i}$ is the shadow value associated with the state variable $% k_{i}$. Each $\lambda _{i}\left( t\right)$ is continuous and satisfies the differential equation \begin{equation} r\lambda _{i}\left( t\right) -\dot{\lambda}_{i}\left( t\right) =\frac{% \partial }{\partial k_{i}}\mathcal{H}\left( \mathbf{k,c},\mathbf{\lambda }% ,t\right) , \end{equation}% except when $c$ is discontinuous. For the transversality conditions, we have \begin{align} k_{i}\left( T\right) & =\bar{k}_{i}\text{ case (1), or} \\ \lambda _{i}\left( T\right) & =0\text{ case (2), or} \\ \lambda _{i}\left( T\right) & \geq 0,\text{ and }\lambda _{i}\left( T\right) \left( k_{i}\left( T\right) -\bar{k}_{i}\right) =0\text{ case (3),} \end{align}% for end-conditions for state variable $i$ belonging to case 1, 2 or 3. \newpage \section{Some numerical methods} %TCIMACRO{\TeXButton{TeX field}{\setcounter{equation}{0}}}% %BeginExpansion \setcounter{equation}{0}% %EndExpansion \subsection{Numerical solution to Bellman equations} When we cannot solve the Bellman equation analytically, there are several methods to approximate a solution numerically. One of the most straightforward methods when the problem is autonomous, is to discretize the state space and iterate on the Bellman equation until it converges. When the Bellman equation is a contraction mapping, strong results make sure that this procedure converges to the correct value function. When we discretize the state space, we restrict the state variable to take values from a finite set \begin{equation} k_{t}\in \left[ k^{1},k^{2},...k^{n}\right] \equiv K \end{equation}% where the superscripts index the elements of $K.$ We then solve for $c$ from the law-of-motion for $k$% \begin{eqnarray} k_{t+1} &=&f\left( k_{t},c_{t}\right) \\ &\implies &c_{t}=\bar{f}\left( k_{t},k_{t+1}\right) \end{eqnarray} We can then write the Bellman equation for the discretized problem as% \begin{equation} V\left( k_{t}\right) =\max_{k_{t+1}\in K}u\left( k_{t},\bar{f}\left( k_{t},k_{t+1}\right) \right) +\beta V\left( k_{t+1}\right) . \end{equation} As you see, this is a Bellman for a constrained problem, i.e., the control variable is constrained relative to the case when $k_{t+1}$ is continuous. Two things should be noted; First, the Bellman equation is the true Bellman equation of the constrained problem and previous results hold, in particular the contraction mapping theorems apply. Second, how important the constraint implied by the discretization depends on how fine the grid is. Many (few) elements of $K$ with small (large) distances between them, imply that the constraint is weak (severe). Denoting an arbitrary initial value function by $V_{0}\left( k_{t}\right) ,$ being $n$ numbers, we update this value function according to% \begin{equation} V_{1}\left( k_{t}\right) =\max_{k_{t+1}\in K}u\left( k_{t},\bar{f}\left( k_{t},k_{t+1}\right) \right) +\beta V_{0}\left( k_{t+1}\right) \end{equation}% giving $V_{1}\left( k_{t}\right) ,$ being a new set of $n$ numbers. We then iterate on the Bellman equation% \begin{equation} V_{s+1}\left( k_{t}\right) =\max_{k_{t+1}\in K}u\left( k_{t},\bar{f}\left( k_{t},k_{t+1}\right) \right) +\beta V_{s}\left( k_{t+1}\right) \end{equation}% until $V_{s+1}\left( k_{t}\right) \approx V_{s}\left( k_{t}\right) .$ For each of the $n$ values of $k_{t}$, we check $u\left( k_{t},\bar{f}\left( k_{t},k_{t+1}\right) \right) +\beta V\left( k_{t+1}\right)$ for all $n$ values of $k_{t+1}$ and choose the $k_{t+1}$ that gives the highest value, giving $V_{s+1}\left( k_{t}\right)$. Therefore, each iteration requires $% n^{2}$ evaluations when the state variable is unidimensional.\footnote{% When the state variable is higher of dimensionality, this method quickly becomes computationally too burdensome.} Lets consider a simple example. \begin{eqnarray} &&\max_{\left\{ k_{t+1},c_{t}\right\} _{t}^{\infty }}\sum_{t=0}^{\infty }\beta ^{t}\ln \left( c_{t}\right) \\ \text{s.t. }k_{t+1} &=&f\left( k_{t},c_{t}\right) \equiv k_{t}^{\alpha }+\left( 1-\delta \right) k_{t}-c_{t} \\ k_{t} &\geq &0\forall t \\ k_{0} &=&\underline{k} \end{eqnarray} First, we solve for \begin{equation} c_{t}=\bar{f}\left( k_{t},k_{t+1}\right) =k_{t}^{\alpha }+\left( 1-\delta \right) k_{t}-k_{t+1}. \end{equation} Then, we note that $k_{t}\in \left[ 0,\delta ^{\frac{1}{\alpha -1}}\right] \implies$ $k_{t+1}\in \left[ 0,\delta ^{\frac{1}{\alpha -1}}\right]$, implying that value function is bounded. If also $\beta <1,$ the Bellman equation% \begin{equation} V\left( k_{t}\right) =\max_{c_{t}}\ln \left( k_{t}^{\alpha }+\left( 1-\delta \right) k_{t}-k_{t+1}\right) +\beta V\left( k_{t+1}\right) \end{equation}% is a contraction mapping. Let us parametrize, setting $\beta =0.9$, $\delta =.2$ and $\alpha =1/2\implies \delta ^{\frac{1}{\alpha -1}}=25$ and discretize the state space by requiring \begin{equation} k_{t}\in \left[ 5,10,15,20,25\right] \equiv K\text{~}\forall t. \end{equation}% Now set an initial value function, for example \begin{equation*} V_{0}\left( k\right) =\ln k~\forall k. \end{equation*} This is then updated in the following way. For each possible $k_{t},k_{t+1}$ we calculate the left hand side of the Bellman equation, and solve the maximization problem.\ So, for $k_{t}=5,$ and all $k_{t+1}\in K$ we have% \begin{eqnarray*} \ln \left( 5^{\alpha }+\left( 1-\delta \right) 5-5\right) +.9\ln 5 &=&1.66 \\ \ln \left( 5^{\alpha }+\left( 1-\delta \right) 5-10\right) +.9\ln 10 &=&-\infty \\ \ln \left( 5^{\alpha }+\left( 1-\delta \right) 5-15\right) +.9\ln 15 &=&-\infty \\ \ln \left( 5^{\alpha }+\left( 1-\delta \right) 5-20\right) +.9\ln 20 &=&-\infty \\ \ln \left( 5^{\alpha }+\left( 1-\delta \right) 5-25\right) +.9\ln 25 &=&-\infty \end{eqnarray*} Implying that the updated value function for $k_{t}=5$ is \begin{equation*} V_{1}\left( 5\right) =1.66. \end{equation*}% For $k_{t}=10,$% \begin{eqnarray*} \ln \left( 10^{\alpha }+\left( 1-\delta \right) 10-5\right) +.9\ln 5 &=&3.27 \\ \ln \left( 10^{\alpha }+\left( 1-\delta \right) 10-10\right) +.9\ln 10 &=&2.22 \\ \ln \left( 10^{\alpha }+\left( 1-\delta \right) 10-15\right) +.9\ln 15 &=&-\infty \\ \ln \left( 10^{\alpha }+\left( 1-\delta \right) 10-20\right) +.9\ln 20 &=&-\infty \\ \ln \left( 10^{\alpha }+\left( 1-\delta \right) 10-25\right) +.9\ln 25 &=&-\infty \end{eqnarray*}% implying \begin{equation*} V_{1}\left( 10\right) =3.27 \end{equation*} In the same way, for $k_{t}=15$ \begin{eqnarray*} \ln \left( 15^{\alpha }+\left( 1-\delta \right) 15-5\right) +.9\ln 5 &=&3.83 \\ \ln \left( 15^{\alpha }+\left( 1-\delta \right) 15-10\right) +.9\ln 10 &=&3.84 \\ \ln \left( 15^{\alpha }+\left( 1-\delta \right) 15-15\right) +.9\ln 15 &=&2.30 \\ \ln \left( 15^{\alpha }+\left( 1-\delta \right) 15-20\right) +.9\ln 20 &=&-\infty \\ \ln \left( 15^{\alpha }+\left( 1-\delta \right) 15-5\right) +.9\ln 25 &=&-\infty \end{eqnarray*} \begin{equation*} V_{1}\left( 15\right) =3.84 \end{equation*} Doing this also for $k_{t}=20$ and $k_{t}=25$ completes the first iteration. Then, we repeat the iterations until we think the process has converged sufficiently, \subsection{Band Matrix Methods for differential equations} Assume we want to solve the differential equation \begin{equation} y^{\prime \prime }\left( t\right) +ay^{\prime }\left( t\right) +by\left( t\right) =g\left( t\right) \end{equation}% over some interval, with initial conditions given. If we want to solve this numerically, we first have to get rid of the abstract infinitely small differences called differentials. We approximate these with finite size forward differences such that \begin{eqnarray} y^{\prime }\left( t\right) &\approx &\frac{y\left( t+\Delta t\right) -y\left( t\right) }{\Delta t}, \\ y^{\prime }\left( t\right) &\approx &\frac{\frac{y\left( t+\Delta t\right) -y\left( t\right) }{\Delta t}-\frac{y\left( t\right) -y\left( t-\Delta t\right) }{\Delta t}}{\Delta t} \\ &=&\frac{y\left( t+\Delta t\right) -2y\left( t\right) +y\left( t-\Delta t\right) }{\Delta t^{2}} \end{eqnarray} Using this we can solve the equation for a finite set of values in the following way. Say we want to solve the equation in the interval $t\in \lbrack p,q]$ and we know $y^{\prime \prime }\left( 0\right) =c_{1}$ and $% y^{\prime }\left( 0\right) =c_{2}.$ We divide the interval for $t$ into $n$ equal parts and use the following notation \begin{eqnarray} t_{k} &=&p+\frac{k\left( q-p\right) }{n}, \\ t_{k}-t_{k-1} &=&\frac{q-p}{n}\equiv \Delta t, \\ y\left( t_{k}\right) &\equiv &y_{k.} \end{eqnarray} This gives us the following equations \begin{eqnarray} \frac{\left( y_{-1}-2y_{0}+y_{1}\right) }{\Delta t^{2}} &=&c_{1} \\ \frac{\left( -y_{0}+y_{1}\right) }{\Delta t} &=&c_{2} \\ \frac{\left( y_{-1}-2y_{0}+y_{1}\right) }{\Delta t^{2}}+a\frac{-y_{0}+y_{1}}{% \Delta t}+by_{0} &=&g\left( t_{0}\right) \\ \frac{\left( y_{0}-2y_{1}+y_{2}\right) }{\Delta t^{2}}+a\frac{-y_{1}+y_{2}}{% \Delta t}+by_{1} &=&g\left( t_{1}\right) \\ &&. \\ &&. \\ \frac{\left( y_{n-1}-2y_{n}+y_{n+1}\right) }{\Delta t^{2}}+a\frac{% -y_{n}+y_{n+1}}{\Delta t}+by_{n} &=&g\left( t_{n}\right) . \end{eqnarray} This provides $n+3$ linear equations for the $n+3$ unknown $y$. Writing this as a system we have \begin{eqnarray} \mathbf{A}\left[ \begin{array}{c} y_{-1} \\ y_{0} \\ y_{1} \\ . \\ y_{n} \\ y_{n+1}% \end{array}% \right] &=&\left[ \begin{array}{c} c_{1} \\ c_{2} \\ g\left( t_{0}\right) \\ g\left( t_{1}\right) \\ . \\ g\left( t_{n}\right)% \end{array}% \right] \\ \left[ \begin{array}{c} y_{-1} \\ y_{0} \\ y_{1} \\ . \\ y_{n} \\ y_{n+1}% \end{array}% \right] &=&\mathbf{A}^{-1}\left[ \begin{array}{c} c_{1} \\ c_{2} \\ g\left( t_{0}\right) \\ g\left( t_{1}\right) \\ . \\ g\left( t_{n}\right)% \end{array}% \right] , \end{eqnarray}% with (setting $n=3$) \begin{gather} \mathbf{A\equiv } \\ _{\left[ \begin{array}{cccccc} _{\Delta t^{-2}} & _{-2\Delta t^{-2}} & _{\Delta t^{-2}} & _{0} & _{0} & _{0} \\ _{0} & _{-\Delta t^{-1}} & _{\Delta t^{-1}} & _{0} & _{0} & _{0} \\ _{\Delta t^{-2}} & _{-2\Delta t^{-2}-a\Delta t^{-1}+b} & _{\Delta t^{-2}+a\Delta t^{-1}} & _{0} & _{0} & _{0} \\ _{0} & _{\Delta t^{-2}} & _{-2\Delta t^{-2}-a\Delta t^{-1}+b} & _{\Delta t^{-2}+a\Delta t^{-1}} & _{0} & _{0} \\ _{0} & _{0} & _{\Delta t^{-2}} & _{-2\Delta t^{-2}-a\Delta t^{-1}+b} & _{\Delta t^{-2}+a\Delta t^{-1}} & _{0} \\ _{0} & _{0} & _{.} & _{\Delta t^{-2}} & _{-2\Delta t^{-2}-a\Delta t^{-1}+b} & _{\Delta t^{-2}+a\Delta t^{-1}}% \end{array}% \right] } \notag \end{gather} To get any accuracy, we should of course set $n$ much larger than $3.$ As we see, the matrix $\mathbf{A}$ contains many zeros, with a \emph{band }of non-zeros around the diagonal. Due to this feature, it is easy for the computer to invert it also if $n$ is in the order of hundreds. \subsection{Newton-Raphson} Suppose we are looking for an optimum of the real valued function $f\left( \mathbf{x}\right) ,\mathbf{x\in R}^{n}$ where $f$ is twice differentiable. A standard way to do this numerically is to apply the Newton-Raphson algorithm. If the optimum is interior, it satisfies the necessary first order conditions that the gradient is zero, i.e., at the optimum, denoted $% \mathbf{x}^{\ast }$, \begin{equation} Df\left( \mathbf{x}\right) =0. \end{equation} Now apply a first order linear approximation to the gradient from some initial point $\mathbf{x}_{0}$ \begin{equation} 0=Df\left( \mathbf{x}^{\ast }\right) \approx Df\left( \mathbf{x}_{0}\right) +D^{2}f\left( \mathbf{x}_{0}\right) \left( \mathbf{x}^{\ast }-\mathbf{x}% _{0}\right) , \end{equation} where $D^{2}f\left( \mathbf{x}_{0}\right)$ is the \emph{Hessian} matrix of second derivatives of $f.$ Provided the Hessian is invertible, we can get an approximation to $\mathbf{x% }^{\ast },$% \begin{equation} \mathbf{x}^{\ast }\approx \mathbf{x}_{0}-\left( D^{2}f\left( \mathbf{x}% _{0}\right) \right) ^{-1}Df\left( \mathbf{x}_{0}\right) . \end{equation} From this we can construct a search algorithm that under some circumstances makes better and better approximations \begin{equation} \mathbf{x}_{s+1}\approx \mathbf{x}_{s}-\left( D^{2}f\left( \mathbf{x}% _{s}\right) \right) ^{-1}Df\left( \mathbf{x}_{s}\right) . \end{equation} If we don't have analytic solutions to the gradient and Hessian we can use numerical approximations, for example the forward difference method; For a small number $\varepsilon$, we have, \begin{eqnarray} \frac{\partial f\left( \mathbf{x}\right) }{\partial x_{1}} &=&\frac{f\left( \mathbf{x}+\left[ \begin{array}{c} \varepsilon \\ \mathbf{0}% \end{array}% \right] \right) -f\left( \mathbf{x}\right) }{\varepsilon } \\ \frac{\partial f^{2}\left( \mathbf{x}\right) }{\partial x_{1}^{2}} &=&\frac{% f\left( \mathbf{x}+\left[ \begin{array}{c} \varepsilon \\ \mathbf{0}% \end{array}% \right] \right) -2f\left( \mathbf{x}\right) +f\left( \mathbf{x-}\left[ \begin{array}{c} \varepsilon \\ \mathbf{0}% \end{array}% \right] \right) }{\varepsilon ^{2}}, \\ \frac{\partial f^{2}\left( \mathbf{x}\right) }{\partial x_{2}\partial x_{1}} &=&\frac{\frac{\partial f\left( \mathbf{x}+\left[ \begin{array}{c} 0 \\ \varepsilon \\ \mathbf{0}% \end{array}% \right] \right) }{\partial x_{1}}-\frac{\partial f\left( \mathbf{x}\right) }{% \partial x_{1}}}{\varepsilon } \\ &=&\frac{f\left( \mathbf{x}+\left[ \begin{array}{c} \varepsilon \\ \varepsilon \\ \mathbf{0}% \end{array}% \right] \right) -f\left( \left[ \begin{array}{c} 0 \\ \varepsilon \\ \mathbf{0}% \end{array}% \right] \right) -f\left( \mathbf{x}+\left[ \begin{array}{c} \varepsilon \\ 0 \\ \mathbf{0}% \end{array}% \right] \right) +f\left( \mathbf{x}\right) }{\varepsilon ^{2}}. \end{eqnarray} One should be very careful with this method since it can only find local optima in the case when $f$ is not globally concave. In well-behaved problems, it is however, easily programmed and fairly quick. \newpage \section{Note on constant hazard models} %TCIMACRO{\TeXButton{TeX field}{\setcounter{equation}{0}} }% %BeginExpansion \setcounter{equation}{0} %EndExpansion Consider as an example a pool of unemployed, with measure (size) at time $t$ given by $x(t).$ Suppose also that there is an outflow of people from the unemployment pool (we are for now disregarding the inflow of new unemployed). The outflow is determined by an assumption that there is a constant probability per unit of time, denoted $h$ to get hired. Using this we can derive a law-of-motion for $x\left( t\right)$. Over a small (infinitesimal) interval of time $dt$, we have \begin{eqnarray} x\left( t+dt\right) &=&x\left( t\right) \left( 1-hdt\right) \\ x\left( t+dt\right) -x\left( t\right) &=&-x\left( t\right) hdt \\ \lim_{dt\rightarrow 0}\frac{x\left( t+dt\right) -x\left( t\right) }{dt} &\equiv &\dot{x}\left( t\right) =-hx\left( t\right) . \end{eqnarray} This is a simple differential equation with a solution \begin{equation} x\left( t\right) =x\left( 0\right) e^{-ht}. \label{poolsize} \end{equation} Now, consider an individual who is unemployed at time $0$ and the random variable $s$ denote the time she will stay unemployed. Let us now derive the probability density function $f\left( s\right)$ and let $F\left( s\right)$ denote the cumulative distribution function, i.e., $F\left( s\right)$ is the probability the unemployment spell is no longer than $s.$ Clearly, \begin{equation} F\left( s\right) =\int_{0}^{s}f\left( t\right) dt. \end{equation} From (\ref{poolsize}) we know that at time $s$, a share \begin{equation} \frac{x\left( t\right) }{x\left( 0\right) }=e^{-hs}, \end{equation} remains unemployed and since hiring is completely random, \begin{equation} 1-\frac{x\left( t\right) }{x\left( 0\right) }=1-e^{-hs}=F\left( s\right) , \end{equation} and consequently \begin{equation} f\left( s\right) =he^{-hs}. \end{equation} Let us also verify that the probability of finding a job per unit of time, conditional on not having found it stays constant. This probability is \begin{equation} \frac{f\left( t\right) }{1-F\left( t\right) }=\frac{he^{-ht}}{e^{-ht}}=h. \end{equation} We can now compute the average spell length as \begin{equation} \int_{0}^{\infty }sf\left( s\right) ds=\int_{0}^{\infty }she^{-hs}ds. \end{equation} Using the formula for integration by parts, \begin{eqnarray} \int_{0}^{\infty }she^{-hs}ds &=&\left[ -se^{-hs}\right] _{0}^{\infty }-\int_{0}^{\infty }-e^{-hs}ds \\ &=&0-0-\left[ \frac{e^{-hs}}{h}\right] _{0}^{\infty }=\frac{1}{h}. \end{eqnarray} Similarly, we can compute the median length $m$, i.e., solving $F(m)=1/2$ from \begin{eqnarray} 1/2 &=&1-e^{-hm} \\ &\rightarrow &e^{-hm}=1/2 \\ -hm &=&-\ln 2 \\ m &=&\frac{\ln 2}{h}\approx \frac{.69}{h} \end{eqnarray} This is sometimes called the \emph{rule of 69; }expressing $h$ in percent per unit of time. The half-life is found by dividing $69$ by $h.$ For example, if the probability of finding job is 5\% per week. It takes $% 69/5\approx 14$ weeks before half the pool of unemployed have found jobs and the average unemployment spell is $1/0.05=20$ weeks. Another example; it is often found that differences in GDP between countries (after controlling for differences is savings and schooling) is closing by $% 3\%$ per year. Then, half the difference is left after $69/3=23$. The rule of 69 also works when something is growing. If a bank account yields 4\% return per year, it takes $69/4\approx 17$ years for it to double. \section{T-mappings\label{Tmap}} Instead of using time as an argument of the value function, let's use time subscripts. We can then write the Bellman equation as \begin{equation*} V_{t}\left( k_{t}\right) =\left\{ \begin{array}{c} \max_{c_{t}}\left( u\left( k_{t},c_{t}\right) +\beta V_{t+1}\left( k_{t+1}\right) \right) \\ s.t.k_{t+1}=f\left( k_{t},c_{t}\right)% \end{array}% \right. \end{equation*} Now, let us define an operator that maps next period value functions$% ,V_{t+1}:K\rightarrow R,$ (where $K$ is the state space, i.e., the set of possible values for the state variable), into functions that provides the current value associated with all $k_{t}$ in the state space. Thus, the operator, which we will call,\thinspace $T,$ maps elements of the space of value functions, call that space $C,$ back into the same space, i.e., $% T:C\rightarrow C$.\footnote{% Later on, we will make assumptions such that we can further restrict the space of possible valute functions.} Formally, we define the $T$ mapping as \begin{equation*} TV_{t+1}:C\rightarrow C\equiv \left\{ \begin{array}{c} \max_{c_{t}}\left( u\left( k_{t},c_{t}\right) +\beta V_{t+1}\left( k_{t+1}\right) \right) \\ s.t.k_{t+1}=f\left( k_{t},c_{t}\right)% \end{array}% \right. . \end{equation*}% When we want to indicate that the mapped function, $TV_{t+1},$ is a function of, e.g., $k_{t},$we append $\left( k_{t}\right) ;$% \begin{equation*} TV_{t+1}=TV_{t+1}\left( k_{t}\right) . \end{equation*} Note that while $V_{t+1}$ is a function of $k_{t+1}$, $TV_{t+1}$ is a function of $k_{t}.$ Let us take an example. Suppose $u\left( k_{t},c_{t}\right) =\ln c_{t}$ and $% f\left( k_{t},c_{t}\right) =k_{t}-c_{t}$ $\forall t.$ Let's us now see what the $T$ operator does. Take a particular element in $C,$ for example $\ln k:R^{+}\rightarrow R.$ So here the state space $K=R^{+}.$ Now, \begin{equation*} T\ln k_{t}=\left\{ \begin{array}{c} \max_{c_{t}}\left( \ln c_{t}+\beta \left( \ln k_{t+1}\right) \right) \\ s.t.k_{t+1}=k_{t}-c_{t}% \end{array}% \right. \end{equation*} The first order condition is \begin{equation*} \frac{1}{c_{t}}=\frac{\beta }{k_{t}-c_{t}}\rightarrow c_{t}=\frac{k_{t}}{% 1+\beta } \end{equation*}% thus, \begin{eqnarray*} T\ln k_{t} &=&\ln \frac{k_{t}}{1+\beta }+\beta \left( \ln \left( k_{t}-\frac{% k_{t}}{1+\beta }\right) \right) \\ &=&\left( \left( 1+\beta \right) \ln k_{t}+\beta \ln \beta -\left( 1+\beta \right) \ln \left( 1+\beta \right) \right) \end{eqnarray*} which is another function $R^{+}\rightarrow R.$ So, we see that the $T$ maps functions into (potentially) other functions, concluding the example. Using the $T$ operator, the Bellman equation in period $t$ can be written \begin{equation*} V_{t}=TV_{t+1}, \end{equation*}% or, equivalently, \begin{equation*} V_{t}\left( k_{t}\right) =TV_{t+1}\left( k_{t}\right) . \end{equation*} Furthemore, in the next period, the next period's Bellman equation is \begin{equation*} V_{t+1}=TV_{t+2}. \end{equation*} Thus, \begin{equation*} V_{t}=TV_{t+1}=T^{2}V_{t+2}. \end{equation*} The meaning of $T^{2}V_{t+2}$ in words is; give me (I am $T^{2}$) a value function that applies in period $t+2$ (you don't need to say anything about what $k_{t+2}$ is going to be), and I give you a function that tells you the value in period $t$ associated with all values of $k_{t}$ in the state space. Formally;% \begin{eqnarray*} T^{2}V_{t+2} &=&\left\{ \begin{array}{c} \max_{c_{t}}\left( u\left( k_{t},c_{t}\right) +\beta \left( TV_{t+2}\left( k_{t+1}\right) )\right) \right) \\ s.t.k_{t+1}=f\left( k_{t},c_{t}\right)% \end{array}% \right. \\ &=&\left\{ \begin{array}{c} \max_{c_{t}}\left( u\left( k_{t},c_{t}\right) +\beta \left( \max_{c_{t+1}}\left( u\left( k_{t+1},c_{t+1}\right) +\beta V_{t+2}\left( k_{t+2}\right) \right) \right) \right) \\ s.t.k_{t+1}=f\left( k_{t},c_{t}\right) ,k_{t+2}=f\left( k_{t+1},c_{t+1}\right)% \end{array}% \right. . \end{eqnarray*}% Now, suppose in the time autonomuous case that we can find a limiting function \begin{equation*} V\left( k_{t}\right) \equiv \lim_{s\rightarrow \infty }T^{s}V_{t+s}\left( k_{t}\right) \end{equation*}% then, as shown in the lecture notes, this function satisfies the Bellman equation \begin{equation*} V\left( k_{t}\right) =TV\left( k_{t}\right) \end{equation*}% i.e., it is a fixed point of the $T$ operator. This means that in the the space $C$, $V$ is an elements such that $T$ maps back onto the same element. If $T$ is a contraction mapping on $C$, this element exists and is unique. \end{document}