$title Putty-SemiPutty Technology
* Technology choice is based on an ex-ante isoquant which forms an envelope
* for ex-post substitution possibilities. That is, once a point on the ex-ante
* isoquant has been established, input-choices can be reoptimized in subsequent periods.
* In this GAMS program, I evaluate the magnitude of approximation errors
* which are introduced when installed technology is selected on the
* basis of weighted or unweighted average factor prices over the life of the capital stock.
set t Time horizon /0*100/;
variables
obj Variable profit -- present value return to a unit of capital
Y(t) Output
L(t) Employment
K(t) Capital stock
E(t) Energy use
Cbar Reference cost
Ebar Reference energy use
Lbar Reference labor use
Ybar Reference output
pKbar Reference price of capital
pLbar Reference price of labor
ALPHA Reference energy value share;
* Benchmark statistics:
parameter
l0 Baseline labor demand /0.45/
k0 Baseline capital stock
e0 Baseline energy price /0.05/
lamda Capital survival share /0.93/
r0 Interest rate /0.05/
sigmaL Long-run elasticity (dual) /0.5/
sigmaS Short-run energy elasticity (dual) /0.2/
sigmaKL Short-run KL elasticity (dual) /0.5/
Kbar Reference capital use
theta Baseline energy value share
rk0 Base year capital price
rhoL Long-run elasticity (primal)
rhoS Short-run energy elasticity (primal)
rhoKL Short-run elasticity (primal)
srvshr Survival share
pref Reference price
kvs Baseline capital value share
py(t) Projected output price
pl(t) Projected wage rate
pe(t) Projected energy price
pinv Price of new vintage capital;
* Assign benchmark parameter values:
rhoL = 1 - 1/sigmaL;
rhoS = 1 - 1/sigmaS;
rhoKL = 1 - 1/sigmaKL;
rk0 = r0+1-lamda;
k0 = (1-e0-l0)/rk0;
kvs = rk0*k0/(rk0*k0+l0);
theta = e0;
pref(t) = 1/(1+r0)**ord(t);
py(t) = pref(t);
pl(t) = pref(t);
pe(t) = pref(t);
pinv = 1;
srvshr(t) = lamda**(ord(t)-1);
* Fix the capital investment:
Kbar = k0;
equations objdef,output,capital,benchmark,alphadef,cbardef,lbardef,kbardef,ebardef;
* Objective function: maximize present value of returns:
objdef.. obj =e= sum(t, py(t)*Y(t) - pl(t)*L(t) - pe(t)*E(t))/Kbar;
* Period by period output is on the ex-post isoquant which is anchored at a reference point
* on the ex-ante isoquant:
output(t).. Y(t) =e= Ybar * ( ALPHA * (E(t)/Ebar)**rhos +
(1-ALPHA) * (kvs * (K(t)/kbar)**rhokl + (1-kvs) * (L(t)/Lbar)**rhokl)**(rhos/rhokl))**(1/rhos);
* Capital depreciates over time:
capital(t).. K(t) =e= Kbar * srvshr(t);
* Reference output level is on the ex-ante isoquant:
benchmark.. (theta * (Ebar/e0)**rhoL + (1-theta) * ((Kbar/k0)**kvs * (Lbar/L0)**(1-kvs))**rhoL)**(1/rhoL) =e= Ybar;
* Benchmark value share, scaling the price of energy to unity:
alphadef.. ALPHA * (Ebar + pLbar*Lbar + rk0*pKbar*Kbar) =e= Ebar;
* Benchmark cost (normalizing the price of energy to be unity):
cbardef.. Cbar =e= (theta + (1-theta)*(pKbar**kvs * pLbar**(1-kvs))**(1-sigmaL))**(1/(1-sigmaL));
Lbardef.. Lbar =e= Ybar * l0 * Cbar**sigmaL * (pKbar**kvs * pLbar**(1-kvs))**(1-sigmaL) / pLbar;
kbardef.. Kbar =e= Ybar * k0 * Cbar**sigmaL * (pKbar**kvs * pLbar**(1-kvs))**(1-sigmaL) / pKbar;
ebardef.. Ebar =e= Ybar * e0 * Cbar**sigmaL;
model puttyputty/all/;
* Assign some initial values:
Y.L(t) = srvshr(t);
L.L(t) = l0 * srvshr(t);
K.L(t) = k0 * srvshr(t);
E.L(t) = e0 * srvshr(t);
Cbar.L = 1;
Ybar.L = 1;
Ebar.L = e0;
Lbar.L = L0;
ALPHA.L = theta;
pLbar.L = 1;
pKbar.L = 1;
obj.L = sum(t, py(t)*Y.L(t) - pl(t)*L.L(t) - pe(t)*E.L(t))/Kbar;
* Install some bounds to avoid bad function calls:
Y.LO(t) = 0.01 * srvshr(t);
L.LO(t) = 0.01 * l0 * srvshr(t);
K.LO(t) = 0.01 * k0 * srvshr(t);
E.LO(t) = 0.01 * e0 * srvshr(t);
Cbar.LO = 0.01 * 1;
Ebar.LO = 0.01 * e0;
Lbar.LO = 0.01 * l0;
ALPHA.LO = 0.01 * theta;
pLbar.LO = 0.01 * 1;
pKbar.LO = 0.01 * 1;
solve puttyputty using nlp maximizing obj;
parameter intensity Energy intensity of output
profit Operating margin;
intensity(t,"Ref") = E.L(t)/(e0*Y.L(t));
profit(t,"Ref") = (py(t)*Y.L(t) - pl(t)*L.L(t) - pe(t)*E.L(t))/(srvshr(t)*pref(t));
* Assume that energy prices on average increase at 3% per year up to a
* factor of 5 above the steady-state level:
loop(t,
pe(t+1) = min(5*pe(t+1),pe(t)*1.03/(1+r0));
);
set tlbl(t) /0,20,40,60,80,100/;
$setglobal batch yes
$setglobal gp_opt3 "set key outside"
$setglobal gp_opt0 "set xlabel 'year'"
$if "%batch%"=="yes" $setglobal gp_opt1 "set term png"
$setglobal domain t
$setglobal labels tlbl
pe(t) = pe(t)/pref(t);
$if "%batch%"=="yes" $setglobal gp_opt2 "set output 'fig1.png'"
$libinclude plot pe
pe(t) = pe(t)*pref(t);
solve puttyputty using nlp maximizing obj;
profit(t,"Dynamic") = (py(t)*Y.L(t) - pl(t)*L.L(t) - pe(t)*E.L(t))/(srvshr(t)*pref(t));
intensity(t,"Dynamic") = E.L(t)/(e0*Y.L(t));
intensity(t,"Dyn_Ref") = Ebar.L/(e0*Ybar.L);
* Now do a calculation in which we use a weighted average of energy prices
* to select a point on the ex-ante isoquant:
pKbar.FX = sum(t,srvshr(t))/sum(t, srvshr(t)*pe(t)/pref(t));
pLbar.FX = sum(t,srvshr(t))/sum(t, srvshr(t)*pe(t)/pref(t));
* After having fixed the point, we then compute the time path of response
* to market prices:
solve puttyputty using nlp maximizing obj;
intensity(t,"Weighted_Ref") = Ebar.L/(e0*Ybar.L);
intensity(t,"Weighted") = E.L(t)/(e0*Y.L(t));
profit(t,"Weighted") = (py(t)*Y.L(t) - pl(t)*L.L(t) - pe(t)*E.L(t))/(srvshr(t)*pref(t));
pKbar.FX = card(t)/sum(t, pe(t)/pref(t));
pLbar.FX = card(t)/sum(t, pe(t)/pref(t));
solve puttyputty using nlp maximizing obj;
intensity(t,"Unweighted_Ref") = Ebar.L/(e0*Ybar.L);
intensity(t,"Unweighted") = E.L(t)/(e0*Y.L(t));
profit(t,"Unweighted") = (py(t)*Y.L(t) - pl(t)*L.L(t) - pe(t)*E.L(t))/(srvshr(t)*pref(t));
$if "%batch%"=="yes" $setglobal gp_opt2 "set output 'fig2.png'"
$libinclude plot intensity
$if "%batch%"=="yes" $setglobal gp_opt2 "set output 'fig3.png'"
$libinclude plot profit
* Repeat the experiment with highly variable energy prices:
pe(t) = uniform(0.5,1.5)*pe(t);
pe(t) = pe(t)/pref(t);
$if "%batch%"=="yes" $setglobal gp_opt2 "set output 'fig4.png'"
$libinclude plot pe
pe(t) = pe(t)*pref(t);
pLbar.LO = 0.01 * 1;
pKbar.LO = 0.01 * 1;
pKbar.UP = +INF;
pLbar.UP = +INF;
solve puttyputty using nlp maximizing obj;
profit(t,"Dynamic") = (py(t)*Y.L(t) - pl(t)*L.L(t) - pe(t)*E.L(t))/(srvshr(t)*pref(t));
intensity(t,"Dynamic") = E.L(t)/(e0*Y.L(t));
intensity(t,"Dyn_Ref") = Ebar.L/(e0*Ybar.L);
* Now do a calculation in which we use a weighted average of energy prices
* to select a point on the ex-ante isoquant:
pKbar.FX = sum(t,srvshr(t))/sum(t, srvshr(t)*pe(t)/pref(t));
pLbar.FX = sum(t,srvshr(t))/sum(t, srvshr(t)*pe(t)/pref(t));
* After having fixed the point, we then compute the time path of response
* to market prices:
solve puttyputty using nlp maximizing obj;
intensity(t,"Weighted_Ref") = Ebar.L/(e0*Ybar.L);
intensity(t,"Weighted") = E.L(t)/(e0*Y.L(t));
profit(t,"Weighted") = (py(t)*Y.L(t) - pl(t)*L.L(t) - pe(t)*E.L(t))/(srvshr(t)*pref(t));
pKbar.FX = card(t)/sum(t, pe(t)/pref(t));
pLbar.FX = card(t)/sum(t, pe(t)/pref(t));
solve puttyputty using nlp maximizing obj;
intensity(t,"Unweighted_Ref") = Ebar.L/(e0*Ybar.L);
intensity(t,"Unweighted") = E.L(t)/(e0*Y.L(t));
profit(t,"Unweighted") = (py(t)*Y.L(t) - pl(t)*L.L(t) - pe(t)*E.L(t))/(srvshr(t)*pref(t));
$if "%batch%"=="yes" $setglobal gp_opt2 "set output 'fig5.png'"
$libinclude plot intensity
$if "%batch%"=="yes" $setglobal gp_opt2 "set output 'fig6.png'"
$libinclude plot profit