Worked Examples in Dynamic Optimization:

Analytic and Numeric Methods

Laurent Cretegny∗

Centre of Policy Studies, Monash University,

Australia

Thomas F. Rutherford†

Department of Economics, University of Colorado

USA

March 29, 2004

Abstract
Economists are accustomed to think about economic growth models in con-
tinuous time. However, applied models require numerical methods because of
the absence of tractable analytical solutions. Since these methods operate by
essence in discrete time, models involve discrete formulation. We demonstrate
the usefulness of two off-the-shelf algorithms to solve these problems : nonlin-
ear programming and mixed complementarity. We then show the advantage
of the latter for approximating infinite-horizon models.

JEL classification: C69; D58; D91

Keywords: Dynamic optimization; Mathematical methods; Infinite-horizon

models

∗ Mailing address : Centre of Policy Studies, PO Box 11E, Monash University, Clayton Vic 3800,

Australia. E-mail: laurent@[Link]. The author gratefully acknowledges financial support

from the Swiss National Science Foundation (post-doctoral research fellowship).
† Mailing address : Department of Economics, University of Colorado, Boulder, USA. E-mail:

rutherford@[Link].
1 Introduction
Dynamic optimization in economics appeared in the 1920s with the work of Hotelling
and Ramsey. In the 1960s dynamic mathematical techniques became then more fa-
miliar to economists mainly due to the work of neoclassical growth theorists. These
techniques involve most of the time formulation of models in continuous time. When
closed form solutions do not exist they are then formulated in discrete time. The
purpose of this document is to provide some sample solutions of a collection of
dynamic optimization problems in two settings, using analytical methods in contin-
uous time and numerical methods in discrete time.
Formulation of infinite-horizon models are not possible with numerical methods.
Therefore approximation issues are crucial in finite-horizon models. We consider two
classes of off-the-shelf algorithms to solve these dynamic models. The first is non-
linear programming (NLP) developed originally for optimal planning models. The
second class is the mixed complementarity problem (MCP) approach. The MCP
formulation is represented by the first-order conditions for nonlinear programming.
Hence any NLP problem can be solved as an MCP formulation, not necessarily as
efficient as using NLP-specific methods.
Approximating infinite-horizon models is illustrated in figure 1. The two inner
circles represent the idea that the finite MCP formulation includes any of the NLP
formulations. These two finite formulations are a subset of the infinite-horizon NLP
formulation. It is then intuitively clear that an MCP formulation should provide a
”better” approximation to infinite-horizon models than an NLP formulation. The
closeness of approximation is informally portrayed by the Euclidian distance in the
figure.

Figure 1: Approximating infinite-horizon models

The outline of the paper is as follows. Starting from the classical mathematical
technique to solve dynamic economizing problems in continuous time, the next
section shows how to derive the NLP and MCP formulation to solve these problems.
Section 3 presents in detail analytical solutions to economic planning problems and
shows how to formulate them in off-the-shelf softwares. The following section moves
on to the neoclassical growth model. The last section explains how to use the
optimal neoclassical growth model in applied economics.

1
2 Mathematical Methods
The dynamic economizing problem may be solved in three different approaches.
The first approach going back up to Bernoulli in the very late 1600s is the calculus
of variations. The second is the maximum principal developed in the 1950s by Pon-
tryagin and his co-workers. The third approach is dynamic programming developed
by Bellman about the same time.
Early applications of dynamic optimization to economics are due to Ramsey
and Hotelling in the 1920s. At that time the mathematical technique used to solve
dynamic problems was the calculus of variations. Therefore in the following section
we first state in a concise way the calculus of variations problem. Then we move
on to the maximum principle which can be considered a dynamic generalization of
the method of Lagrange multiplier. This method is well-known among economists
and is especially suited to the formulation in discrete time. Regarding dynamic
programming it is usually applied to stochastic models and then will not be covered
here.

2.1 Continuous time approach

The classical calculus of variations problem may be written as
Z t1
max J = I(x(t), x0 (t), t)dt
{x(t)} t0
subject to various initial and endpoint conditions

where these conditions are defined as follow:

a. Euler equation: Fx = dFx0 /dt, t0 ≤ t ≤ t1 .
b. Legendre condition: Fx0 x0 ≤ 0, t0 ≤ t ≤ t1 .
c. Boundary conditions:
- Initial conditions always apply: x(t0 ) = x0 .
- The terminal time and terminal value may be fixed exogenously or free.
d. Transversality conditions apply when the terminal value and time are free:
- If only the terminal value is free, then Fx0 = 0 at t1 .
- If only the terminal time is free, then F − x0 Fx0 = 0 at t1 .
- If both the terminal value and time are free, then F = 0 and Fx0 = 0 at t1 .
These necessary conditions of the calculus of variations can be derived from
the maximum principle. Intuitively it remains to let the rates of change of the
state variables to be the control variables in the maximum principle, which means
u(t) = x0 (t). Assuming that the terminal time value is fixed, which is always the
case in numerical problems, the corresponding maximum principle may be defined
as
Z t1
max J = I(x(t), u(t), t)dt + F (x (t1 ) , t1 )
{u(t)} t0
subject to x0 (t) = f (x(t), u(t), t)
t0 , t1 and x (t0 ) = x0 fixed
x (t1 ) = g (x (t1 ) , t1 ) or free

2
where I(·) is the intermediate function, F (·) is the final function, f (·) is the state
equation function and g(·) is the terminal constraint function.
In a concise way the maximum principle technique involves adding costate vari-
ables λ(t) to the problem, defining a new function called the Hamiltonian,

H(x(t), u(t), λ(t), t) = I(x, u, t) + λ(t)f (x, u, t)

and solving for trajectories {u(t)}, {λ(t)}, and {x(t)} satisfying the following con-
ditions
∂H
optimality condition = Iu + λfu = 0
∂u
∂H
costate equation λ0 = − = − (Ix + λfx )
∂x
∂H
state equation x0 = = f with x (t0 ) = x0
∂λ
∂F
terminal conditions · x (t1 ) ≥ 0 ⊥ λ (t1 ) ≥
∂x
∂F ∂g
· λ (t1 ) = + λ̃
∂x ∂x
with λ̃ ≥ 0 ⊥ x (t1 ) = g (x (t1 ) , t1 )

which are necessary for a local maximum.

2.2 Discrete time formulation

The formulation of the discrete time version of the maximum principle is straight-
forward. Forming the Hamiltonian,

H(xt , ut , λt+1 , t) = I(x, u, t) + λt+1 f (x, u, t)

the necessary conditions are as follow:

∂H
optimality condition = Iu + λt+1 fu = 0
∂ut
∂H
costate equation λt+1 − λt = − = − (Ix + λt+1 fx )
∂xt
∂H
state equation xt+1 − xt = = f with xt0 = x0
∂λt+1
∂F
terminal conditions · xt1 +1 ≥ 0 ⊥ λt1 +1 ≥
∂x
∂F ∂g
· λt1 +1 = + λ̃
∂xt1 +1 ∂xt1 +1
with λ̃ ≥ 0 ⊥ xt1 +1 = g (xt1 +1 , t1 + 1)

As mentioned earlier, the maximum principle can be considered the extension

of the method of Lagrange multipliers to dynamic optimization problems. This
method allows us to state problems in the same way they would be written in off-
the-shelf softwares. Write L for the Lagrangian of the full intertemporal problem.
The NLP formulation is then
t1
X
max L= I (xt , ut , t) + λt+1 [xt + f (xt , ut , t) − xt+1 ]
{u(t)},{λ(t)},{x(t)}
t=t0

+ λt0 [xt0 − x0 ] + λ̃t1 +1 [g (xt1 +1 , t1 + 1) − xt1 +1 ] + F (xt1 +1 , t1 + 1)

3
and the MCP formulation follows from the first-order conditions
∂L
= Iu (xt , ut , t) + λt+1 fu (xt , ut , t) = 0 t = t 0 , . . . , t1
∂ut
∂L
= Ix (xt , ut , t) + λt+1 fx (xt , ut , t) + λt+1 − λt = 0 t = t 0 , . . . , t1
∂xt
∂L
= −λt1 +1 + λ̃t1 +1 (gx − 1) + Fx = 0
∂xt1 +1
∂L
= f (xt , ut , t) − (xt+1 − xt ) = 0 t = t 0 , . . . , t1
∂λt+1
∂L
= xt0 − x0 = 0
∂λt0
∂L
= g (xt1 +1 , t1 + 1) − xt1 +1 = 0
∂ λ̃t1 +1

which are the necessary conditions of the maximum principle.

The method of Lagrange multipliers shows clearly that when the system has a
fixed final state, as here, there are two constraints for the terminal period t1 +1: the
state variable xt1 +1 must satisfy the terminal constraint and still satisfies the state
equation. This explains the two Lagrange multipliers associated with xt1 +1 : λt1 +1
for the state equation, and λ̃t1 +1 for the terminal constraint. When the system
has a free final state, which means that the terminal constraint is not specified, the
Lagrange multiplier λt1 +1 is equal to zero if the value of the final function is zero.

3 Economic Planning Models

3.1 Optimal consumption plan
1
Find the consumption plan C(t), 0 ≤ t ≤ T , over a fixed period to maximize the
discounted utility stream
Z T
e−rt C a (t)dt subject to C(t) = iK(t) − K 0 (t), K(0) = K0 , K(T ) = 0
0

where 0 < a < 1 and K represents the capital stock.

Analytic solution
We have a variational problem based on the function:

F (t, k, k 0 ) = e−rt U (ik − k 0 )

Taking derivatives, we have:

Fk = ie−rt U 0
and
Fk0 = −e−rt U 0
The Euler equation is then:
d £ −rt 0 ¤
−e U = ie−rt U 0
dt
1 Problem 4.5 in Kamien and Schwartz (2000).

4
Although the underlying problem is defined in terms of the capital stock, it is
convenient at this point to use consumption as the decision variable, when we form
the time derivative we have:
re−rt U 0 − e−rt U 00 c0 = ie−rt U 0
which reduces to:
U 00 c0
− =i−r
U0
In the case of constant elasticity utility,
U (c) = ca ,
the Euler equation becomes:
c0
(1 − a)=i−r
c
and, integrating we determine the growth rate of consumption:
i−r
c = c0 e 1−a t
In this expression, initial consumption level is a constant of integration. In order
to determine the consumption level, we need to focus on the initial and terminal
conditions for the capital stock. To determine the time path of the capital stock, we
consider the equation which relates consumption to capital earnings and investment:
i−r
ik − k 0 = c0 e 1−a t
In order to solve an equation of this form, it is necessary to use a standard method
for solving this sort of an equation, multiplying by an integrating factor:
e−it [k 0 − ik] = −c0 eθt
in which we define:
i−r ai − r
θ= −i=
1−a 1−a
This equation can be written:
£ ¤
d e−it k(t) = −c0 eθt dt
which integrates to:
eθt
e−it k(t) = γ − c0
θ
so:
c0 1−a
i−r
t
k(t) = γeit −
e
θ
We then have two boundary conditions to determine the constants of integration:
k(0) = k0 , k(T ) = 0
The initial condition produces:
c0
γ = k0 +
θ
Substituting into the terminal condition, we have:
θk0
c0 =
eθT − 1
Finally, substitute the integrating constants back into the expression for the capital
stock to obtain: · θT ¸
e − eθt
k(t) = k0 eit
eθT − 1

5
 Numeric solution
Working in discrete time, the following code sets up the model as a nonlinear op-
timization problem. The first solution is used to compare results from the analytic
and numeric models. The second set of solutions evaluate the qualitative properties
of the consumption path for alternative elasticities parameters, a. The final calcu-
lation in this program presents an alternative representation of the choice problem
as budget-constrained welfare maximization.
1 $title Kamien and Schwartz, problem 4.5 - NLP formulation
2
3 sets t time periods / 0*60 /
4 decade(t) decades / 10, 20, 30, 40, 50 /
5 tfirst(t) first period of time
6 tlast(t) last period of time;
7
8 tfirst(t) = yes$(ord(t) eq 1);
9 tlast(t) = yes$(ord(t) eq card(t));
10
11 scalars r discount rate / 0.03 /
12 i interest rate / 0.04 /
13 a utility coefficient / 0.5 /;
14
15 variables c(t) consumption level
16 k(t) capital stock
17 kt terminal capital stock
18 u utility function;
19
20 equations market(t) market clearance in period t
21 market_t terminal market clearance
22 const_kt terminal capital constraint
23 utility objective function definition;
24
25 market(t).. k(t) - k(t-1) =e= 1$tfirst(t) + i*k(t-1) - c(t-1);
26
27 market_t.. (kt - sum(tlast, k(tlast))) =e= sum(tlast, i*k(tlast) - c(tlast));
28
29 const_kt.. kt =e= 0;
30
31 utility.. u =e= sum(t, (1/(1+r))**((ord(t)-1)) * c(t)**a);
32
33 model ramsey / all /;
34
35 * Lower bound to avoid domain errors
36
37 [Link](t) = 0.001;
38
39 * Numeric solution
40
41 parameters compare comparison of analytic and numeric solution;
42
43 solve ramsey using nlp maximizing u;
44
45 compare(t,’numeric c’) = c.l(t);
46 compare(t,’numeric k’) = k.l(t);
47
48 * Analytic solution
49
50 scalars theta, c0;
51
52 theta = (a * i - r) / (1 - a);
53 c0 = theta / (exp(theta * (card(t)-1)) - 1);
54
55 compare(t,’analytic c’) = c0 * exp((theta+i)*(ord(t)-1));
56 compare(t,’analytic k’) = exp(i*(ord(t)-1)) * (exp(theta*(card(t)-1))
57 -exp(theta*(ord(t)-1))) / (exp(theta*(card(t)-1))-1);
58

6
 59 * Alternative elasticities parameters
60
61 sets elasval alternative elasticity values / ’0.3’,’0.6’,’0.9’/;
62
63 parameters consum consumption path for alternative elasticities;
64
65 a = 0;
66 loop(elasval,
67 a = 0.3 + a;
68 solve ramsey using nlp maximizing u;
69 consum(t,elasval) = c.l(t);
70 );
71
72 * Alternative representation of the choice problem
73
74 parameter p(t) present value of consumption in period t;
75
76 p(t) = (1/(1 + i))**(ord(t)-1);
77
78 equations budget present-value budget constraint;
79
80 budget.. sum(t, p(t) * c(t)) =e= 1 + i;
81
82 model altmodel /utility, budget/;
83
84 a = 0.5;
85 solve altmodel using nlp maximizing u;
86
87 compare(t,"altmodel c") = c.l(t);
88
89 $if %batch%==yes $setglobal batch yes
90 $if %batch%==yes $setglobal gp_opt1 "set term postscript eps monochrome ’Times-Roman’ 20"
91 $if %batch%==yes $setglobal gp_opt2 "set title"
92
93 $setglobal domain t
94 $setglobal labels decade
95
96 $if %batch%==yes $setglobal gp_opt3 "set output ’[Link]’"
97 $libinclude plot compare
98
99 $if %batch%==yes $setglobal gp_opt3 "set output ’[Link]’"
100 $setglobal gp_opt4 "set key left"
101 $libinclude plot consum

7
 Figure 2: Analytic and Numeric Solutions
1.2
numeric c
numeric k
analytic c
1 analytic k
altmodel c

0.8

0.6

0.4

0.2

0
10 20 30 40 50

Figure 3: Consumption Paths

0.7
0.3
0.6
0.6 0.9

0.5

0.4

0.3

0.2

0.1

0
10 20 30 40 50

8
 Figure 2 compare the state variable time path between analytic and numeric
models. They don’t seem to be identical especially the value of capital at the final
period. The reason comes from the treatment of time being a discrete succession of
periods. It results that, the final period in continuous-time model corresponds to
the end of the final period in discrete-time model, which is the beginning of period
t1 + 1. Since the continuous time is the limit of discrete periods shrinking to zero,
differences between the two approaches are reduced when smaller periods of time
are considered. An illustration is given below.
1 $title Kamien and Schwartz, problem 4.5 - Smaller time periods
2
3 sets t time periods / 0*120 /
4 m(t) main time periods / 0 /
5 tfirst(t) first period of time
6 tlast(t) last period of time;
7
8 tfirst(t) = yes$(ord(t) eq 1);
9 tlast(t) = yes$(ord(t) eq card(t));
10
11 scalars r discount rate / 0.03 /
12 i interest rate / 0.04 /
13 a utility coefficient / 0.5 /
14 dt increment of time subperiod / 2 /;
15
16 loop(t$m(t), m(t+dt)=yes; );
17
18 variables c(t) consumption level
19 k(t) capital stock
20 kt terminal capital stock
21 u utility function;
22
23 equations market(t) market clearance in period t
24 market_t terminal market clearance
25 const_kt terminal capital constraint
26 utility objective function definition;
27
28 market(t).. dt * (k(t) - k(t-1)) =e= (1*dt)$tfirst(t) + i*k(t-1) - c(t-1);
29
30 market_t.. (kt - sum(tlast, k(tlast))) =e= sum(tlast, i*k(tlast) - c(tlast));
31
32 const_kt.. kt =e= 0;
33
34 utility.. u =e= sum(t, (1/(1+r))**((ord(t)-1)/dt) * c(t)**a / dt);
35
36 model ramsey / all /;
37
38 * Lower bound to avoid domain errors
39
40 [Link](t) = 0.001;
41
42 * Do a comparison of numeric and analytic solutions
43
44 parameters compare comparison of analytic and numeric solution;
45
46 solve ramsey using nlp maximizing u;
47
48 compare(m,’numeric c’) = c.l(m);
49 compare(m,’numeric k’) = k.l(m);
50
51 scalars theta, c0;
52
53 theta = (a * i - r) / (1 - a);
54 c0 = theta / (exp(theta * (card(t)-1)/dt) - 1);
55
56 compare(m(t),’analytic c’) = c0 * exp((theta+i)*(ord(t)-1)/dt);
57 compare(m(t),’analytic k’) = exp(i*(ord(t)-1)/dt) * (exp(theta*(card(t)-1)/dt)

9
58 -exp(theta*(ord(t)-1)/dt)) / (exp(theta*(card(t)-1)/dt)-1);
59
60 option decimals = 8;
61 display compare;
62
63 sets decade / 20 ’10’, 40 ’20’, 60 ’30’, 80 ’40’, 100 ’50’ /;
64
65 $if %batch%==yes $setglobal batch yes
66 $if %batch%==yes $setglobal gp_opt1 "set term postscript eps monochrome ’Times-Roman’ 20"
67 $if %batch%==yes $setglobal gp_opt2 "set title"
68
69 $setglobal labels decade
70
71 $if %batch%==yes $setglobal gp_opt3 "set output ’[Link]’"
72 $libinclude plot compare m

Figure 4: Solution for 60 years with 120 time periods

1.2
numeric c
numeric k
analytic c
1 analytic k

0.8

0.6

0.4

0.2

0
10 20 30 40 50

MCP formulation
It may be useful to represent prices explicitly in the model. Below are therefore the
MPSGE ((Rutherford, 1999)) and algebraic models using the MCP formulation.
1 $title Kamien and Schwartz, problem 4.5 - MPSGE and MCP formulation
2
3 sets t time periods / 0*60 /
4 tfirst(t) first period of time
5 tlast(t) last period of time;
6
7 tfirst(t) = yes$(ord(t) eq 1);
8 tlast(t) = yes$(ord(t) eq card(t));
9
10 scalars r discount rate / 0.03 /
11 i interest rate / 0.04 /;
12
13 variables c(t) consumption level
14 u utility function;
15
16 positive variables k(t) capital stock
17 kt terminal capital stock;
18
19 equations market(t) market clearance in period t

10
20 market_t terminal market clearance
21 const_kt terminal capital constraint
22 utility objective function definition;
23
24 market(t).. k(t) - k(t-1) =e= 1$tfirst(t) + i*k(t-1) - c(t-1);
25
26 market_t.. (kt - sum(tlast, k(tlast))) =e= sum(tlast, i*k(tlast) - c(tlast));
27
28 const_kt.. kt =e= 0;
29
30 utility.. u =e= sum(t, (1/(1+r))**((ord(t)-1)) * log(c(t)));
31
32 model ks_nlp / all /;
33
34 * Lower bound to avoid domain errors
35
36 [Link](t) = 0.001;
37
38 * NLP solution
39
40 solve ks_nlp using nlp maximizing u;
41
42 * MPSGE formulation
43
44 alias (t,t1);
45
46 parameters theta budget share over time
47 epsilon budget share over time (lagged)
48 k_nlp capital value from NLP solution;
49
50 theta(t) = (1/(1+r))**(ord(t)-1)/sum(t1,(1/(1+r))**(ord(t1)-1));
51 epsilon(t+1) = theta(t);
52 k_nlp(t) = k.l(t);
53
54 $ontext
55
56 $model:ks_mge
57
58 $sectors:
59 k(t) ! capital stock
60
61 $commodities:
62 pk(t) ! price of capital stock
63 pkt ! price of terminal capital stock
64
65 $consumers:
66 ra ! representative agent
67
68 $auxiliary:
69 kt ! terminal capital stock
70 pktc ! price of constraint terminal capital stock
71
72 $prod:k(t)
73 o:pk(t+1) q:(1+i)
74 o:pkt$tlast(t) q:(1+i)
75 i:pk(t) q:1
76
77 $demand:ra s:1.0
78 e:pk(t)$tfirst(t) q:1
79 d:pk(t)$(not tfirst(t)) q:epsilon(t)
80 d:pkt q:(sum(tlast,theta(tlast)))
81
82 $constraint:kt
83 pkt =e= pktc;
84
85 $constraint:pktc
86 kt =e= 0;

11
 87
88 $report:
89 v:cons(t) d:pk(t) demand:ra
90 v:cons_t d:pkt demand:ra
91
92 $offtext
93 $sysinclude mpsgeset ks_mge
94
95 * NLP values to initialize the model
96
97 k.l(t) = k_nlp(t);
98 pk.l(t) = market.m(t);
99 pkt.l = market_t.m;
100 pktc.l = -const_kt.m;
101
102 * MPSGE solution
103
104 ks_mge.iterlim = 0;
105 $include ks_mge.gen
106 solve ks_mge using mcp;
107
108 * MCP formulation
109
110 equations pr_k(t) zero profit condition for capital stock
111 pr_kt zero profit condition for terminal capital stock
112 demand(t) demand function;
113
114 pr_k(t).. (1+i) * (pk(t+1) + pkt$tlast(t)) =e= pk(t);
115
116 pr_kt.. pkt =e= pktc;
117
118 demand(t).. c(t) * (pk(t+1) + pkt$tlast(t)) =e=
119 theta(t) * sum(tfirst, pk(tfirst)*k(tfirst));
120
121 MODEL ks_mcp / pr_k.k, pr_kt.kt, [Link], market_t.pkt, const_kt.pktc, demand.c /;
122
123 * MCP solution
124
125 ks_mcp.iterlim = 0;
126 solve ks_mcp using mcp;

3.2 The monopolist

2
The demand function for a monopolist depends on both the product price and the
rate of change of the product price, according to:

x = a0 p + b0 + c0 p0

Assume that the cost of production at rate x is:

C(x) = a1 x2 + b1 x + c1

Given the initial price, p(0) = p0 , and the required ending price, p(T ) = pT , find
the price policy over 0 ≤ t ≤ T which maximizes profits:
Z T
[px − C(x)]dt
0

Analytic solution
Notice that because the time period is fixed, the fixed term in the cost function,
c1 , is irrelevant if the firm is committed to produce; so we will ignore that term to
2 Problem 5.4 in Kamien and Schwartz (2000).

12
conserve on algebra.
Substituting with the demand function, we see that this problem corresponds to a
calculus of variations problem in which the function depends only on the price and
the gradient of price, but not on the time path, i.e.

F (p, p0 ) = px(p, p0 ) − C(p, p0 )

Neglecting constants, this reduces to:

F (p, p0 ) = a0 (1 − a0 a1 )p2 − a1 c0 2p02 + c0(1 − 2a0 a1 )pp0 + (b0 − 2a0 a1 b0 − a0 b1 )p

Some differentiation:

Fp = 2a0 (1 − a0 a1 )p + c0 (1 − 2a0 a1 )p0 + b0 − 2a0 a1 b0 − a0 b1

and
dFp0
= c0 (1 − 2a0 a1 )p0 − 2a1 c0 2p00
dt
The Euler equation is then a second-order, linear differential equation with constant
coefficients:
p00 + Bp = R
where
a0 (1 − a0 a1 )
B=
a1 c0 2
and
a0 b1 + 2a0 a1 b0 − b0
R=
2a1 c0 2
Notice that in steady-state, where p00 = 0, the Euler condition implies that:
R a0 b1 + 2a0 a1 b0 − b0
p∗ = =
B 2a0 (1 − a0 a1 )

which is equivalent the optimal monopoly price in the static equilibrium.

If we are to assume that the static equilibrium model is based on a downward sloping
demand function and a convex technology, then:

a0 < 0, b0 > 0, a1 > 0, and b1 > 0

Hence, we have may conclude:

p∗ > 0, B<0

In order to solve the differential equation, we begin with the adjacent homogeneous
system:
p00 + Bp = R
We know that the solution of this equation has the form:

p(t) = cert

Hence:
cert (r2 + B) = 0
Defining: s
a0 (a0 a1 − 1)
r̂ =
a1 c0 2

13
 The solution to the non-homogeneous equation therefore has the form:

p(t) = c1 er̂t + c2 e−r̂t + c3

and follows from the definition of r̂ that

R
c3 = = p∗
B
And boundary conditions determine c1 and c2 as solutions to the following system
of equations:
c1 + c2 = p0 − p∗ , c1 er̂T + c2 e−r̂T = pT − p∗
When the initial and final prices are both equal to the static monopoly price, c1 =
c2 = 0 and the optimal policy is to keep the price fixed over the time horizon.
If the terminal price equals the optimal static value, then over the horizon the
price moves monotonically from the initial value to the terminal value (when the
terminal price equals the p∗ , then c1 and c2 are of opposite sign).

Numeric solution
Working in discrete time, we can formulation this model as a nonlinear optimization
problem and solve it using GAMS/MINOS, as illustrated in the following code:
1 $title Kamien and Schwartz, problem 5.4 - NLP formulation
2
3 set t /1*100/,
4 decade(t) /10,20,30,40,50,60,70,80,90/;
5
6 scalar sigma elasticity of demand /4/
7 eta elasticity of supply /0.25/
8 c0 multiplier /-20/
9 a0,a1,b1, r,coef1,coef2;
10
11 parameter
12 compare comparison of numerical approximation with analytic solution,
13 pricepath approach path for prices from various starting points,
14 turnpike illustrating turnpike property of the optimal price path;
15
16 * Impute a1 and b1 so that we have a steady-state with
17 * the price and quantity both equal to unity:
18
19 b1 = (1 - 1/sigma) * (1 - 1/eta);
20 a1 = (1/2) * (1 - 1/sigma - b1);
21 a0 = -sigma;
22
23 * Declare the model:
24
25 variables p(t) price
26 x(t) quantity
27 c(t) cost
28 profit maximand;
29
30 equations demand, cost, objdef;
31
32 * By declaring equations over t+1, we omit equations for the
33 * first period in which the price is fixed exogenously:
34
35 demand(t).. x(t) =e= (1+sigma) - sigma * p(t) + c0 * (p(t)-p(t-1));
36
37 cost(t).. c(t) =e= a1 * x(t)*x(t) + b1*x(t);
38
39 objdef.. profit =e= sum(t, x(t) * p(t) - c(t));
40
41 * Create a model with all of these equations:

14
42
43 model dynamic /all/;
44
45 * Fix terminal period price at the equilibrium price:
46
47 [Link]("100") = 1;
48
49 * Fix initial period values:
50
51 [Link]("1") = 0.5;
52 solve dynamic using nlp maximizing profit;
53
54 compare(t,"numeric") = p.l(t);
55
56 r = sqrt( a0 * (a0 * a1 - 1) / (a1*c0*c0) );
57 coef2 = (p.l("1") - 1) / (1 - exp(-2 * r * 99));
58 coef1 = - coef2 * exp(-2 * r * 99);
59 compare(t,"analytic") = 1 + coef1 * exp(r * (ord(t)-1) ) + coef2 * exp(-r * (ord(t)-1));
60
61 * Create a set defined by either the initial or terminal period price:
62
63 set p0 /"0.1","0.3","0.5","0.7","0.9"/;
64
65 loop(p0,
66 [Link]("1") = 0.1 + 0.2 * (ord(p0)-1);
67 solve dynamic using nlp maximizing profit;
68 pricepath(t,p0) = p.l(t);
69 );
70
71 * Now illustrate the turnpike property:
72
73 [Link]("1") = 0.6;
74 loop(p0,
75 [Link]("100") = 0.1 + 0.2 * (ord(p0)-1);
76 solve dynamic using nlp maximizing profit;
77 turnpike(t,p0) = p.l(t);
78 );
79
80 * Display the results using GNUPLOT:
81
82 $if %batch%==yes $setglobal batch yes
83 $if %batch%==yes $setglobal gp_opt1 "set term postscript eps monochrome ’Times-Roman’ 20"
84 $if %batch%==yes $setglobal gp_opt2 "set title"
85
86 $setglobal domain t
87 $setglobal labels decade
88
89 $if %batch%==yes $setglobal gp_opt3 "set output ’[Link]’"
90 $libinclude plot compare
91
92 $if %batch%==yes $setglobal gp_opt3 "set output ’[Link]’"
93 $libinclude plot pricepath
94
95 $if %batch%==yes $setglobal gp_opt3 "set output ’[Link]’"
96 $setglobal gp_opt4 "set key bottom left"
97 $libinclude plot turnpike

15
 Figure 5: Analytic and Numeric Solutions
1
numeric
0.95 analytic

0.9

0.85

0.8

0.75

0.7

0.65

0.6

0.55

0.5
10 20 30 40 50 60 70 80 90

Figure 6: Approach Paths for Various Starting Points

1
0.1
0.3
0.9 0.5
0.7
0.8 0.9

0.7

0.6

0.5

0.4

0.3

0.2

0.1
10 20 30 40 50 60 70 80 90

16
 Figure 7: Turnpike Property for Optimal Price Path
1

0.9

0.8

0.7

0.6

0.5

0.4

0.3 0.1
0.3
0.2 0.5
0.7
0.9
0.1
10 20 30 40 50 60 70 80 90

3.3 Non-renewable resource

3
Suppose a mine contains an amount B of a mineral resource (like coal, copper or
oil). The profit rate that can be earned from selling the resource at rate x is ln x.
Find the rate at which the resource should be sold over the fixed period [0,T] to
maximize the present value of profits from the mine. Assume the discount rate a
constant r. Assume the resource has no value beyond time T .

Analytic solution
Following the hint, define y(t) as the cumulative sales by time t. Then y 0 (t) is the
sales rate at time t Find y(t) to:
Z T
max e−rt ln y 0 (t)dt
0

subject to:
y(0) = 0, y(T ) = B
We therefore have:
µ ¶
dFy0 −e−rt y 00
F (t, y, y 0 ) = e−rt ln y 0 (t), Fy = 0 and = 0 +r
dt y (t) y0
The Euler equation then gives us the differential equation:
y 00
= −r
y0
Integrating, we have:
y(t) = c1 e−rt + c2
Then applying the boundary conditions, we have:
e−rt − 1
y(t) = B
e−rT − 1
3 Problem 5.5 in Kamien and Schwartz (2000).

17
 Numeric solution
Working in discrete time, we can formulation this model as a nonlinear optimization
problem and solve it using GAMS/MINOS, as illustrated in the following code:
1 $title Kamien and Schwartz, problem 5.5 - NLP formulation
2
3 set t /0*100/,
4 decade(t) /10,20,30,40,50,60,70,80,90/;
5
6 scalar r interest rate /0.05/;
7
8 variables profit present value of extraction
9 x(t) production at time t;
10
11 equations objdef defines profit
12 supply defines cumulative extraction;
13
14 objdef.. profit =e= sum(t, exp(-r * (ord(t)-1)) * log(x(t)));
15
16 supply.. sum(t, x(t)) =e= 1;
17
18 model hotelling /all/;
19
20 [Link](t) = 0.00001;
21 x.l(t) = 1/card(t);
22
23 solve hotelling using nlp maximizing profit;
24
25 parameter compare comparison of numeric and analytic solutions
26 y(t) cumulative extraction in the analytic solution;
27
28 y(t) = (exp(-r * (ord(t)-1)) - 1) / (exp(-r * 100) - 1);
29
30 compare(t,"numeric") = x.l(t);
31 compare(t,"analytic") = y(t+1) - y(t);
32 compare("100","analytic") = 0;
33 display compare;
34
35 $if %batch%==yes $setglobal batch yes
36 $if %batch%==yes $setglobal gp_opt1 "set term postscript eps monochrome ’Times-Roman’ 20"
37 $if %batch%==yes $setglobal gp_opt2 "set title"
38
39 $setglobal domain t
40 $setglobal labels decade
41
42 $if %batch%==yes $setglobal gp_opt3 "set output ’[Link]’"
43 $libinclude plot compare

18
 Figure 8: Analytic and Numeric Solutions
0.05
numeric
0.045 analytic

0.04

0.035

0.03

0.025

0.02

0.015

0.01

0.005

0
10 20 30 40 50 60 70 80 90

3.4 Non-renewable resource (general case)

4
Reconsider the previous problem but suppose that the profit rate is P (x) when the
resource is sold at rate x, where P 0 (0) > 0 and P 00 < 0.

1. Show that the present value of the marginal profit from extraction is constant
over the planning period (otherwise it would be worthwhile to shift the time
of sale of a unit of the resource from a less profitable moment to a more
profitable one). Marginal profit, P 0 (t) therefore grows exponentially at the
discount rate r.
2. Show that the optimal extraction rate declines through time.

Analytic solution
We have a calculus of variations problem in which:

F (t, y, y 0 ) = e−rt P (y 0 (t)), Fy = 0 and Fy0 = e−rt P 0 (y 0 ).

The Euler condition therefore implies:

dFy0 de−rt P 0 (y 0 )
= =0
dt dt
or, in answer to question 1:

e−rt P 0 (y 0 ) = constant

Then if P 00 < 0, then only way that P 0 (y 0 ) increases at an exponential rate r over
time is that the extraction rate, y 0 , must be declining through time.
4 Problem 5.6 in Kamien and Schwartz (2000).

19
 Numeric solution
As the demand curve becomes more elasticity, the production profile must decline
at a faster rate so that the present value of the marginal from extraction remains
constant over the planning period.
1 $title Kamien and Schwartz, problem 5.6 - NLP formulation
2
3 set t /0*100/,
4 decade(t) /10,20,30,40,50,60,70,80,90/;
5
6 scalar r interest rate /0.05/,
7 sigma elasticity of demand /0.5/;
8
9 variables profit present value of extraction
10 x(t) production at time t;
11
12 equations objdef defines profit
13 supply defines cumulative extraction;
14
15
16 objdef.. profit =e= sum(t, exp(-r * (ord(t)-1)) * x(t)**(sigma-1)/sigma );
17
18 supply.. sum(t, x(t)) =e= 1;
19
20 model hotelling /all/;
21
22 [Link](t) = 0.00001;
23 x.l(t) = 1/card(t);
24
25 set sigval /"1.0","1.2","1.4"/
26
27 parameter extract Extraction profile over time;
28
29 loop(sigval,
30 sigma = 0.81 + 0.2 * ord(sigval);
31 solve hotelling using nlp maximizing profit;
32 extract(t,sigval) = x.l(t);
33 );
34
35 $if %batch%==yes $setglobal batch yes
36 $if %batch%==yes $setglobal gp_opt1 "set term postscript eps monochrome ’Times-Roman’ 20"
37 $if %batch%==yes $setglobal gp_opt2 "set title"
38
39 $setglobal domain t
40 $setglobal labels decade
41
42 $if %batch%==yes $setglobal gp_opt3 "set output ’[Link]’"
43 $libinclude plot extract

20
 Figure 9: Extraction Values for Alternative Values of σ
0.09
1.0
1.2
0.08 1.4

0.07

0.06

0.05

0.04

0.03

0.02

0.01

0
10 20 30 40 50 60 70 80 90

3.5 Pollution control

Utility U (C, X) increases with the consumption rate C and decreases with the stock
of pollution, X. For C > 0, P > 0,

UC > 0, UCC < 0, lim UC = ∞;

C→0
UX < 0, UXX < 0, lim UX = 0; UCX = 0.
X→0

The constant rate of output Y is to be divided between consumption and pollution

control. Consumption contributes to pollution, while pollution control reduces it;
Z(C) is the net contribution to the pollution flow, with Z 0 > 0, Z 00 > 0. For
small C, little pollution is created and much abated; thus net pollution declines:
Z(C) < 0. But for large C, considerable pollution is created and few resources
remain for pollution control, therefore on net pollution increases: Z(C) > 0. Let
C ∗ be the consumption rate that satisfies Z(C ∗ ) = 0. In addition, the environment
absorbs pollution at a constant proportionate rate b. Characterize the consumption
path C(t) that maximizes the discounted utility stream:
Z ∞
e−rt U (C, X)dt
0

subject to

X 0 = Z(C) − bX, X(0) = X0 , 0 ≤ C ≤ Y, 0≤X

Also characterize the corresponding optimal pollution path and the steady state.
This kind of problems are typically the ones which are much more convenient
to solve with numerical methods rather than analytically.
1 $title Kamien and Schwartz, problem II.8.5 - NLP formulation
2
3 sets t time periods / 0*10 /
4 tfirst(t) first period of time
5 tlast(t) last period of time;
6
7 tfirst(t) = yes$(ord(t) eq 1);

21
 8 tlast(t) = yes$(ord(t) eq card(t));
9
10 scalars r discount rate / 0.03 /
11 b rate of pollution decay / 0.05 /
12 psi disutility rate of pollution / 0.15 /
13 alpha pollution control parameter / 5 /
14 beta function curvature parameter / 16 /
15 y constant rate of output / 8 /
16 x0 initial stock of pollution / 30 /;
17
18 variables u utility function;
19
20 positive variables c(t) consumption level
21 x(t) pollution stock
22 xt terminal capital stock;
23
24 equations steq_x(t) state equation of pollution
25 steq_xt state equation of terminal pollution
26 appr_xt approximation of terminal pollution
27 utility objective function definition;
28
29 steq_x(t).. x(t) - x(t-1) =e= x0$tfirst(t)
30 + (-alpha + beta / (y - c(t-1)) - b * x(t-1))$(not tfirst(t));
31
32 steq_xt.. (xt - sum(tlast, x(tlast))) =e=
33 - alpha + sum(tlast, beta / (y - c(tlast)) - b * x(tlast));
34
35 appr_xt.. xt =e= 150;
36
37 utility.. u =e= sum(t, (1/(1+r))**((ord(t)-1)) * (log(c(t)) - psi * log(x(t))));
38
39 model pollution / all /;
40
41 * Lower bound to avoid domain errors
42
43 [Link](t) = 0.001;
44 [Link](t) = y-0.001;
45 [Link](t) = 0.001;
46
47 * NLP solution
48
49 solve pollution using nlp maximizing u;

4 The Neoclassical Growth Model

4.1 Factor shares
5
For a neoclassical function, show that each factor of production earns its marginal
product. Show that if owners if capital save all their income and workers consume
all their income, then the economy reaches the golden rule of capital accumulation.
Explain the results.

Analytic solution
The neoclassical function and its properties:

Y = F (K, L)

Non-negative and diminishing marginal products:

FK ≥ 0, FKK < 0, FL ≥ 0, FLL < 0

5 Problem 1.5 in Barro and Sala-I-Martin (2004).

22
Constant returns to scale:

F (λK, λL) = λF (K, L)

Inada conditions assuring an interior solution:

lim FK = ∞, lim FK = 0
K→0 K→∞

lim FL = ∞, lim FL = 0
L→0 L→∞

The production function may then be expressed in intensive form:

Y = LF (K/L, 1) ≡ Lf (k)

where k = K/L, or y = f (k) where y = Y /L. The marginal production of capital:

∂Y ∂yL ∂y
= = = f 0 (k)
∂K ∂kL ∂k
The marginal product of labour:
∂Y ∂yL ∂y ∂k
= =y+L = f (k) + Lf 0 (k) = f (k) − f 0 (k)(K/L) = f (k) − kf 0 (k)
∂L ∂L ∂L ∂L
The firm’s objective is to maximize profits defined as:

max Π ≡ F (K, L) − wL − rK

Dividing this expression by L, we have:

max π ≡ f (k) − w − rk

The first order condition for k is:

f 0 (k) = r

Under constant returns to scale, all revenue is returned to capital and labour:

f (k ∗ ) = w + rk ∗

Substituting for r, we determine the wage rates which results in zero profit:

w∗ = f (k ∗ ) − k ∗ f 0 (k ∗ )

We see that this wage is precisely the marginal product of labour. Hence, when
firms maximize profits constant returns to scale assures that profits are driven to
zero.
Assume now that all capital income is fully reinvested, so:

I ∗ = r∗ K ∗

Also assume that all labour is consumed:

C ∗ = w∗ L

We therefore have:
Ċ L̇
= =n
C L

23
 If we define c = C/L, then we have:
ċ
=0
c
The laws of motion for capital in the Solow-Swann model are defined as:

K̇ = I − δK = r∗ K − δK

so it follows that:
k̇ = r∗ k − (δ + n)k
On a steady-state growth path, we have:

K̇ Ċ L̇
= = =n
K C L
hence, k̇ = 0, and
r∗ − (δ + n) = 0
Substituting for the marginal product of capital, we recover the Golden Rule con-
dition:
f 0 (k) = δ + n

Numeric solution
1 $title Barro and Sala-i-Margin, problem 1.5 - Solow-Swann Growth Model
2
3 * Declare the time horizon here. It is important to
4 * declare this set as ordered (using the "*") so that
5 * we can reference "t+1" from inside the loop over t:
6
7 set t /0*200/, t0 /0/;
8
9 * Data describing the base year are given here:
10
11 scalar alpha base year capital value share /0.6/
12 r0 base year gross return to capital /0.12/
13 delta capital depreciation rate /0.07/
14 n labor growth rate /0.02/
15 s_L savings rate of workers /0.10/
16 s_K savings rate of capital owners /0.90/
17
18 * The elasticity of substitution between labor and capital is a free
19 * parameter. For the initial simulation, we set it to unity (almost).
20 * (We use 1.01 rather than 1 in order to avoid having to change the
21 * functional form from Cobb-Douglas to CES in the model):
22
23 sigma elasticity of substitution /1.01/
24
25 * Calibrated parameters:
26
27 c0 base year consumption (calibrated)
28 rho primal CES exponent
29 k0 base year capital stock (calibrated)
30 l0 base year labor supply (calibrated);
31
32 * The following parameters hold an equilibrium time path:
33
34 parameter k Time path of capital stock
35 l Time path of labor
36 c Time path of total consumption
37 y Time path of output
38 r Time path of return to capital
39 w Time path of wage rate

24
 40
41 * The following parameters hold output to be plotted:
42
43 timepath Time path of per-capita variables,
44 output Time path of output (alternative sigma values),
45 return Time path of return (alternative sigma values),
46 consum Time path of consumption (alternative sigma values),
47 wage Time path of wage (alternative sigma values);
48
49 * Compute the primal elasticity exponent:
50
51 rho = 1 - 1/sigma;
52
53 * Calibrate the base year capital stock and labor supply
54 * in efficiency units, taking base year output equal to unity
55 * and measuring labor in efficiency units:
56
57 k0 = alpha / r0;
58 l0 = (1-alpha);
59
60 * Base year consumption is based on capital and labor earnings
61 * shares and the marginal propensity to save out of those income
62 * sources:
63
64 c0 = alpha * (1-s_K) + (1-alpha) * (1-s_L);
65
66 * Initialize base year (time 0) output, capital and labor stock:
67
68 y(t0) = 1;
69 k(t0) = k0;
70 l(t0) = l0;
71
72 * Do an initial simulation with the specified value of simga
73 * (1.01 = Cobb Douglas).
74
75 loop(t,
76
77 * Entering period t the values of capital and labor are known, so the
78 * output is known:
79
80 y(t) = ( alpha * (k(t)/k0)**rho + (1-alpha) * (l(t)/l0)**rho)**(1/rho);
81
82 * The return to capital and labor are computed as marginal products:
83
84 r(t) = (y(t)*k0/k(t))**(1/sigma) * alpha / k0;
85 w(t) = (y(t)*l0/l(t))**(1/sigma) * (1-alpha) / l0;
86
87 * Consumption is the sum of consumption levels by capital owners and
88 * workers:
89
90 c(t) = (1-s_K) * r(t) * k(t) + (1-s_L) * w(t) * l(t);
91
92 * Capital evolves through depreciation and investment:
93
94 k(t+1) = k(t) * (1 - delta) + y(t) - c(t);
95
96 * Labor growth is exogenous at rate n:
97
98 l(t+1) = l(t) * (1 + n);
99 );
100
101 * Store the time path of key values for plotting:
102
103 timepath(t,"output") = y(t) * l0 / l(t);
104 timepath(t,"consum") = (c(t)/c0) * l0 / l(t);
105 timepath(t,"return") = r(t) * k0 / alpha;
106 timepath(t,"wage") = w(t) * l0 / (1-alpha);

25
107
108 * Declare a set over values of the elasticity of substitution (sigma)
109 * to be compared:
110
111 set sigval /"0.5", "1.0", "2.0" /;
112
113 parameter sigvalue(sigval) / "0.5" 0.5, "1.0" 1.01, "2.0" 2.0 /;
114
115 loop(sigval,
116
117 * Assign the elasticity:
118
119 sigma = sigvalue(sigval);
120 rho = 1 - 1/sigma;
121
122 * Compute the equilibrium time path (period 0 values are the same in all
123 * simulations):
124
125 loop(t,
126 y(t) = ( alpha * (k(t)/k0)**rho + (1-alpha) * (l(t)/l0)**rho)**(1/rho);
127 r(t) = (y(t)*k0/k(t))**(1/sigma) * alpha / k0;
128 w(t) = (y(t)*l0/l(t))**(1/sigma) * (1-alpha) / l0;
129 c(t) = (1-s_K) * r(t) * k(t) + (1-s_L) * w(t) * l(t);
130 k(t+1) = k(t) * (1 - delta) + y(t) - c(t);
131 l(t+1) = l(t) * (1 + n);
132 );
133
134 * Save some values to plot comparisons:
135
136 consum(t,sigval) = (c(t)/c0) * l0 / l(t);
137 output(t,sigval) = y(t) * l0 / l(t);
138 return(t,sigval) = r(t) * k0 / alpha;
139 wage(t,sigval) = w(t) * l0 / (1-alpha);
140
141 );
142
143 * Generate some labeled plots:
144
145 set tics(t) / 0, 25, 50, 75, 100, 125, 150, 175, 200 /
146
147 $if %batch%==yes $setglobal batch yes
148 $if %batch%==yes $setglobal gp_opt1 "set term postscript eps monochrome ’Times-Roman’ 20"
149 $if %batch%==yes $setglobal gp_opt2 "set title"
150
151 $setglobal gp_xl tics
152 $setglobal gp_xlabel years
153 $setglobal domain t
154 $setglobal labels tics
155
156 $if %batch%==yes $setglobal gp_opt3 "set output ’[Link]’"
157 $setglobal gp_opt4 "set key outside"
158 $setglobal gp_opt5 "set xlabel ’years’"
159 $setglobal gp_opt6 "set ylabel ’% change’"
160 $libinclude plot timepath
161
162 $if %batch%==yes $setglobal gp_opt3 "set output ’[Link]’"
163 $libinclude plot output
164
165 $if %batch%==yes $setglobal gp_opt3 "set output ’[Link]’"
166 $libinclude plot consum
167
168 $if %batch%==yes $setglobal gp_opt3 "set output ’[Link]’"
169 $libinclude plot return
170
171 $if %batch%==yes $setglobal gp_opt3 "set output ’[Link]’"
172 $libinclude plot wage

26
 Figure 10: Key Variables
1.5
output
consum
1.4 return
wage
1.3

1.2
% change

1.1

0.9

0.8

0.7
0 25 50 75 100 125 150 175 200
years

Figure 11: Output

2.4
0.5
1.0
2.2 2.0

2
% change

1.8

1.6

1.4

1.2

1
0 25 50 75 100 125 150 175 200
years

27
 Figure 12: Return to Capital
1
0.5
1.0
2.0
0.95

0.9
% change

0.85

0.8

0.75
0 25 50 75 100 125 150 175 200
years

Figure 13: Wage Rate

1.6
0.5
1.0
2.0
1.5

1.4
% change

1.3

1.2

1.1

1
0 25 50 75 100 125 150 175 200
years

28
4.2 Distortions in the Solow-Swan model
6
Assume that output is produced by the CES production function,

Y = [(aF KFη + aI KIη )φ/η + aG KG

φ 1/φ
]

where Y is output; KF is formal capital, which is subject to taxation; KI is informal

capital, which evades taxation; KG is public capital, provided by government and
used freely by all producers; aF , aI , aG > 0; η < 1, and φ < 1. Installed formal
and informal capital differ in their location and form of ownership and, therefore,
in their productivity.
Output can be used on a one-for-one basis for consumption or gross investment
in the three types of capital. All three types of capital depreciate at the rate δ.
Population is constant, and technology progress is nil.
Formal capital is subject to tax at the rate τ at the moment of its installation.
Thus, the price of formal capital (in units of output) is 1 + τ . The price of a unit of
informal capital is one. Gross investment in public capital is the fixed fraction sG
of tax revenues. Any unused tax receipts are rebated to households in a lump-sum
manner. The sum of investment in the the two forms of private capital is the faction
s of income net of taxes and transfers. Existing private capital can be converted on
a one-to-one basis in either direction between formal and informal capital.

a. Derive the ratio of informal to formal capital used by profit-maximizing produc-

ers.
b. In the steady-state, the three forms of capital grow at the same rate. What is
the ratio of output to formal capital in the steady-state?
c. What is the steady-state growth rate of the economy?
d. Numerical simulations show that, for reasonable parameter values, the graph of
the growth rate against the tax rate, τ , initially increases rapidly, then reaches a
peak, and finally decreases steadily. Explain this nonmonotonic relation between
the growth rate and the tax rate.

Analytic solution
The problem states that output is given by the CES production function
y = f (kF , kI , kG )
h i1/ψ
ψ/η
= (aF kFη + aI kIη ) ψ
+ aG kG
·³ ³ ´η ´ψ/η ³ ´ψ ¸1/ψ
kI
= kF aF + aI kF + aG kkG
F
,

where k denotes capital and subscripts F , I and G denote formal, informal and
government, respectively; population growth is constant and technological progress
is nil; depreciation is the same for all forms of capital, implying that

k̇F = iF − δkF ,
k̇I = iI − δkI , and
k̇G = iG − δkG ,
where i denotes the investment at time t for the kind of capital specified by the
subscript; taxes are collected as a fixed fraction of formal investment,

T = τ iF ;
6 Problem 1.7 in Barro and Sala-I-Martin (2004) based on Easterly (1993).

29
gross investment in public capital is a fixed fraction of taxes,

iG = sG T = sG τ iF ;

unused taxes,
TU = (1 − sG )T = (1 − sG )τ iF ,
are rebated to households in a lump-sum manner; prices (in units of output) are

PF = (1 + τ ) and PI = 1 ;

the sum of private investments in formal and informal capital is given by

iF + iI = s(Y − T + TU ) = s(Y − sG τ iF ) ,

which implies
iI = s(Y − (1 + sG τ )iF ) ;
Treating public capital as an externality, a profit maximizing producer chooses
formal and informal investments to solve

max {Py y − PF iF − PI iI } ≡ max {y − (1 + τ )iF − iI }

iF ,iI iF ,iI

subject to
def k̇F iF
gF = = − δ , and
kF kF
def k̇I iI
gI = = −δ ,
kI kI
or equivalently,
iF iI
kF = , and kI = .
gF + δ gI + δ
The first order conditions for this problem are
∂y dkF ∂y dkI
= (1 + τ ) and =1,
∂kF diF ∂kI diI
which implies
µ ¶µ ¶−1 µ ¶η−1 µ ¶
∂y dkI ∂y dkF kI aI (gF + δ) 1
= = .
∂kI diI ∂kF diF kF aF (gI + δ) 1+τ

a. From the first order conditions, the ratio of informal to formal capital used by
profit-maximizing producers can be computed to be
· ¸1/(η−1) · ¸1/(1−η)
kI aF (gI + δ) (1 + τ )aI (gF + δ)
= = .
kF (1 + τ )aI (gF + δ) aF (gI + δ)

b. If the three forms of capital grow at the same rate so that

gF = gI = gG ,

then
· µ ¶¸1/(1−η)
kI aI kG iG
= (1 + τ ) and = = sG τ ,
kF aF kF iF

30
 and thus
·³ ´ψ/η ¸1/ψ
η/(1−η) ψ
y = kF βF + βI (1 + τ ) + βG τ ,

where
µ ¶η/(1−η)
def def aI def
βF = aF , βI = aI and βG = aG sψ
G .
aF
The ratio of output to formal capital in the steady state is therefore
·³ ´ψ/η ¸1/ψ
y η/(1−η) ψ
= βF + βI (1 + τ ) + βG τ .
kF

c. Since the ratio of output to formal capital in the steady state is constant, it must
be the case that the growth rate of the economy is equal to the growth rate of
formal capital, which is given by
µ ¶
y
g = s −δ .
kF
The growth rate for the economy is therefore
g(τ ) = sG(τ )1/ψ − δ
where
·³ ´ψ/η ¸
def
G(τ ) = βF + βI (1 + τ )η/(1−η) + βG τ ψ > 0.

d. An economic explanation for simulated nonmonotonic behaviour for the change

in growth as a function of taxes is as follows:
As taxes grow from zero, the public good becomes available but capital
is also moved from the formal to the informal sector. For reasonable
parameter values, the increased productivity due to the public good
dominates the loss of productivity due to the transfer from formal
capital to informal capital, and so the economy grows.
At some point, the growth with respect to taxes is maximized, indicat-
ing that the loss of productivity due to movement away from formal
capital to informal capital is exactly offset by the increase in the pro-
ductivity due to public capital.
As the tax rate increases beyond the maximal point and approaches
unity, the steady state growth will decrease as the low productivity
informal capital dominates the increase in public sector productivity.
However, µ ¶
s
g 0 (τ ) = G(τ )(−1+1/ψ) G0 (τ )
ψ
where
·µ ¶³ ´−1+ψ/η ¸
βI
G0 (τ ) = ψ βF + βI (1 + τ )η/(1−η) (1 + τ )−1+η/(1−η) + βG τ ψ−1 .
1−η
The sign of g 0 is always positive, since the sign of G is always positive and the
sign of G0 is the same as the ratio of savings rate to the elasticity of substitution
between private and public capital. This analytic result thus indicates that
growth is always increasing in the tax rate, thereby contradicting the simulated
nonmonotonic behaviour discussed above.

31
 Numeric solution
1 $title Barro and Sala-i-Margin, problem 1.7 - Easterly model
2
3 * This programs illustrates how to calibrate the Easterly
4 * model, evaluate how the steady-state growth rate depends on
5 * the tax rate, and then evaluates the transition path for a
6 * change in the tax rate
7
8
9 set taxrate Alternative tax rates to evaluate (%) / 1*200 /,
10 taxlabel(taxrate) /10,30,50,70,90,110,130,150,170,190 /,
11 t Years to simulate /1997*2041/,
12 return Assumed returns to public capital /low, medium, high/,
13 tplot Time periods to plot /1997*2040/,
14 decade(tplot) / 2000, 2010, 2020, 2030, 2040/;
15
16 scalar
17
18 *===============================================================================
19 * Base year data are specified here:
20
21 tau0 benchmark tax rate on formal sector investment /0.50/
22 tau_s tax rate on formal sector in simulated adjustment /0.90/
23 r_g benchmark relative return to public sector capital /1.0/
24 g baseline growth rate /0.03/
25 s_g public savings rate /0.75/
26 delta depreciation rate /0.07/
27 iratio ratio of informal to formal capital /0.15/
28 sigmag substitution elasticity between private and public capital /0.5/
29 sigmaf substitution elasticity between formal and informal capital /4.0/
30
31 *===============================================================================
32 * Calibrated or temporary parameters:
33
34 s private savings rate
35 thetag implicit benchmark value share of public capital
36 thetaf share of private capital in the formal sector
37 tau tax rate in counter-factual
38 k_f baseline formal capital
39 k_i baseline informal capital
40 k_g baseline public capital
41 k_p private sector capital stock (state variable)
42 y0 scale parameter in benchmarking
43 a_f CES share parameter for formal capital
44 a_i CES share parameter for informal capital
45 a_g CES share parameter for public capital
46 eta CES exponent (inner nest)
47 psi CES exponent (outer nest)
48
49 ki_ratio ratio of informal to formal capital
50 yf_ratio ratio of output to formal capital;
51
52 parameter
53
54 r0(return) Implicit baseyear relative return to public capital
55
56 growth Steady-state growth rate (sensitivity to return on public capital)
57
58 y Output level in simulated transition
59 kg Public capital in simulated transition
60 kf Formal capital in simulated transition
61 ki Informal capital in simulated transition
62
63 transition Growth rates through the transition;
64
65
66 *===============================================================================

32
 67 * Benchmarking steps are explained here:
68
69 eta = 1 - 1/sigmaf;
70 psi = 1 - 1/sigmag;
71
72 * For purpose of deriving coefficients, set magnitude of formal
73 * capital to unity:
74
75 k_f = 1;
76 k_i = iratio;
77
78 thetaf = k_f / (k_f + k_i);
79
80 * Given formal capital, we know the public capital stock from the
81 * tax rate and public sector savings rate:
82
83 k_g = s_g * tau0 * k_f;
84
85 * Value share of public capital depends on assumed shadow return
86 * to public sector capital:
87
88 thetag = r_g * k_g / (k_f + (1+tau0)*k_i + r_g * k_g);
89
90 * Calibrate public capital coefficient from the base year quantity,
91 * the value share and the elastiicty:
92
93 a_g = k_g**(-psi) * thetag / (1 - thetag);
94
95 * Calibrate relative size of formal and informal coefficients,
96 * based on relative size of the capital stock and the elasticity:
97
98 a_f = 1;
99 a_i = (1 / (1 + tau0)) * (k_i/k_f)**(1-eta);
100
101 * Now compute the implicit output level and rescale capital stocks to
102 * be consistent with benchmark output equal to unity:
103
104 y0 = ( (a_f * k_f**eta + a_i * k_i**eta)**(psi/eta) + a_g * k_g**psi )**(1/psi);
105
106 k_f = k_f / y0;
107 k_i = k_i / y0;
108 k_g = k_g / y0;
109
110 * Calibrate private savings to be consistent with steady-state:
111
112 s = (g + delta) * (1 + tau0 + iratio) / (1/k_f);
113
114
115 *===============================================================================
116 * 1) Simulate the transitional dynamics associated with a change
117 * in the tax rate.
118
119 * Apply the new tax rate:
120
121 tau = tau_s;
122
123 * Given the tax rate, we know the formal share of private capital use:
124
125 thetaf = 1 / ((a_i * (1 + tau) / a_f)**(1/(1-eta)) + 1);
126
127 * Initialize the state variable:
128
129 k_p = k_f + k_i;
130
131 loop(t,
132
133 * Record current capital stocks:

33
134
135 kg(t) = k_g;
136 kf(t) = k_f;
137 ki(t) = k_i;
138
139 * Compute output:
140
141 y(t) = ( (a_f * k_f**eta + a_i * k_i**eta)**(psi/eta)
142 + a_g * k_g**psi )**(1/psi);
143
144 * Update capital stocks for the next period.
145
146 * Note: the price of a unit of new capital is given by a share-weighted
147 * average of the informal price (1) and the formal sector price (1+tau).
148 * This explains the termin the the denominator of the k_p expression:
149
150 k_p = k_p * (1 - delta) + s * y(t) / (1 - thetaf + (1+tau) * thetaf);
151 k_g = k_g * (1 - delta) + s_g * tau * thetaf * s * y(t) / (1 - thetaf
152 + (1+tau) * thetaf);
153 k_f = k_p * thetaf;
154 k_i = k_p * (1 - thetaf);
155 );
156
157 * Compute growth rates:
158
159 transition(t,"g0")$y(t+1) = 100 * g;
160 transition(t,"y")$y(t+1) = 100 * (y(t+1)-y(t)) / y(t);
161 transition(t,"kg")$kg(t+1) = 100 * (kg(t+1)-kg(t)) / kg(t);
162 transition(t,"kf")$kf(t+1) = 100 * (kf(t+1)-kf(t)) / kf(t);
163 transition(t,"ki")$ki(t+1) = 100 * (ki(t+1)-ki(t)) / ki(t);
164
165 display transition;
166
167
168 *===============================================================================
169 * 2) Perform a sensitivity analysis: growth as a function of the
170 * tax rate, accounting for alternative assumptions regarding
171 * the base year shadow price on public capital:
172
173 * Assign base year relative returns to public capital
174 * (base year return to informal capital = 1):
175
176 r0("low") = 0.5;
177 r0("medium") = r_g;
178 r0("high") = 1 + tau0 + 0.5;
179
180 * For each alternative base year return to public capital,
181 * recalibrate the model:
182
183 loop(return,
184
185 r_g = r0(return);
186 eta = 1 - 1/sigmaf;
187 psi = 1 - 1/sigmag;
188 k_f = 1;
189 k_i = iratio;
190 k_g = s_g * tau0 * k_f;
191 a_f = 1;
192 a_i = (1 / (1 + tau0)) * (k_i/k_f)**(1-eta);
193 thetag = r_g * k_g / (k_f + (1+tau0)*k_i + r_g * k_g);
194 a_g = k_g**(-psi) * thetag / (1 - thetag);
195 y0 = ( (a_f * k_f**eta + a_i * k_i**eta)**(psi/eta) + a_g * k_g**psi )
196 **(1/psi);
197 k_f = k_f / y0;
198 k_i = k_i / y0;
199 k_g = k_g / y0;
200 s = (g + delta) * (1 + tau0 + iratio) / (1/k_f);

34
201 display s;
202
203 * Then for each model, evaluate how the growth rate
204 * changes with the tax rate:
205
206 loop(taxrate,
207 tau = 0.01 * ord(taxrate);
208 ki_ratio = (a_i * (1 + tau) / a_f)**(1/(1-eta));
209 yf_ratio = ( (a_f + a_i * (ki_ratio)**eta )**(psi/eta)
210 + a_g * (s_g * tau)**psi )**(1/psi);
211 growth(taxrate,return) = 100 * (s * yf_ratio
212 / (1 + tau + ki_ratio) - delta);
213 );
214
215 );
216
217 display growth;
218
219
220 *===============================================================================
221 * Generate some plots:
222
223 $if %batch%==yes $setglobal batch yes
224 $if %batch%==yes $setglobal gp_opt1 "set term postscript eps monochrome ’Times-Roman’ 20"
225 $if %batch%==yes $setglobal gp_opt2 "set title"
226
227 $if %batch%==yes $setglobal gp_opt3 "set output ’[Link]’"
228 $setglobal gp_opt4 "set key outside width 4"
229 $setglobal gp_opt5 "set xlabel ’Tax rate on formal capital (%)’"
230 $setglobal gp_opt6 "set ylabel ’Economic growth rate (%)’"
231 $setglobal gp_opt7 "set grid"
232 $setglobal gp_opt8 "set yrange [-2:4]"
233
234 $setglobal domain taxrate
235 $setglobal labels taxlabel
236 $libinclude plot growth
237
238 $if %batch%==yes $setglobal gp_opt3 "set output ’[Link]’"
239 $setglobal gp_opt5 "set xlabel ’Year’"
240 $setglobal gp_opt6 "set ylabel ’Economic growth rate (%)’"
241 $setglobal gp_opt8 "set yrange [-1:5]"
242
243 $setglobal domain tplot
244 $setglobal labels decade
245 $libinclude plot transition

35
 Figure 14: Capital Taxes and Steady-State Growth
4
low
medium
high
3
Economic growth rate (%)

-1

-2
10 30 50 70 90 110 130 150 170 190
Tax rate on formal capital (%)

Figure 15: Growth Rates through the Transition

5
g0
y
kg
4 kf
ki
Economic growth rate (%)

-1
2000 2010 2020 2030 2040
Year

36
5 The Neoclassical Optimal Growth Model
This section lays down the basics for developing applied dynamic CGE models. We
begin by going through the logic of the Ramsey model which is often presented as
a dynamic optimization problem7 :
X∞ µ ¶t 1−θ
1 Ct − 1
max
t=0
1 + ρ 1−θ

s.t.
Ct = f (Kt ) − It
Kt+1 = (1 − δ)Kt + It
K0 = K̄0
The maximand in this problem is often called constant-elasticity-of-intertemporal
substitution (CEIS) utility function. As will be shown below, it simply represents
a monotonic transformation of conventional CES utility function.
Here, as in many macroeconomics textbooks, aggregate output is expressed as
a function of the capital stock alone, i.e.:

Yt = f (Kt )

In the MPSGE representation of the Ramsey model8 , it is convenient to work

with a constant-returns production function in which we have inputs of both labour
and capital:
Yt = F (L̄t , Kt )
When labour is in fixed supply, the production function exhibits diminishing
returns to capital. There is therefore no loss of generality by formulating the model
on the basis of a constant returns to scale technology.
In writing down a model it is helpful to employ the unit cost function associated
with the production function F (·):9

c(pL K L K
t , rt ) ≡ min pt aL + rt aK

s.t.
F (aL , aK ) = 1
Shephard’s lemma tells us that the compensated demand functions for labour
and capital are the partial derivatives of the unit cost function:

∂c(pL K
t , rt )
aK (rK , pL ) = K
∂rt
and
∂c(pL K
t , rt )
aL (rK , pL ) =
∂pLt
The representative agent model can be formulated as a general equilibrium
model which is completely routine, apart from the fact that there are an infinite-
number of variables. Following the conventional GAMS/MPSGE framework, equi-
librium in the model is characterized by three classes of equations:
7 Forsimplicity it is assumed that there is no population growth.
8 Seethe last part in section 3.1 for a simple model represented in three equivalent formulations,
i.e. NLP, algebraic MCP and MPSGE.
9 Note that the lower-case function c(·) represents unit cost, while the upper case C represents
t
consumption in year t. In the equilibrium model Ct (p, M ) represents the demand for output in
year t as a function of output prices and aggregate present value of income.

37
 1. Market clearance conditions and associated market prices are as follows:10
- Output market (market price pt ):

Yt = Ct (p, M ) + It

- Labour market (wage rate pL

t ):

L̄t = aL (rtK , pL
t ) Yt

- Market for capital services (capital rental rate rtK ):

Kt = aK (rtK , pL
t ) Yt

- Capital stock (capital purchase price pK

t ):

Kt+1 = (1 − δ)Kt + It

2. Zero profit conditions and associated activities are:11

- Output (Yt ):
pt = c(pL K
t , rt )

- Investment (It ≥ 0):

pt ≥ pK
t+1

- Capital stock (Kt ):

pK K K
t = rt + (1 − δ)pt+1

3. Income balance:
∞
X
M = pK
0 K̄0 + pL
t L̄t
t=0

Two questions might arise for an MPSGE modeler looking at this equilibrium
model. First, the careful observer might note that the demand functions, Ct (p, M ),
have not been specified, and because these arise from CIES preferences so there
may be some details to work out. This problem is considerably easier than the
second issue, namely how do we solve an infinite-dimensional system of nonlinear
equations. Let’s first look at this latter issue. The issue of CEIS preferences will be
considered in the calibration section below.
In order to solve a finite approximation of the model with a T -period model
horizon, we need to decompose the consumer’s problem. Consider the infinite-
horizon problem of the representative agent in Ramsey’s model:
X∞ µ ¶t
1
max u(ct )
t=0
1+ρ

s.t.
∞
X ∞
X
pc ct = pK
0 K̄0 + pL
t L̄t
t=0 t=0
10 The demand functions employed in this model assure that all prices will be nonzero in equilib-

rium. There is no formal need, therefore, to associated prices with market clearance conditions, as
would be required in a conventional complementarity problem. We provide an associated here in
order to help understand how the model might be extended with demand functions which would
admit zero prices.
11 The only activity level which could possibly fall to zero would be investment, and that would

only happen in a policy scenario which resulted in a substantial reduction in the return to capital.

38
 c1−θ −1
in which u(c) = 1−θ , and define a value of terminal assets to be:
∞
X ¡ ¢
A∗T = pc c∗t − pL
t L̄t
t=T +1

Then consider the equivalent model:

XT µ ¶t ∞ µ
X ¶t
1 1
max u(ct ) + u(ct )
t=0
1+ρ 1+ρ
t=T +t

s.t.
T
X T
X
pc ct = pK
0 K̄0 + pL
t L̄t − AT
t=0 t=0
∞
X ∞
X
pc ct = AT + pL
t L̄t
t=T +1 t=T +1

If AT is fixed then this can be posed as two separate optimization problems, one
running through time period T and another for the post-terminal period. When
terminal assets are assigned a value of A∗T , corresponding to the infinite-horizon
solution, then the finite horizon model will then produce consumption levels for
years 0 through T which are identical to the ∞-horizon model. The question is how
do we find A∗T ?
Terminal assets in the closed economy model are simply equal to the value of
the capital stock at the start of period T + 1. The model running through year
T then produces a good approximation to the consumer problem when we have a
good approximation to the terminal capital stock. The key insight provided by Lau,
Pahlke, and Rutherford (2002) is that the state variable KT +1 can be determined as
part of the equilibrium calculation by targeting the associated control variable, IT .
In the present model this could be based on any of the following primal constraints:
- Terminal investment growth rate set equal to the long-run steady-state growth
rate:
IT /IT −1 = 1 + g

- Terminal investment growth rate set equal to the growth rate of aggregate
output:
IT /IT −1 = YT /YT −1

- Terminal investment growth rate set equal to the growth rate of consumption:
IT /IT −1 = CT /CT −1

State-variable targeting provides a very compact means of determining the ter-

minal capital stock. In models with multiple consumers living beyond period T , it
would be necessary to account for which of these agents owns the assets. Note that
some agents may have negative asset positions at the end of the model – particularly
in overlapping generations models where young households accumulate debt which
is repaid in middle age.
The final detail involved in implementing a dynamic model in MPSGE is cali-
bration. The simplest approach is to set up the model along a steady-state growth
rate in which the interest rate (r̄) and growth rate (ḡ) are given. The first thing to
work out is to determine the structure of the benchmark equilibrium.
Here are the steps involved in sorting out the steady-state conditions which
related investment and capital earnings in a static data set which is consistent with
a steady-state growth path:

39
 1. The zero-profit condition for It reveals the price level for capital:

pK
t
pK
t+1 = = pt
1 + r̄
hence
pK
t = (1 + r̄)pt

The base year price of capital is then:

p̄K = 1 + r̄

2. The zero profit condition for Kt determines the price level for rtK :

pK K K
t = rt + (1 − δ)pt+1

Substituting the values of pK K

t and pt+1 reveals that the base year rental price
of capital is sufficient to cover interest plus depreciation:

r̄K = r̄ + δ

3. The main challenge involved in calibrating a dynamic model centers on the

reconciliation of base year capital earnings, investment, the steady-state in-
terest rate and the capital depreciation rate. To see how this works, consider
the market clearance condition for capital in the first period:

K1 = K̄0 (1 − δ) + I¯ = (1 + ḡ)K̄0

This implies that base year investment can be calculated on the basis of growth
and depreciation of the base year capital stock:

I¯ = K̄0 (ḡ + δ)

Finally, we can use r̄K to determine K̄0 on the basis of the value of capital
earnings in the base year, V K, hence:
ḡ + δ
I¯ = V K
r̄ + δ

The problem that arises in applied models is that I¯ and V K will not satisfy
this relation for arbitrary values of ḡ, r̄ and δ. Something typically has to be
adjusted to match up the dataset with the baseline growth path.
The second issue to work out is the representation of CEIS preferences in a
MPSGE model. Consider the following equivalent representations of intertemporal
preferences:
1. Additively separable utility:
X∞ µ ¶t 1−θ
1 Ct − 1
U (C) =
t=0
1 + ρ 1−θ

2. Linearly homogeneous utility:

"∞ µ ¶t # 1−θ
1
X 1
Û (C) = Ct1−θ
t=0
1 + ρ

40
 It is possible to determine the equivalence of U and Û by recalling that a mono-
tonic transformation of utility does not alter the underlying preference ordering.
Observe that:
1/a
Û = V (U ) = [aU + κ]
where
X∞ µ ¶t
1 1+ρ
κ= = ,
t=0
1+ρ ρ
and
a = 1 − θ.
V (·) is a monotonic transformation (V 0 > 0), hence optimization of U and Û yield
identical demand functions.
Alternatively, recall that preference orderings are defined by the marginal rate
of substitution. In both of these models we have:
µ ¶θ
∂U/∂Ct+1 1 Ct
=
∂U/∂Ct 1+ρ Ct+1

There are several advantages associated with the use of linearly homogeneous
representation. First of all, these preferences can be represented in MSPGE. Second,
the reporting of welfare changes as Hicksian-equivalent variations is trivial with Û
: a 1% change in Û corresponds to a 1% equivalent variation in income.
CEIS preferences over a finite horizon can be represented in MPSGE as follows
(lines 110 to 112 in the code given below):
$PROD:U s:sigma
O:PU Q:(c0*sum(t, pref(t)*qref(t)))
I:P(t) Q:(qref(t)*C0) P:pref(t)

Intertemporal preferences in an MPSGE model are typically based on the fol-

lowing parameters:

- c0 is the base year consumption level,

- qref(t) = (1+g0)**(ord(t)-1) is the baseline equilibrium index of economic
activity, calculated on the basis of a steady-state growth rate equal to g0,
- pref(t) = (1/(1+r0))**(ord(t)-1) is the baseline present value price path,
calculated on the basis of a steady-state interest rate equal to r0, and
- sigma is the intertemporal elasticity of substitution.
Figure 5 shows how the utility function is calibrated using these parameters.
Benchmark quantities determine an anchor point for the set of indifference curves.
Benchmark prices fix the slope of the indifference curve at that point, and the
elasticity describes the curvature of the indifference curve.
The MPSGE representation includes a discount rate which is define implicitly
as:
1+r
ρ= −1
(1 + g)θ

Numeric implementation
The code on the following pages presents a GAMS/MPSGE model which has been
formulated following these ideas. Lines 1 to 72 reads base year data describing a
steady-state equilibrium. Investment levels are imputed from the base year capital
stock which is in-turn inferred from the assumed capital value share. Lines 73 to 150

41
 Figure 16: Calibrated intertemporal preferences

Ct+1 6

1+g

pt
pt+1
=1+r
-
1 Ct

declares the GAMS/MPSGE model, assigns steady-state values for activity levels
and price, and then checks consistency of the resulting model. Lines 151 to 157
runs a policy experiment. It assigns a tax on capital earnings beginning in year 6.
The resulting equilibrium is computed assuming that economic agents anticipate
the application of the tax, resulting in a sharp response in investment and other
economic variables to the new economic environment. Over time the tax leads to
a reduction in the steady-state capital stock and the real wage. Finally, lines 158
to the end show how to present output in graphs using GNUPLOT, both to the
windows screen and as encapsulated postscript.
1 $TITLE Ramsey Model - MPSGE formulation
2
3 $ontext
4
5 Calibrate to the steady-state condition:
6
7 I0 = KD0 * (g + delta) / (r + delta)
8
9 where g=2, delta=7, r=5, so
10
11 I0 = 48 * 9 / 12 = 36
12
13 Y I FD
14 P 100 -36 -64
15 PL -52 52
16 RK -48 48
17 PS 36 -36
18
19 $offtext
20
21 SET tt Time horizon (including the first year of the post-terminal period)
22 /2004*2081/,
23 t(tt) Time period over the model horizon
24 /2004*2080/;
25
26 SET t0(t), tl(t), tterm(tt);
27
28 PARAMETER g Growth rate /0.02/
29 r Interest rate /0.05/
30 delta Depreciation rate /0.07/
31 kvs Capital value share /0.48/
32 sigma Elasticity of substitution /1.00/
33
34
35 y0 Base year output

42
 36 kd0 Base year rental value of capital
37
38 k0 Base year capital stock
39 i0 Base year investment
40 c0 Base year consumption
41 l0 Base year labor input
42
43 kstock Base year capital stock multiplier /1/
44
45 taxk(t) Capital tax rate in period T
46
47 qref(t) Reference quantity path
48 pref(tt) Reference price path;
49
50
51 * Use the GAMS ORD (ordinality) and CARD (cardinality)
52 * functions to automate the identification of the first
53 * and last periods of the model horizon:
54
55 t0(t) = yes$(ord(t) eq 1);
56 tl(t) = yes$(ord(t) eq card(t));
57 tterm(tt) = yes$(ord(tt) eq card(tt));
58
59 * Calibrate the model to the baseline growth path:
60
61 y0 = 100;
62 kd0 = kvs * y0;
63 l0 = y0 - kd0;
64 k0 = kd0 / (r + delta);
65 i0 = (g + delta) * k0;
66 c0 = y0 - i0;
67 taxk(t) = 0;
68 qref(t) = (1+g)**(ord(t)-1);
69 pref(tt) = (1/(1+r))**(ord(tt)-1);
70
71 DISPLAY y0, kd0, l0, k0, i0, c0, g, r, delta;
72
73 $ONTEXT
74
75 $MODEL:RAMSEY
76
77 $SECTORS:
78 U ! Intertemporal utility index
79 Y(t) ! Output
80 I(t) ! Investment
81 K(t) ! Capital stock
82
83 $COMMODITIES:
84 PU ! Intertemporal utility price index
85 P(t) ! Output price
86 RK(t) ! Return to capital
87 PK(tt) ! Capital price
88 PL(t) ! Wage rate
89
90 $CONSUMERS:
91 RA ! Representative agent
92
93 $AUXILIARY:
94 TK ! Post-terminal capital stock
95
96 $PROD:Y(t) s:1
97 O:P(t) Q:Y0
98 I:PL(t) Q:L0
99 I:RK(t) Q:KD0 A:RA T:TaxK(t)
100
101 $PROD:K(tt)$T(tt)
102 O:PK(TT+1) Q:(K0*(1-DELTA))

43
103 O:RK(tt) Q:KD0
104 I:PK(tt) Q:K0
105
106 $PROD:I(tt)$T(tt)
107 O:PK(TT+1) Q:I0
108 I:P(tt) Q:I0
109
110 $PROD:U s:sigma
111 O:PU Q:(c0*sum(t, pref(t)*qref(t)))
112 I:P(t) Q:(qref(t)*c0) P:pref(t)
113
114 $DEMAND:RA
115 D:PU
116 E:PL(t) Q:(L0*qref(t))
117 E:PK(T0) Q:(K0*KSTOCK)
118 E:PK(TTERM) Q:-1 R:TK
119
120 $REPORT:
121 V:C(t) I:P(t) PROD:U
122 V:W W:RA
123
124 $CONSTRAINT:TK
125 SUM(T$TL(T+1), I(T+1)/I(t) - Y(T+1)/Y(t)) =E= 0;
126
127 $OFFTEXT
128 $SYSINCLUDE mpsgeset RAMSEY
129
130 * Assign steady-state equilibrium values for quantities and prices:
131
132 Y.L(t) = qref(t);
133 I.L(t) = qref(t);
134 K.L(t) = qref(t);
135
136 P.L(t) = pref(t);
137 RK.L(t) = pref(t);
138 PL.L(t) = pref(t);
139
140 * The steady-state price of capital is the output price
141 * times one plus the interest rate:
142
143 PK.L(tt) = (1+r) * pref(tt);
144 TK.L = k0 * (1+g)**card(t);
145
146 [Link] = 0;
147 $INCLUDE [Link]
148 SOLVE RAMSEY USING MCP;
149 [Link] = 1000;
150
151 * Apply a tax on capital inputs of 25% beginning in year 6:
152
153 TAXK(t)$(ORD(t) > 5) = 0.25;
154
155 $INCLUDE [Link]
156 SOLVE RAMSEY USING MCP;
157
158 * Generate some reports with graphs:
159
160 PARAMETER indices "Consumption, Investment and Capital Stock Indices";
161
162 indices(t,"C") = C.L(t)/(c0*qref(t));
163 indices(t,"I") = I.L(t)/qref(t);
164 indices(t,"K") = K.L(t)/qref(t);
165
166 DISPLAY INDICES;
167
168 * Define the domain over which the X-axis will be defined:
169

44
170 $setglobal domain t
171
172 * Define the labels to be printed along the X axis:
173
174 set tlbl(t) Time periods to be labelled in output plots /2010,2030,2050,2070/;
175
176 $setglobal labels tlbl
177
178 * Place the key to the figures outside the graph (see GNUPLOT 3.7 Help File):
179
180 $setglobal gp_opt0 "set key outside"
181
182 * Plot the graph with horizonal and vertical grid lines (see GNUPLOT 3.7 Help File):
183
184 $setglobal gp_opt1 "set grid"
185
186 * Generate the plot to the screen -- it can subsequently be copied to the clipboard
187 * using a right-click of the mouse, and then pasted into a separate program for
188 * publication:
189
190 $if %batch%==yes $setglobal batch yes
191 $if %batch%==yes $setglobal gp_opt1 "set term postscript eps monochrome ’Times-Roman’ 20"
192 $if %batch%==yes $setglobal gp_opt2 "set title"
193
194 $if %batch%==yes $setglobal gp_opt3 "set output ’[Link]’"
195 $setglobal gp_opt4 "set key outside width 4"
196
197 $libinclude plot indices
198
199 * Repeat the report generation process a couple more times:
200
201 parameter price "Capital prices and wage rate";
202 price(t,"RK") = RK.L(t)/P.L(t);
203 price(t,"PK") = PK.L(t)/((1+r)*P.L(t));
204 price(t,"PL") = PL.L(t)/P.L(t);
205
206 $if %batch%==yes $setglobal gp_opt3 "set output ’[Link]’"
207 $libinclude plot price
208
209 parameter grates(t,*) Growth rates through the transition;
210
211 grates(t,"c") = 100 * (C.l(t+1)/C.l(t) - 1);
212 grates(t,"i") = 100 * (I.l(t+1)/I.l(t) - 1);
213 grates(t,"k") = 100 * (K.l(t+1)/K.l(t) - 1);
214
215 set gs /c,i,k/;
216 grates(tl,gs) = na;
217
218 $if %batch%==yes $setglobal gp_opt3 "set output ’[Link]’"
219 $libinclude plot grates

45
Figure 17: Consumption, Investment and Capital Stock Indices
1.15
C
1.1 I
K
1.05

0.95

0.9

0.85

0.8

0.75

0.7

0.65

0.6
2010 2030 2050 2070

Figure 18: Capital Prices and Wage Rates

1.05
RK
PK
PL
1

0.95

0.9

0.85

0.8
2010 2030 2050 2070

46
 Figure 19: Growth Rates through the Transition
3
C
I
2 K

-1

-2

-3

-4

-5
2010 2030 2050 2070

References
Barro, R. J. and X. Sala-I-Martin (2004). Economic Growth. New York: McGraw-
Hill.
Easterly, W. (1993, November). How much do distortions affect growth? Journal
of Monetary Economics 32, 187–212.
Kamien, M. I. and N. L. Schwartz (2000). Dynamic Optimization. Amsterdam:
Elsevier Science.
Lau, M. I., A. Pahlke, and T. F. Rutherford (2002). Approximating infinite-
horizon models in a complementarity format: A primer in dynamic general
equilibrium analysis. Journal of Economic Dynamics & Control 26, 577–609.
Rutherford, T. F. (1999). Applied general equilibrium modeling with MPSGE
as a GAMS subsystem: An overview of the modeling framework and syntax.
Computational Economics 14, 1–46.

47