This is my model. Equations condition_1 , condition_2, condition_3 accordingly refer to the bigM method. So,
- when my reference_actual(s) value is more than 1, binary variable u1(s)=1 and u2(s)=0 and u3(s)=0.
- when my reference_actual(s) value is between 0.97 and 1, binary variable u2(s) =1, one and u1(s)=0 and u3(s)=0.
- Otherwise,when my reference_actual(s) value is less than 0.97, binary variable u3(s) =1 and u2(s)=0 and u1(s)=0.
I don't know why but the model gives me sometimes even other values than zero or one for the binary variable.
$title Minimize_cost for random load profile
SETS
t #set of time slots
/1*10/
s #set of different scenaria for DSM activation
/1*10/
$setglobal path "C:\Users\Komnas\Desktop\";
*##data for load profile##########
*#################################
$call =xls2gms r=b2:c11 i=%path%case_1.xls o=%path%DEMAND.inc
* #Consumption of the producer KW
PARAMETER DEMAND(t) /
$include %path%DEMAND.inc
/
;
display DEMAND;
*##data for electricity price##########
*######################################
$call =xls2gms r=b2:b11;d2:d11 i=%path%case_1.xls o=%path%EL_PRICE.inc
* #Electricity price input data
PARAMETER EL_PRICE(t) /
$include %path%EL_PRICE.inc
/
;
display EL_PRICE;
*##data for different senaria####
*################################
$call =xls2gms r=b2:b11;f2:f11 i=%path%case_1.xls o=%path%A_SCENARIO.inc
* #Probability of Scenario to be examined
PARAMETER A_SCENARIO(s) /
$include %path%A_SCENARIO.inc
/
;
display A_SCENARIO;
*##data for activation of DSM####
*################################
TABLE
BINARY_ACTIVATION(t, s)
1 2 3 4 5 6 7 8 9 10
1 1 0 1 0 0 0 1 0 0 0
2 1 1 0 0 0 0 0 0 0 0
3 0 0 0 1 0 1 0 0 0 0
4 0 0 0 0 1 0 0 0 0 0
5 0 0 1 0 0 0 1 0 1 0
6 1 0 0 0 0 0 1 1 0 0
7 0 0 0 0 0 0 0 0 0 0
8 0 0 0 1 1 0 0 0 0 0
9 1 1 0 1 0 0 0 0 1 0
10 1 0 0 0 0 1 1 1 0 0
;
display BINARY_ACTIVATION ;
PARAMETERS
CONSUMPTION_PER_UNIT #Electricity consumption per 1 unit tonne production of chlorine /3.4/
PRODUCT_PRICE #The price that industry shell their product /850/
UNSHEDABLE_LOAD #The capacity that cannot be offered to grid MW /30/
PENALTY_PRICE #Price that producer has to pay for penalty /100/
REMUNERATION_RESERVATION #The price that ELIA remunarates producer for offering flexibility /1.41/
REMUNERATION_ACTIVATION #The price that ELIA remunarates producer for activation of DSM €\MWh /78/
HOURS_ACTIVE_DSM #Total time in Hours that DSM is activated- defined in contract /3/
TIME #Total time of contracted period
SHEDDIND_LIMIT #Marginal level real demand under this level when DSM activated /60/
ACTIV
BIG /15/
;
TIME = card(t);
ACTIV = smax(t, DEMAND(t));
POSITIVE VARIABLES
dsm(t,s) #Demand response for each time interval
real_load(t,s) #Actual load consumed demand response is applied
units(t,s) #Units produced per time interval by selling the product
total_DSM #Total activation of DSM during the whole period
revenue_sales(t,s) #Revenue the T time slot earning by selling each unit
reservetion_advance_revenue #Revenue come from DSM reservation price €
activ_remun_revenue(t,s) #Revenue come from each MW activated as DSM
contract_volume
advance
reference_actual(s)
extra_remuneration(s)
extra_penal_1(s)
extra_penal_2(s)
activation(s)
activ2(s)
dsmreal(s)
;
BINARY VARIABLE
u1
u2
u3
;
VARIABLES
cost #Total cost of industry from production procces and DSM activation
divergence_contract(s) #The difference between the signed contract and the sumation of DSM all over the period
payback_or_revenue(s) #Extra remuneration or penalty according to the reference per actual provision of reserve
;
*************
**Equations**
*************
EQUATIONS
ob_funct #Objective function (maximize profit for producer)
real_demand(t,s) #Calculation of the actual consumption (after modification of DSM)
n_units(t,s) #Number of units (Depends on the hourly actual consumption)
revenues(t,s) #Calculation of revenue the T time slot came from sell of each unit
least_unshedable_load(t,s) #The unshedable load level of industry
advance_calc #Calculation of reservation revenue advance that elia pays to the customer for flexibility
activation_remuneration(t,s) #Calculation of activation revenue
diverg_contract(s) #Residue from the agreed contract volume
condition(t,s) #Condition to ensure that product respects the limits set by ELIA = 5 MW
contract_volume_max
interruptible_DSM(t,s)
ref_per_act(s) #Reference per actual provision of reserve
* extra_remun_calc(s)
* extra_pen_calc_1(s)
* extra_pen_calc_2(s)
* payback_revenue(s) #Calculation of extra remuneration or penalty
condition_1(s) #condition to aply the bigM conditions
condition_2(s)
condition_3(s)
* condition_5(s)
;
*#################################
*######problem formulation########
*#################################
ob_funct..
cost
=e=
sum(s,
A_SCENARIO(s) * (
+ sum(t, real_load(t,s) * EL_PRICE(t))
- sum(t, revenue_sales(t,s))
- sum(t, activ_remun_revenue(t,s))
- reservetion_advance_revenue
* + payback_or_revenue(s))
)
);
contract_volume_max..
contract_volume =l= sum(t, DEMAND(t))/TIME - SHEDDIND_LIMIT ;
* contract_volume =l= smax(t, DEMAND(t));
*###########################
*####CONSTRAINTS############
*###########################
*#To ensure that product respects the limits
condition(t,s)..
dsm(t,s) =l= ACTIV - UNSHEDABLE_LOAD;
* dsm(t,s)
* =l= dsm(t,s) + BIG * (1 - BINARY_ACTIVATION(t,s));
*#In this ICH product target is fixed the shedding limit
interruptible_DSM(t,s)..
dsm(t,s)
=g=
(DEMAND(t) - SHEDDIND_LIMIT)* BINARY_ACTIVATION(t,s);
*#Actual Offtake energy from grid [MW]
real_demand(t,s)..
real_load(t,s) =e= DEMAND(t) - dsm(t,s)* BINARY_ACTIVATION(t,s) ;
*#Number of units produced per time slot
n_units(t,s)..
units(t,s) =e= real_load(t,s) / CONSUMPTION_PER_UNIT ;
*#Calculation of revenue the T time slot, came from sell of each unit
revenues(t,s)..
revenue_sales(t,s) =e= units(t,s) * PRODUCT_PRICE ;
*#To ensure that we don't overpass the unshedable load level of industry
least_unshedable_load(t,s)..
real_load(t,s) =g= UNSHEDABLE_LOAD ;
*#Calculation of activation remuneration revenue in €
activation_remuneration(t,s)..
activ_remun_revenue(t,s)
=e= dsm(t,s) * BINARY_ACTIVATION(t,s) * REMUNERATION_ACTIVATION ;
*#Calculation of reservation revenue,advances
advance_calc..
reservetion_advance_revenue
=e= TIME * REMUNERATION_RESERVATION * contract_volume ;
* hours €/ÌW/ h MW
*#calculation if there is a residue from the agreed contract volume
diverg_contract(s)..
divergence_contract(s)
=e= sum(t, BINARY_ACTIVATION(t,s) * (contract_volume - dsm(t,s)) ) ;
*#reference per actual provision of reserve
ref_per_act(s)..
reference_actual(s)
=e= sum(t, BINARY_ACTIVATION(t,s) * dsm(t,s)) /
( sum(t,BINARY_ACTIVATION(t,s) * contract_volume) + 0.000001 ) ;
*#we use these binary variables to take the different cases of reference per actual provision of reserve
*condition_1(s)..
u1(s) + u2(s) + u3(s) =e= 1 ;
*#if the average power reserve actually provided is between 97-100% of the reference power reserve
condition_2(s)..
reference_actual(s) =l= 0.97 + BIG * u2(s) ;
*#if the average power reserve actually provided is more than 100% of the reference power reserve
*#should be changed
condition_3(s)..
reference_actual(s) =l= 1 + BIG * u1(s) ;
*extra_remun_calc(s)..
* extra_remuneration(s)
* =e=
* u1(s) * ( - divergence_contract(s) * 1.2 * REMUNERATION_RESERVATION ) ;
*extra_pen_calc_1(s)..
* extra_penal_1(s)
* =e=
* u2(s) * reservetion_advance_revenue ;
*extra_pen_calc_2(s)..
* extra_penal_2(s)
* =e=
* u3(s) * (
* + divergence_contract(s) * 1.2 * REMUNERATION_RESERVATION
* + reservetion_advance_revenue ) ;
*#calculiation of extra remuneration or penalty to be paid in ELIA
*payback_revenue(s)..
* payback_or_revenue(s)
* =e=
* + extra_remuneration(s)
* - extra_penal_1(s)
* - extra_penal_2(s) ;
;
MODEL Minimization_cost /all/;
SOLVE Minimization_cost USING minlp minimization cost ;