Solving Euler Equations : Classical Methods and the C 1 Contraction Mapping Method Revisited ∗

In this paper we describe the classical methods used to solve the Euler equations. Special attention is devoted to the iterative method based on a contraction mapping derived from these equations in Maldonado and Moreira (2003). We test the numerical robustness of this method when it is used in models with sensitiveness to initial conditions. Finally we extend the method to the case of stochastic Euler equations.


INTRODUCTION
One of the main difficulties of numerical methods solving intertemporal economic models is to find accurate estimates for stationary solutions.For example, in dynamic programming problems, the Bellman equation approach provides a contraction mapping with the value function as its fixed Wilfredo L. Maldonado, Humberto Moreira point.With a value function approximation given by that method, policy function approximations are calculated in such a way that they can be used to describe the optimal path.The precision of the estimation is guaranteed by the asymptotic behavior of the sequence of policy function approximations (see Stokey and Lucas (1989) for details).
Numerical methods based on the Bellman iterations are provided by Taylor and Uhlig (1990) and the estimation error for the approximate policy can be found in Santos and Vigo-Aguiar (1998) and Maldonado and Svaiter (2001).However, Bellman's method can only be applied when the economic model admits a representative agent.Moreover, the speed of convergence of this method is slow.
When the model has no representative agent, the usual way to state the optimality and equilibrium conditions is through the Euler equations.Of course, these equations can only be used if the solution is interior and all the functions involved are differentiable.Numerical methods implemented from these equations have been more efficient than those based on the Bellman equations at least in two aspects: they can be used even when there is no representative agent and they are empirically faster when an adequate approximation scheme is performed (see Judd (1992)).
In this paper we describe the classical methods to estimate the stationary solution of the Euler equations associated to an intertemporal economic model.Among these methods we will describe those in Baxter et al. (1990), Coleman (1990) and Coleman (1991), essentially based on the value function approximation; the projection method of Judd (1992); and the parameterizing-expectations method of Marcet (1988) and Marcet and Lorenzoni (1998).There is also a method described by Li (1998) to solve a specific monetary model.We also revisite the method provided by Maldonado and Moreira (2003) 1 to find approximations to stationary solutions of the Euler equations.This method is defined from a contraction mapping in a functional space which has the stationary solution as its fixed point.The major contribution of this technique is its capability to approximate stationary solutions in the uniform topology of continuously differentiable functions.
In this work we also test the numerical robustness of Maldonado and Moreira's method by applying it to Euler equations with the logistic map as its stationary solution.The test shows that even though the stationary solution has sensitiveness to initial conditions, the contractive method also works efficiently.Finally we extend the C 1 contraction mapping method of Maldonado and Moreira to the case of stochastic Euler equations.This is an important issue since most of the applied models in economic dynamics deal with stochastic shocks in the fundamentals of the model.We illustrate the method by applying it to the neoclassical stochastic growth model.
As highlighted by Maldonado and Moreira (2003) their contractive method is not equivalent to Bellman's one and that it does not need a discretization of the state variable domain (i.e., we can solve it point pointwisely).Moreover, the sufficient condition for the convergence of that method is slightly stronger than the determinacy of the steady state.
The paper is divided as follows.Section 2 presents the deterministic Euler equations and the classical methods to estimate the stationary solution.Section 3 extends the contraction method given by Maldonado and Moreira (2003) to the case of stochastic Euler equations.In section 4 we provide examples that illustrate the method, its robustness and its extension to the stochastic case.In section 5 concluding remarks are given.In the Appendix we present the proof of the stochastic dynamic programming case.
Solving Euler Equations: Classical Methods and the C 1 Contraction Mapping Method Revisited

The Euler equations and assumptions
In intertemporal economic models the equilibrium paths are usually defined by a set of equations that embody optimality and market clearing conditions.These equations, in their simplest form, depend on the current and lagged values and expectations of the state variables of the model.Roughly speaking, let X ⊂ R n be the state space, D ⊂ X × X × X represent the intertemporal feasibility set and E : D → R n be the function whose zeros define the optimality and market clearing conditions.The model is defined by the triple (X, D, E).
An equilibrium path from x 0 ∈ X is a sequence Equations ( 1) are called the Euler equations of the model.In many cases it is possible to construct an equilibrium path from a single function.A stationary solution for (X, D, E) is a function g : X → X such that for all x ∈ X: E(x, g(x), g (2) (x)) = 0. (2) We say that x ∈ R n is a steady state if g(x) = x.Let |.| be one of the equivalent norms of R n .This norm induces a norm in the space of real square matrices of order n (M n ) denoted by . . 2 Finally, let B r (x) = {y ∈ R n ; |y − x| < r}.
We will assume that E is a twice continuously differentiable and that E 2 is negative definite nmatrix on the interior of D.3 We will make the following: Assumption D. There exists an interior steady state x such that: where all the derivatives are evaluated at (x, x, x).
Remark 2.1.When n = 1, it is not difficult to show that (i) implies (ii).
Remark 2.2.The condition for local uniqueness of stationary solution is that the steady state must be a saddle point of the linearization of E = 0. Namely, the characteristic polynomial of the linearization must have n roots with modulus lower than one and n roots with modulus greater than one.A necessary and sufficient condition for this is: which is weaker than assumption D.
As an example of this structure, let us consider the deterministic dynamic programming problem.We have the following elements: the state variable set X ⊂ R n , the correspondence of feasible pairs Γ : X → X, the one-period return function F : Graph(Γ) → R and the discount factor β ∈ (0, 1).So, the dynamic programming problem is to find a path (x t ) t≥0 such that: Wilfredo L. Maldonado, Humberto Moreira In this case, X is the state space, D = {(x, y, z) ∈ X 3 ; (x, y) ∈ Graph(Γ) and (y, z) ∈ Graph(Γ)} and if (x t ) t≥0 is an interior optimal path, then For this structure the stationary solution g corresponds to the optimal policy function of the problem.It is easy to verify that for this model the assumption D-(i) is equivalent to the dominant diagonal condition at the steady state.

Some methods to solve the Euler equations
Now we will present some methods used in the literature for solving the Euler equations.

The policy function iteration method
This method consists in defining an application T from the set of feasible maps h : X → X (h(x) ∈ Γ(x), ∀x ∈ X) into itself such that the following holds: (3) The stationary solution g is a fixed point of T and, under some technical conditions, the sequence (T n (h 0 )) n converges to g pointwisely for a convenient choice of h 0 .It is important to note that, when restricted to the dynamic programming problem, the algorithm given in ( 3) is the same as the Bellman iteration method.Indeed, define the following sequence of functions: From the first order condition and the Envelope Theorem: where h n (x) = argmax {y;(x,y)∈A} F (x, y) + βv n−1 (y).Using the definition of T and this last equation, it is easy to see that h n = T n (h 0 ).Hence, Bellman's method implies that h n converges to g.In particular, this method is equivalent to Bellman's one.
To solve numerically for the value function iterations we can use a discretization of the domain and/or linear-quadratic approximations of the approximated value functions (see Christiano (1990) and McGrattan (1990)).
The method defined in (3) is the same as the Bellman iteration method.Indeed, define the following sequence of functions: for all n ≥ 1, where v 0 (x) = F (x, h 0 (x)).From the first order condition and the Envelope Theorem: Using the definition of T and this last equation, it is easy to see that h n = T n (h 0 ).Hence, Bellman's method implies that h n converges to g.In particular, this method is equivalent to Bellman's one.
To solve numerically for the value function iterations we can use a discretization of the domain and/or linear-quadratic approximations of the approximated value functions (see Christiano (1990) and McGrattan (1990)).
The method defined in (3) was used by Baxter et al. (1990) and Coleman (1990) to solve stochastic growth models.Both papers use a discretization of the state space and linear interpolations to find each element of the sequence of approximated solutions.Coleman (1991) proved the monotonicity of the operator defined in (3) and the uniform convergence of (T n (h 0 )) n to g.
Solving Euler Equations: Classical Methods and the C 1 Contraction Mapping Method Revisited

Projection method
The second group of methods to solve equation ( 2) is based on projection methods.They consist in approximating the stationary solution by finite dimensional families of functions.Namely, suppose that the stationary solution g, defined in a compact domain [m, M ], is approximated by a polynomial: and T i 's are the Chebyshev polynomials4 defined on [0, 1] recursively by T n+1 (x) = 2xT n (x) − T n−1 (x), T 0 (x) = 1 and T 1 (x) = x.Then, the vector a = (a 1 , ..., a n ) is calculated in order to minimize the residuals: R(x; a) = E(x, g(x; a), g(g(x; a); a)). (4) There are at least three ways to minimize the errors in ( 4).The first one consists in selecting some points (x 1 , . . ., x n ) in the domain [m, M ] and solving the system of non-linear equations: R(x i ; a) = E(x i , g(x i ; a), g(g(x i ; a); a)) = 0.
The second one is to minimize the L 2 −norm of those residuals, i.e., Minimize Finally, the following system of non-linear equations might be solved m j=0 R(x j ; a)ψ j (x j ) = 0, where (x 1 , ..., x m ) is the m zeroes of ψ m+1 .Examples of this method are given in Judd (1992).
Similar ideas were developed for stochastic models using "parameterized expectations" and finding the parameters that minimize the residuals (see Marcet (1988) and Marcet and Lorenzoni (1998) for details of the method and examples).

Li's contractive operator defined from the Euler equation
In a specific monetary model, Li (1998) defined the following map: The stationary solution of the model satisfies Φ(g) = g.In that model suitable conditions on the set of parameters guarantee that Φ is a contraction mapping, hence (Φ n (h 0 )) n converges to g in the uniform norm.

Maldonado and Moreira's C
1 contractive method defined from the Euler equation Maldonado and Moreira (2003) provided an iterative method based on the Euler equations.The method consists in defining implicitly an operator from the space of feasible stationary functions to itself.Then, they showed that the operator is a C 1 contraction mapping with the stationary solution as its fixed point.Next we summarize the method and the resulting algorithm from it.
where C 2 (B r (x)) is the space of twice continuously differentiable functions from B r (x) into itself.
Define the norm Dh(x) , for all h ∈ C 2 (B r (x).
Let Hr be the closure of H r with respect to this norm.Therefore, ( Hr , • 1 ) is a complete metric space.Observe that, by the definition of the metric (uniform convergence in the first derivative), it is easy to see that Hr is a subset of C 1 (B r (x)).
Theorem 2.2 claims that the map ϕ : Hr → Hr is a η−contraction and consequently has a fixed point in H r .It is easy to see that such a fixed point is a stationary solution of (X, D, E).
An extension of Theorem 2.2 can be obtained if assumption D is satisfied in a neighborhood of the steady state (see details and proofs in Maldonado and Moreira (2003)).
At this point it is important to highlight that although Li's method and Maldonado and Moreira's method are based on contraction mappings, they are structurally different.Li's contraction mapping is defined explicitly by Φ(h)(x) = E(x, h(x), h (2) (x)) + h(x) whereas Maldonado and Moreira's one is defined implicitly by E(x, ϕ(h)(x), h (2) (x)) = 0. Maldonado and Moreira (2003) showed two advantages of their method over the former: the convergence in the C 1 norm and the set of parameters where convergence holds contains the Li's set of parameters.
Using the contraction mapping ϕ we can define the following algorithm to approximate the stationary solution of the Euler equation.
Let h 0 : X → X be a constant function (for instance, h 0 ≡ x, the steady state).Given x ∈ X, solve E(x, y, h 0 (h 0 (x))) = 0 and denote its solution by h 1 (x).To compute h 2 (x) we have to solve So, first we need to find y = h 1 (h 1 (x)) such that: Then, plugging h 1 (h 1 (x)) in the previous equation we find h 2 (x).In general, we have to proceed as follows: first step given h n , using the Euler equation, find the following sequence where x was dropped in the notation.
Solving Euler Equations: Classical Methods and the C 1 Contraction Mapping Method Revisited Observe that it is not necessary to make any partition of the state space in order to calculate the value of the functions evaluated at x in the sequence above.This is a remarkable property of the method since interpolation generates a source of increasing errors in any algorithm.
Remark 2.3.Since this method is based on a contractive mapping, the error is uniformly decreasing in each step (with respect to the C 1 topology).

THE C 1 CONTRACTION MAPPING METHOD FOR THE STOCHASTIC EULER EQUATIONS
In this section we provide an extension of the C 1 contraction method described in section 2 to the stochastic case.To facilitate the exposition we will restrict ourselves to the stochastic dynamic programming problem.
Let (X, X ) be the set of the state variable values and (Z, Z) the set of shock values.Let both sets be measurable spaces.Assume that the shocks are given by a transition function The feasible set is defined by Ω ⊂ X × X × Z and the one-period return function is F : Ω → R n .As in the deterministic case, β ∈ (0, 1) represents the discount factor.Therefore, the stochastic dynamic programming problem is defined by (X, Z, Q, Ω, F, β).
Under standard assumptions, 6 we can obtain the existence of an optimal policy function g : X × Z → X.In addition, under differentiability and interiority of solution hypotheses the optimal policy function must satisfy the stochastic Euler equation: for all k ∈ X and z ∈ Z.
We assume the following assumption:
(ii) There exists The first assumption says that the optimal path is strongly interior.assumption S (ii.a) is the stochastic version of the diagonal dominant condition.Finally, Assumption S (ii.b) is a technical one which is the extension of the assumption D (ii). Let The following lemma defines the map ϕ in analogous way of the deterministic case.
It is worth noting that S is a condition on the optimal policy function of the problem (which is unknown).Therefore, it has to be verified a posteriori, namely, we have to apply the iterative method to find the sequence h n = ϕ(h n−1 ) and for n large enough to verify if S is satisfied on h n .

SOME EXAMPLES
In this section we illustrate the C 1 contraction mapping method given by Maldonado and Moreira (2003) to approximate the stationary solution of the Euler equation.The first one shows how it is applied in linear models.The second example shows the numerical robustness of the method even though the stationary solution has sensitiveness to the initial condition.The last example illustrates how the method runs in stochastic models.

The linear model
The first example is the simple linear model mostly used as a first approximation of non-linear models.It also works when objective functions are quadratic.Consider the homogeneous linear equations: M x t−1 + N x t + P x t+1 = 0, where M, N, P ∈ R n×n .If N is nonsingular, then the equation above can be written as: If ||A|| + ||B|| < 1 and ||B|| < 1/2, then D is satisfied and, so it is the saddle stability point of x = 0.The iterative method runs as follows: given a differentiable function h n with ||Dh n || < 1 we define h n+1 by so that the dynamics will evolve through linear functions.Finally, the algorithm is: given any matrix M 0 ∈ R n×n , with ||M 0 || < 1, the sequence: will converge to M ∈ R n×n ,with ||M || < 1, and the stationary solution will be g(x) = M x. 7 7 It is important to remark that the iterative method only approximates solutions with norm lower than one.For instance, if A = 0 and B = −0.25,then the Euler equation has two solutions g 1 (x) = 0 and g 2 (x) = 4x.Obviously the method will select g 1 .
Solving Euler Equations: Classical Methods and the C 1 Contraction Mapping Method Revisited

An example with chaotic stationary solution
In this example we consider a deterministic dynamic programming problem such the optimal policy function is g(x) = 4x(1 − x) in the interval [0, 1].This example is important in order to test the robustness of our iterative method.Since g exhibits sensitiveness to initial conditions (positive Liapounov exponent), then the calculated iterations could be far from the real iteration.
According to Boldrin and Montrucchio (1986) we can explicitly define a deterministic model that has g as its stationary solution.In this case we must choose a discount factor small enough.Specifically, taking E(x, y, z) = x − x 2 + δz − (0.25 + 4δ + 2δz)y + 12δy 2 − 8δy 3 as the Euler Equation, it is easy to check that g is the stationary solution of this equation for all δ > 0 sufficiently small.The figure below shows the evolution of the sequence of policy functions, where the optimal one is the continuous line and the iterations h 1 to h 10 (which is the closest to g) are the dotted lines; the distance between h 10 and g is 7.863 × 10 −4 .The discount factor we chose was δ = 0.1.In this example the assumption D is not satisfied, since the policy function has an interior unstable steady state.Nevertheless, our method still works and the convergence occurs in all the domain.It suggests that the method might work even if D is not satisfied.This may be a subject for future research.

The stochastic growth model
In this last example we apply the iterative method generated from Theorem 3.2 to the neoclassical stochastic growth model.To do that, recall that the stochastic Euler equation for the stochastic dynamic Wilfredo L. Maldonado, Humberto Moreira programming problem is: The specifications for the growth model are: technological shocks are i.i.d, the technology is given by the production function f and the utility function is u.So, F (x, y, z) = u(zf (x) − y).
The algorithm is defined by the map ϕ : Given h 0 ∈ C 1 (X × Z), define the sequence h n+1 = ϕ hn which may converge to optimal policy function of the problem.To state explicitly the map ϕ, let us consider the following specification: f (x) = x and u(c) = c − (b/2)c 2 .Thus, the (contraction) mapping ϕ is: Knowing the parameter values β ∈ (0, 1), b > 0 and the distribution of the shock z we can perform the iterates of ϕ from h 0 ∈ C(X × Z).By Theorem 3.2 those iterations will converge to the optimal policy function of the problem.

CONCLUSIONS
The use of accurate numerical methods in economic dynamic models has an increasing demand among applied researchers.The particular characteristic of those models (when the optimal path is found by iterating the stationary solution function) requires a strong control of the errors in each iteration.
This paper provides a survey on the classical numerical methods to estimate the stationary solution of the Euler equations.Among them, Bellman's iteration method, projection methods and contraction methods provide the most popular numerical algorithms to solve those equations.We also test the robustness of the method defined by Maldonado and Moreira (2003) by applying it to solve the dynamic programming problem which has the logistic map as the optimal policy function.This test shows that even though the optimal policy has sensitiveness to initial conditions the method has a good performance to approximate the solution.Finally, we extend that method to the stochastic Euler equations showing that under hypothesis that can be verified a posteriori we can obtain the convergence to the stationary solution of those equations.
As mentioned in Maldonado and Moreira (2003), the C 1 contraction method has two main advantages over the others in the literature: The convergence is in the C 1 topology (allowing in some cases for accurate comparative dynamics exercises) and it is not necessary any grid or interpolation in each iteration.For the stochastic case the same properties are satisfied.
For this method to work, the required hypothesis in the deterministic case is an open condition related with the first derivatives of the structural equations evaluated at the steady state.This hypothesis is slightly stronger than the determinacy of the steady state condition.For the stochastic case the hypothesis is a natural extension of that and can be verified ex-post.