Hy,
I have some trouble solving a nash-cournot game using emp. I get a
solution for some values of a parameter I want to vary but not for all
values.
I have tried:
-Making sure the problem is well scaled (tricky part here is that I
get 0 if I scale too much).
-use initial values that are plausibly close to the solution.
-trying to solve with different initial values.
-setting the maximum iterations to a higher level.
Can someone give me other tips?
Here are my results… I put the maximum level of iterations to 11.
alpha is: 0.10–objfunc: -3.25–i invests in
coal: 0.00–i invests in CCGT: 0.00 --i invests
in OGGT: 2.15–found after 11.00–
alpha is: 0.20–objfunc:-7.71103E+10–i invests in coal:
18776.23–i invests in CCGT: 60161.16–i invests in OGGT:
53552.12–found after 11.00–
alpha is: 0.30–objfunc:-7.33490E+10–i invests in coal:
20823.27–i invests in CCGT: 59427.37–i invests in OGGT:
49994.43–found after 3.00–
alpha is: 0.40–objfunc:-7.20880E+10–i invests in coal:
20098.60–i invests in CCGT: 60791.08–i invests in OGGT:
48999.15–found after 11.00–
alpha is: 0.50–objfunc:-7.49596E+10–i invests in coal:
0.00–i invests in CCGT: 81393.48–i invests in OGGT:
51030.34–found after 3.00–
alpha is: 0.60–objfunc:-6.09310E+10–i invests in coal:
0.00–i invests in CCGT: 80871.89–i invests in OGGT:
39687.42–found after 11.00–
alpha is: 0.70–objfunc:-7.72888E+10–i invests in coal:
0.00–i invests in CCGT: 83264.56–i invests in OGGT:
51786.41–found after 1.00–
alpha is: 0.80–objfunc:-7.82225E+10–i invests in coal:
0.00–i invests in CCGT: 83747.55–i invests in OGGT:
52354.27–found after 1.00–
alpha is: 0.90–objfunc:-6.84628E+10–i invests in coal:
0.00–i invests in CCGT: 83375.87–i invests in OGGT:
45310.62–found after 11.00–
alpha is: 1.00–objfunc:-8.02378E+10–i invests in coal:
0.00–i invests in CCGT: 84957.44–i invests in OGGT:
53351.17–found after 1.00–
I have included my program below. If you run it look for a .dat file
in you workspace called investvariable.dat.
$title expansion capacity model in electricity market with risk averse
investors.
$ontext
this is a model with fixed demand.
$offtext
$eolcom //
set S /s1s8/
loopy /l1l10/
TIME time / t1, t2, t3, t4, t5, t6/
PLANT the type of plant /coal, CCGT, OGCT /
P possible events /p1*p10/
Risk the risk factors /load, fuel, CO2price/;
option limrow = 0, limcol = 0 ;
$offlisting
*the following parameters give initial values to the variables.
parameter uppercoal /100/;
parameter lowercoal /0/;
parameter upper /100/;
parameter lower /0/;
parameter step /1000/;
parameter perturbation/0/;
*this is the end of the initialization parameters
*This gives you the maximum amount that path tries to solve the
optimization problem for a specific alpha.
parameter interationlimit /10/;
*the initial probabilities.
parameter pbase(P) /
p1 1
p2 1
p3 1
p4 1
p5 1
p6 1
p7 1
p8 1
p9 1
p10 1
/;
\
- alpha determines the level of risk aversion. This program solves it
for different values of alpha.
parameter alpha /0/;
parameter totCost(P);
parameter counter;
alias(D,P);
pbase(P)= pbase(P)/sum(D,pbase(D));
- Q=A-B*p. A and B change with time.
Parameter A(TIME) the demand of electricity
/ t1 1.29
t2 1.245
t3 1.20
t4 0.90
t5 0.60
t6 0.30 /;
Parameter B(TIME) the demand of electricity
/ t1 283
t2 273
t3 613
t4 459
t5 320
t6 160 /;
*This is needed to obtain the inverse demand function
B(TIME)=1/(B(TIME));
A(TIME)=A(TIME)/(B(TIME));
Parameter tau(TIME) the duration of the demand
/ t1 10
t2 40
t3 310
t4 4400
t5 3000
t6 1000 /;
Parameter c(PLANT) The cost of operating a plant
/ coal 25
CCGT 45
OGGT 80 /;
Parameter I(PLANT) Investment cost to a plant to add one
additional capacity
/ coal 140000
CCGT 80000
OGGT 60000 /;
Parameter e(PLANT) The amount of polution emitted by a plant
/ coal 1
CCGT 0.35
OGGT 0.6
/;
parameter intermediatex /1/;
parameter intermediatey /3/;
parameter thepreviousworked /0/;
parameter previousx(PLANT) /coal 0, ogct 0/;
parameter previousy(PLANT,TIME,P);
parameter previousq(P);
previousq(P)=0;
previousy(PLANT,TIME,P)=0;
Scalar PC /300/ ;
Scalar iter /0/;
Scalar CO2 the price of emitting one tonne of CO2 /20/ ;
Table SCENARIO(P,RISK)
load Fuel CO2price
p1 .9 .7 0.5
p2 .9 1 0.5
p3 .9 1.3 0.5
p4 1 .7 0.5
p5 1 1 1.5
p6 1 1.3 1.5
p7 1.1 .7 1.5
p8 1.1 1 1.5
p9 1.1 1.3 1.5
p10 0.8 .7 1.5;
variables h1,h3;
positive variables x(PLANT),y(PLANT,TIME,P);
positive variables q(P);
equations defh1,
capacityconstraint(PLANT,TIME,P),priceconstraint(TIME,P);
equations defh3, a0, a1(P);
alias(planty,plant);
defh1… h1 =E=
sum(PLANT,I(PLANT)x(PLANT))+sum(P,q(P)( +sum(PLANT,sum( TIME,
( -( A(TIME)-B(TIME)*sum(PLANTY,y(PLANTY,TIME,P)) )+
c(PLANT)SCENARIO(P,“fuel”) +
CO2e(PLANT)*SCENARIO(P,“co2price”) )*tau(TIME)y(PLANT,TIME,P) ) ) ));
capacityconstraint(PLANT,TIME,P)… x(PLANT)-y(PLANT,TIME,P) =g=
0;
priceconstraint(TIME,P)… A(TIME)-
B(TIME)(sum(PLANT,y(PLANT,TIME,P)))=G=0;
- h3 is the objective function of the risk agent. a0 and a1(S)
guarantees that q is a probability measure that is in the risk set of
CVaRalpha.
defh3… h3 =e= sum(P,q(P)*( +sum(PLANT,sum( TIME,( -
( A(TIME)-B(TIME)*sum(PLANTY,y(PLANTY,TIME,P)) )+
c(PLANT)SCENARIO(P,“fuel”) +
CO2e(PLANT)*SCENARIO(P,“co2price”) )*tau(TIME)*y(PLANT,TIME,P) ) ) ));
a0… sum(P,q(P)) =E= 1;
a1(P)… -q(P)+pbase(P)/(alpha) =G= 0;
model nash /all/;
file myinfo / ‘%emp.info%’ /;
put myinfo ‘* Cournot-Nash equilibrium with two players and a least
risk averse player’;
put / ‘equilibrium’
put / ‘min h1 x y defh1 capacityconstraint priceconstraint’;
put / ‘max h3 q defh3 a0 a1’;
putclose;
file investments /investvariable.dat/;
alias(loopy1,loopy)
*Basically I solve for alpha=0.1…1. Every alpha I try to solve it a
number of times with different initial values and then I give up if it
hasn’t found a solution.
*It then moves to the next alpha.
loop(loopy1,
alpha=alpha+0.1;
iter=iter+1;
counter=0;
while(nash.modelstat2,
counter=counter+1;
lowercoal=lowercoal+step;
uppercoal=uppercoal+step;
lower=lower+step;
upper=upper+step;
q.l(P)= uniform(pbase(P)/alpha,pbase(P)/alpha);
while( intermediatey > intermediatex,
intermediatex=uniform(lowercoal,uppercoal);
intermediatey=uniform(lower,upper);
x.l(PLANT)=intermediatex;
y.l(PLANT,TIME,P)=intermediatey;
h1.l = sum(PLANT,I(PLANT)x.l(PLANT))
+sum(P,q.l(P)( +sum(PLANT,sum( TIME,( -( A(TIME)-
B(TIME)*sum(PLANTY,y.l(PLANTY,TIME,P)) )+
c(PLANT)*SCENARIO(P,“fuel”) )*tau(TIME)*y.l(PLANT,TIME,P) ) ) ));
)
if(thepreviousworked=1,
perturbation=uniform(0,1000);
x.l(PLANT)=previousx(PLANT)-
perturbation;
y.l(PLANT,TIME,P)=previousy(PLANT,TIME,P)-perturbation;
q.l(P)=(previousq(P)-0.0001)*(alpha
+0.1)/(alpha);
h1.l = sum(PLANT,I(PLANT)x.l(PLANT))
+sum(P,q.l(P)( +sum(PLANT,sum( TIME,( -( A(TIME)-
B(TIME)*sum(PLANTY,y.l(PLANTY,TIME,P)) )+
c(PLANT)*SCENARIO(P,“fuel”) )*tau(TIME)*y.l(PLANT,TIME,P) ) ) ));
)
option iterlim = 20000;
solve nash using emp;
if(nash.modelstat=2,
thepreviousworked=1;
previousx(PLANT)=x.l(PLANT);
previousy(PLANT,TIME,P)=y.l(PLANT,TIME,P);
previousq(P)=q.l(P);
)
if(counter>interationlimit,
nash.modelstat=2;
thepreviousworked=0;
)
)
*totCost(P)=sum(TIME,tau(TIME)*z.l(TIME,P)*PC)+sum(PLANT,sum( TIME,
(c(PLANT)SCENARIO(P,“fuel”)
+CO2e(PLANT))*tau(TIME)*y.l(PLANT,TIME,P) ) );
-
display totCost;
display x.l,h1.l;
display nash.MODELSTAT,nash.SOLVESTAT;
put investments;
put @1#iter, ‘alpha is:’, alpha;
put ‘–objfunc:’,h1.l,‘–’;
put ‘i invests in coal:’,x.l(“COAL”),‘–’;
put ‘i invests in CCGT:’,x.l(“CCGT”),‘–’;
put ‘i invests in OGCT:’,x.l(“OGCT”),‘–’;
put ‘found after ‘,counter,’–’;
nash.modelstat=0;
)
putclose
\