$TITLE DIAM with recursive utility
$OFFUPPER
*
* This GAMS model by Minh Ha Duong, CIRED and Nicolas Treich, CIRANO-INRA
* Correspondence to : haduong@centre-cired.fr
*
* Reference: Minh Ha Duong and Nicolas Treich (1999),
*  Uncertainty and intergenerational equity in the global warming problem :
*  A recursive preference approach, FEEM Nota di Lavoro 92-99
*  Available online at           http://www.feem.it/web/resun/wp/92-99.html
*  (last access 2000/12/05)
*
* This is version: $Id: ultral.gms,v 1.28 2003/01/25 06:48:43 haduong Exp $
*
*    This program is free software; you can redistribute it and/or modify
*    it under the terms of the GNU General Public License as published by
*    the Free Software Foundation; either version 2 of the License, or
*    (at your option) any later version.
*

SETS

T            Periods                 /1*4/
I            Climate change impact   /Low,Big/

SCALARS

Dura         Period Length in years                       / 30.0  /
RES          Inter-period resistance to substitution      /  2.0  /
RRA          Relative Risk-Aversion                       /  1.5  /
Disc         Inter period utility Discount factor         /  0.74 /
DepCap       Inter period capital depreciation            /  0.215/
InvAccu      Inter period investment accumulation         / 16.14 /
SCap         Share of Capital in Cobb Douglas prod fun    /  0.36 /
SEne         Share of Energy                              /  0.04 /
Y1           Production 1990 no damage (k dol.90 percap)  /  3.82 /
E1           Energy 1990 no X (toe per year per capita)   /  1.21 /
K1           Initial capital (k dol.90 per capita)        / 11.46 /
Proba        Probability of Big                           /  0.25 / 
TSplit       Period of strategic choice                   /  3    / 

PARAMETERS
*SRESA2
A(T)         Productivity progress /1 1.00 , 2 1.22 , 3 1.55 , 4 1.62 /
P(T)         Relative energy cost  /1 1.00 , 2 1.07 , 3 1.33 , 4 1.79 /
L(T)         Population (Giga)     /1 5.26 , 2 8.19 , 3 11.30, 4 12.35 /
S(T)         Saving ratio          /1 0.30 , 2 0.30 , 3 0.25 , 4 0.25 /
CContent(T)  (tC per toe)          /1 0.97 , 2 0.91 , 3 0.84 , 4 1.17 /
OOrigin(T)   (ppmv per GtC)        /1  371 , 2  390 , 3  297 , 4 337 /
OSlope(T)    (ppmv)                /1 0.154, 2 0.129, 3 0.278, 4 0.168/

Impact(I)    Productivity loss     /Low  0.08, Big  0.30/

VARIABLES

X(I,T)       Energy cost multiplier (control variable)
Y(I,T)       Production (PIB per capita per year in k dol. 90)
K(I,T)       Capital stock per capita (k dol. 90)
E(I,T)       Energy (toe per capita per year)
C(I,T)       Consumption per capita per year (k dol. 90)
M(I,T)       End-of-period cuMulative carbon emissions (global GtC)
O(I,T)       End-of-period CO2 mixing ratio (ppmv)
F(I,T)       Climate damage (a multiplier in Y production function)
U(I,T)       Per capita cardinal Utility (k. dol. 90)
V(I,T)       Certainty equiValent of nextgen utility (k dol. 90 p.c.)
GOAL;

* We defining a range and an initial value for all variables to help solver.
* Model bounds (real constraints): X.LO, F.FX
* Algorithmic bounds, necessary  : E-K-U-V.LO, X.UP (arguments in power function) 
* Algorithmic constraints, useful: all other.
* Provide a sensible starting point hopefully can't hurt.
* No explicit scaling: units are so that hopefully the model scales well.
* In particular, variables F and Z should be around unity.
*
X.LO(I,T) = 1    ;  X.UP(I,T) = 4      ;  X.L(I,T) = 1;
Y.LO(I,T) = Y1/2 ;  Y.UP(I,T) = 30 * Y1;  Y.L(I,T) = Y1 * 1.01**(30*(ORD(T)-1));
K.LO(I,T) = K1   ;  K.UP(I,T) = 20 * K1;  K.L(I,T) = K1 / Y1 * Y.L(I,T);
E.LO(I,T) = E1/10;  E.UP(I,T) = 20 * E1;  E.L(I,T) = SEne / P(T) * Y.L(I,T);
C.LO(I,T) = Y1/10;  C.UP(I,T) = 20 * Y1;  C.L(I,T) = (1-S(T))* Y.L(I,T);
M.LO(I,T) = 100  ;  M.UP(I,T) = 1500   ;  M.L(I,T) = 1 + 200 * (ORD(T)-1);
O.LO(I,T) = 354  ;  O.UP(I,T) = 750    ;  O.L(I,T) = 354 + 40*ORD(T);
F.LO(I,T) = 0.5  ;  F.UP(I,T) = 1.5    ;  F.L(I,T) = 1;
U.LO(I,T) = 1    ;  U.UP(I,T) = 200    ;  U.L(I,T) = Y.L(I,T);
V.LO(I,T) = 1    ;  V.UP(I,T) = 200    ;  V.L(I,T) = Y.L(I,T);

* Initial conditions, should always appear after other bounds assignments
*
K.FX(I,'1') = K1;

EQUATIONS

EQY(I,T)     Cobb-Douglas production function
EQC(I,T)     Consumption at constant saving rate (net of adjustment costs)
EQE(I,T)     Constant share of income to energy sector
EQK(I,T)     Capital accumulation
EQM(I,T)     Emission accumulation
EQO(I,T)     Carbon cycle model
EQF(I,T)     Climate Change Impact
EQU(I,T)     Recursive CES felicity
EQV(I,T)     Expected CRRA next generation utility
EQFINALV(I)  Final condition on V (steady state)
EQFINALO(I)  Final condition on O (concentration stabilisation)
EQUNX(T)     Control unicity before TSplit
EQGOAL       Objective;

* Variables are per-capita, so the coefficients in CD function add up to 1-SLabor.
* Productivity progress is relative to 1990 : A(1) = 1.
* Assuming investment rate is preserved, capital installation costs are taken off consumption.
*
EQY(I,T)..       Y(I,T)/Y1 =E= A(T) * F(I,T) * (K(I,T)/K1)**SCap * (E(I,T)/E1)**SEne;

EQC(I,T)..          C(I,T) =E= (1 - S(T)) * Y(I,T);

* EQE stems implicitly from a comportmental assumption:
* At each period, the energy sector maximizes its profit myopically.
* Here the price of output is 1, the price of energy is P X.
* Then a well-known property of Cobb-Douglas production function is that
* the share of energy sector's income in production is its exponent (SEne)
* Consequently, at any period Y / (P X E) = SEne, therefore:

EQE(I,T)..       Y(I,T)/Y1 =E= P(T) * X(I,T) * E(I,T)/E1;

* GAMS specific: Indexing equations by T+1 generates only Card(T)-1 equations :
* Because EQK was declared as (I,T), GAMS will not assign EQK(I,5) (which is what we want).

EQK(I,T+1)..      K(I,T+1) =E= (DepCap*K(I,T) + InvAccu*S(T)*Y(I,T)) * L(T)/L(T+1);

* Using cumultative emissions M, a simple index used by IPCC.
* Note:  M is since 1990, not since preindustrial era.
* GAMS knows that M(I, 0) should be zero

EQM(I,T)..      M(I,T) =E= M(I,T-1) + Dura * L(T) * CContent(T) * E(I,T);

* OOrigin and OSlope were obtained by interpolating
* on seven major carbon cycle models, see file ccycle.xls.

EQO(I,T)..      O(I,T) =E= OOrigin(T) + OSlope(T) * M(I,T);

* The stochastic damage
* Assuming well known that 275 ppmv is the pre-industrial concentration level.

EQF(I,T)..      F(I,T) =E= 1 - Impact(I) * ((O(I,T)-275)/275) ;

  
* The recursive utility function.
* It seems harder to solve with U(I,T)**(1-RES) =E= when RES is large

EQU(I,T).. U(I,T)**(1-RES) =E= (1 - Disc) * (L(T)/L('4'))**RES * C(I,T)**(1-RES) + Disc * V(I,T)**(1-RES);

*EQU(I,T).. U(I,T) =E= ((1 - Disc) * (L(T)/L('4'))**RES * C(I,T)**(1-RES) + Disc * V(I,T)**(1-RES))**(1/(1-RES));

* Code for V and UNX is specific to binary uncertainty resolved at once:
* In clear:
* If T <> TSplit, then there is no uncertainty, consequently L(T) V(T) = U(T+1) L(T+1)
* If T =  TSplit, then we take the CRRA expected value of the above equation

EQV(I,T+1)..
             L(T) * V(I,T) =E= 
     (L(T+1) * U(I,T+1)
     )$(ORD(T) <> TSplit) +
     (( Proba    *(L(T+1)*U('Big',T+1))**(1-RRA) +
       (1-Proba) *(L(T+1)*U('Low' ,T+1))**(1-RRA)   ) ** (1/(1-RRA))
     )$(ORD(T) = TSplit);

* The 5th generation consumption is approximated as C(5) = C(4) * K(5)/K(4)
* Where K(5) is computed as in EQK.

EQFINALV(I)..     V(I,'4') =E= C(I,'4') * (DepCap*K(I,'4') + InvAccu*S('4')*Y(I,'4')) / K(I,'4');

EQUNX(T)$(ORD(T) LT TSplit+1)..
             X('Big',T) =E= X('Low',T);

* Last generation stabilizes CO2 concentration
*
EQFINALO(I)..   O(I,'4') =E= O(I,'3');

* Assertion: U('Low','1') equals U('Big','1').
*
EQGOAL..            GOAL =E= U('Low', '1');

MODEL ULTRAL /ALL/;

ULTRAL.optfile=1;
OPTION NLP=CONOPT;

SOLVE ULTRAL MAXIMIZING GOAL USING NLP