$ontext
DICE 2007 after canceling extra variables and equations and simplification
$offtext
SETS T Time periods /1*60/ ;
parameter deltat interval length of each period (number of years);
deltat = 10;
SCALARS
** Preferences
B_ELASMU Elasticity of marginal utility of consumption / 2.0 /
B_PRSTP Initial rate of social time preference per year / .015 /
** Population and technology
POP0 2005 world population millions /6514 /
GPOP0 Growth rate of population per year /.035 /
POPASYM Asymptotic population / 8600 /
A0 Initial level of total factor productivity /.02722 /
GA0 Initial growth rate for technology per year /.0092 /
DELA Decline rate of technol change per year /.001 /
DK Depreciation rate on capital per year /.100 /
GAMA Capital elasticity in production function /.300 /
K0 2005 value capital trill 2005 US dollars /137. /
** Emissions
SIG0 CO2-equivalent emissions-GNP ratio 2005 /.13418 /
GSIGMA Initial growth of sigma per year /-.00730 /
DSIG Decline rate of decarbonization per year /.003 /
ELAND0 Carbon emissions from land 2005(GtC per year) / 1.1000 /
** Carbon cycle
MAT2000 Concentration in atmosphere 2005 (GtC) /808.9 /
MU2000 Concentration in upper strata 2005 (GtC) /1255 /
ML2000 Concentration in lower strata 2005 (GtC) /18365 /
b12 Carbon cycle transition matrix per year /0.0189288 /
b23 Carbon cycle transition matrix per year /0.005 /
** Climate model
T2XCO2 Equilibrium temp impact of CO2 doubling oC / 3 /
FEX0 Estimate of 2000 forcings of non-CO2 GHG / -.06 /
FEX1 Estimate of 2100 forcings of non-CO2 GHG / 0.30 /
TOCEAN0 2000 lower strat. temp change (C) from 1900 /.0068 /
TATM0 2000 atmospheric temp change (C)from 1900 /.7307 /
C1 Climate-equation coefficient for upper level per year /.0220 /
C3 Transfer coeffic upper to lower stratum per year /.300 /
C4 Transfer coeffic for lower level per year /.0050 /
FCO22X Estimated forcings of equilibrium co2 doubling /3.8 /
** Climate damage parameters calibrated for quadratic at 2.5 C for 2105
A1 Damage intercept / 0.00000 /
A2 Damage quadratic term / 0.0028388 /
A3 Damage exponent / 2.00 /
** Abatement cost
EXPCOST2 Exponent of control cost function /2.8 /
PBACK Cost of backstop 2005 000$ per tC 2005 /1.17 /
BACKRAT Ratio initial to final backstop cost / 2 /
GBACK Initial cost decline backstop pc per year /.005 /
LIMMIU Upper limit on control rate / 1 /
** Participation
PARTFRACT1 Fraction of emissions under control regime 2005 /1 /
PARTFRACT2 Fraction of emissions under control regime 2015 /1 /
PARTFRACT21 Fraction of emissions under control regime 2205 /1 /
DPARTFRACT Decline rate of participation per year /0 /
** Availability of fossil fuels
FOSSLIM Maximum cumulative extraction fossil fuels / 6000 /
;
SET TLAST(T);
TLAST(T) = YES$(ORD(T) EQ CARD(T));
parameters
bb11 Carbon cycle transition matrix
bb12 Carbon cycle transition matrix
bb21 Carbon cycle transition matrix
bb22 Carbon cycle transition matrix
bb23 Carbon cycle transition matrix
bb32 Carbon cycle transition matrix
bb33 Carbon cycle transition matrix
CC1 Climate-equation coefficient for upper level
CC3 Transfer coeffic upper to lower stratum
CC4 Transfer coeffic for lower level
GGPOP0 Growth rate of population
GGA0 Initial growth rate for technology
DDELA Decline rate of technol change
GGSIGMA Initial growth of sigma
DDSIG Decline rate of decarbonization
EELAND0 Carbon emissions from land 2005(GtC)
GGBACK Initial cost decline backstop pc
DDPARTFRACT Decline rate of participation
scale1 Scaling coefficient in the objective function;
bb12 = b12 * deltat;
bb23 = b23 * deltat;
bb11 = 1 – bb12;
bb21 = 587.473*bb12/1143.894;
bb22 = 1 – bb21 – bb23;
bb32 = 1143.894*bb23/18340;
bb33 = 1 – bb32 ;
CC1 = C1*deltat;
CC4 = C4*deltat;
CC3 = C3;
GGPOP0 = GPOP0*deltat;
GGA0 = GA0*deltat;
DDELA = DELA*deltat;
GGSIGMA = GSIGMA*deltat;
DDSIG = DSIG*deltat;
EELAND0 = ELAND0*deltat;
GGBACK = GBACK*deltat;
DDPARTFRACT = DPARTFRACT*deltat;
scale1 = POPASYM*deltat;
PARAMETERS
L(T) Level of population and labor
AL(T) Level of total factor productivity
SIGMA(T) CO2-equivalent-emissions output ratio
RR(T) Average utility social discount factor
beta one-period discount factor
GA(T) Growth rate of productivity from 0 to T
FORCOTH(T) Exogenous forcing for other greenhouse gases
GCOST1 Growth of cost factor
GSIG(T) Cumulative improvement of energy efficiency
ETREE(T) Emissions from deforestation
COST1(t) Adjusted cost for backstop
PARTFRACT(T) Fraction of emissions in control regime
LAM Climate model parameter
Gfacpop(T) Growth factor population ;
* Important parameters for the model
LAM = FCO22X/ T2XCO2;
Gfacpop(T) = (exp(GGPOP0*(ORD(T)-1))-1)/exp(GGPOP0*(ORD(T)-1));
L(T)=POP0* (1- Gfacpop(T))+Gfacpop(T)*popasym;
ga(T)=GGA0*EXP(-DDELA*(ORD(T)-1));
al(“1”) = a0;
LOOP(T, al(T+1)=al(T)/((1-ga(T))););
gsig(T)=GGSIGMA*EXP(-DDSIG*(ORD(T)-1));sigma(“1”)=sig0;LOOP(T,sigma(T+1)=(sigma(T)/((1-gsig(T+1)))););
cost1(T) = (PBACK*SIGMA(T)/EXPCOST2)* ( (BACKRAT-1+ EXP (-GGBACK* (ORD(T)-1) ) )/BACKRAT);
ETREE(T) = EELAND0*(1-0.01*deltat)**(ord(T)-1);
RR(t)=1/((1+B_PRSTP)**(deltat*(ord(T)-1)));
beta = 1 / ((1+B_PRSTP)**deltat);
FORCOTH(T)= FEX0+ (FEX1-FEX0)/(120/deltat-2)*(ORD(T)-1)$(ORD(T) LT 120/deltat)+ 0.36$(ORD(T) GE 120/deltat);
partfract(t) = partfract21;
PARTFRACT(T)$(ord(T)<250/deltat) = Partfract21 + (PARTFRACT2-Partfract21)*exp(-DDPARTFRACT*(ORD(T)-2));
*partfract(“1”)= PARTFRACT1;
partfract(“1”)= 0.25372;
*****************************************************************
VARIABLES
K(T) Capital stock trillions US dollars
MAT(T) Carbon concentration in atmosphere GtC
MU(T) Carbon concentration in shallow oceans Gtc
ML(T) Carbon concentration in lower oceans GtC
TATM(T) Temperature of atmosphere in degrees C
TOCEAN(T) Temperature of lower oceans degrees C
E(T) CO2-equivalent emissions GtC
C(T) Consumption trillions US dollars
MIU(T) Emission control rate GHGs
CCA total industrial carbon emissions GTC
UTILITY;
POSITIVE VARIABLES K, MAT, MU, ML, TATM, TOCEAN, C, MIU, E;
EQUATIONS
KK0 Initial condition for capital
MMAT0 Starting atmospheric concentration
MMU0 Initial shallow ocean concentration
MML0 Initial lower ocean concentration
TATM0EQ Initial condition for atmospheric temperature
TOCEAN0EQ Initial condition for lower ocean temperature
KK(T) Capital balance equation
MMAT(T) Atmospheric concentration equation
MMU(T) Shallow ocean concentration
MML(T) Lower ocean concentration
TATMEQ(T) Temperature-climate equation for atmosphere
TOCEANEQ(T) Temperature-climate equation for lower oceans
EE(T) Emissions equation
FIXSEQ(T) Savings rate equation
CCACCA Cumulative carbon emissions
KendCond(T) terminal condition for capital
UTIL Objective function
;
** Equations of the model
KK0.. K(‘1’) =E= K0;
MMAT0.. MAT(‘1’) =E= MAT2000;
MMU0.. MU(‘1’) =E= MU2000;
MML0.. ML(‘1’) =E= ML2000;
TATM0EQ.. TATM(‘1’) =E= TATM0;
TOCEAN0EQ.. TOCEAN(‘1’) =E= TOCEAN0;
KK(T).. K(T+1) =L= (1-DK)**deltat *K(T)+deltat*( (1-((PARTFRACT(T)**(1-expcost2))*(cost1(t)*(MIU(T)**EXPcost2))))*
(AL(T)*L(T)**(1-GAMA)*K(T)**GAMA) / (1+A1*TATM(T)+ A2*TATM(T)**A3) – C(T) );
MMAT(T+1).. MAT(T+1) =E= MAT(T)*bb11+MU(T)*bb21 + E(T);
MML(T+1).. ML(T+1) =E= ML(T)*bb33+bb23*MU(T);
MMU(T+1).. MU(T+1) =E= MAT(T)*bb12+MU(T)*bb22+ML(T)*bb32;
TATMEQ(T+1).. TATM(T+1) =E= TATM(t)+CC1*(FCO22X*((log(((MAT(T+1)+MAT(T+2))/2+0.000001)/596.4)/log(2)))+FORCOTH(T) –
LAM*TATM(t)-CC3*(TATM(t)-TOCEAN(t)));
TOCEANEQ(T+1).. TOCEAN(T+1) =E= TOCEAN(T)+CC4*(TATM(T)-TOCEAN(T));
EE(T).. E(T) =E= deltat*SIGMA(T)*(1-MIU(T))* (AL(T)*L(T)**(1-GAMA)*K(T)**GAMA) + ETREE(T);
FIXSEQ(T).. (1- ((PARTFRACT(T)**(1-expcost2))*(cost1(t)*(MIU(T)**EXPcost2))))*
(AL(T)*L(T)**(1-GAMA)*K(T)**GAMA) / (1+A1*TATM(T)+ A2*TATM(T)**A3) – C(T) =E= 0.22*( .001 + (1-
((PARTFRACT(T)**(1-expcost2))*(cost1(t)*(MIU(T)**EXPcost2))))*
(AL(T)*L(T)**(1-GAMA)*K(T)**GAMA) /
(1+A1*TATM(T)+ A2*TATM(T)**A3) );
CCACCA.. CCA =E= sum(T$(ord(T)<card(T)), E(T));
KendCond(TLAST).. 0.02*K(TLAST) =L= (1-((PARTFRACT(TLAST)**(1-expcost2))*(cost1(TLAST)*(MIU(TLAST)**EXPcost2))))*
(AL(TLAST)*L(TLAST)**(1-GAMA)*K(TLAST)**GAMA) / (1+A1*TATM(TLAST)+ A2*TATM(TLAST)**A3) – C(TLAST);
UTIL.. UTILITY =E= sum( T, (((C(T)/L(T))**(1-B_ELASMU)-1)/(1-B_ELASMU)) * (RR(T)*L(T)*deltat/194) ) + 381800 ;
** Upper and Lower Bounds: General conditions for stability
K.lo(T) = 100;
MAT.lo(T) = 10;
MU.lo(t) = 100;
ML.lo(t) = 1000;
TOCEAN.up(T) = 20;
TOCEAN.lo(T) = -1;
TATM.lo(t) = 0;
TATM.up(t) = 20;
C.lo(T) = 20;
miu.up(t) = LIMMIU;
* First period predetermined by Kyoto Protocol
miu.fx(“1”) = 0.005;
** Cumulative limits on carbon use at 6000 GtC
CCA.up = FOSSLIM;
**********************************
** Solution options
option iterlim = 99900;
option reslim = 99999;
option solprint = on;
option limrow = 0;
option limcol = 0;
option nlp = conopt;
model CO2 /all/;
solve CO2 maximizing UTILITY using nlp ;
solve CO2 maximizing UTILITY using nlp ;
solve CO2 maximizing UTILITY using nlp ;
solve CO2 maximizing UTILITY using nlp ;
solve CO2 maximizing UTILITY using nlp ;
solve CO2 maximizing UTILITY using nlp ;
*****************************
Parameters
opt_tax(t)
opt_mcemis(t);
opt_tax(t)=-1*ee.m(t)*1000/(kk.m(t)+.00000000001);
opt_mcemis(t)= expcost2*cost1(t)*miu.l(t)**(expcost2-1)/sigma(t)*1000;
parameters
Kend terminal capital
MATend terminal carbon in atmosphere
MUend terminal carbon in upper ocean
MLend terminal carbon in lower ocean
TATMend terminal temperature in atmosphere
TOCEANend terminal temperature in ocean
;
Kend = sum(TLAST, (1-DK)**deltat *K.l(TLAST)+deltat*( (1- ((PARTFRACT(TLAST)**(1-expcost2))*
(cost1(TLAST)*(MIU.l(TLAST)**EXPcost2))))* (AL(TLAST)*L(TLAST)**(1-GAMA)*K.l(TLAST)**GAMA)/
(1+A1*TATM.l(TLAST)+ A2*TATM.l(TLAST)**A3) – C.l(TLAST)));
MATend = sum(TLAST, MAT.l(TLAST)*bb11+MU.l(TLAST)*bb21 + E.l(TLAST));
MUend = sum(TLAST, MAT.l(TLAST)*bb12+MU.l(TLAST)*bb22+ML.l(TLAST)*bb32);
MLend = sum(TLAST, ML.l(TLAST)*bb33+bb23*MU.l(TLAST));
TATMend = sum(TLAST, TATM.l(TLAST)+CC1*(FCO22X*((log((MAT.l(TLAST))/596.4)/log(2)))+FORCOTH(TLAST) – LAM*TATM.l(TLAST)-CC3*(TATM.l(TLAST)-TOCEAN.l(TLAST))));
TOCEANend = sum(TLAST, TOCEAN.l(TLAST)+CC4*(TATM.l(TLAST)-TOCEAN.l(TLAST)));
display Kend, MATend, MUend, MLend, TATMend, TOCEANend, opt_tax, opt_mcemis;