# Infeasibilities with economic dispatch model

**URL:** <https://discourse.julialang.org/t/infeasibilities-with-economic-dispatch-model/82571>\
**Category:** Optimization (Mathematical)\
**Tags:** question, jump, optimization, dispatch\
**Created:** [June 10, 2022, 4:14pm UTC](https://discourse.julialang.org/t/infeasibilities-with-economic-dispatch-model/82571 "2022-06-10T16:14:44Z")\
**Posts on this page:** 5\
**Page:** 1

<div class="post-metadata">

**Author:** ![mbah](https://avatars.discourse-cdn.com/v4/letter/m/bcef8e/32.png) [@mbah](https://discourse.julialang.org/u/mbah)\
**Post date:** [June 10, 2022, 4:14pm UTC](https://discourse.julialang.org/t/infeasibilities-with-economic-dispatch-model/82571/1 "2022-06-10T16:14:44Z")

</div>

Hello,

I am currently working on a network economic dispatch model. The basic transport model below works with only conventional generators. However, whenever I add renewable injection constraints with hourly capacity factors, the model cannot solve and I get an infeasibility message. I have tried various options to no avail. Any assistance would be highly appreciated.

```julia
#using DataFrame, JuMP

generators = DataFrame(index = ["conv_gas", "conv_nuc", "conv_oil", "ror_plant"], 
                      node = ["zone_1", "zone_1", "zone_2", "zone_2"],
                      R_ID = 1:4, node_id = [1,1,2,2], 
                                g_max = [100,250,150,60], cost = [50,10,20,0], conventional = [1,1,1,0], 
                                ROR = [0,0,0,1])

demand = DataFrame(load_z1 = [100,150,80], load_z2 = [150,165,170])

variability = DataFrame(conv_gas = [1,1,1], conv_nuc = [1,1,1], conv_oil = [1,1,1], ror_plant = [0.6,0.6,0.6])

network_zones = ["z1","z2"]

lines = DataFrame(network_lines = 1, z1 = 1, z2 = -1)

```

```julia
# Sets
G = generators.R_ID #all generators
T = convert(Array{Int64}, 1:3)
Z = convert(Array{Int64}, 1:2) #zones
L = lines.network_lines 

#SUBSETS

#C = generators.R_ID[generators.conventional.==1] #conventional generators
#ROR = intersect(generators.R_ID[generators.ROR.==1], G) #Run-Of-River

```

```julia
#Dispatch_Model = Model(Clp.Optimizer)

#@variable(Dispatch_Model, vGEN[T, c in C] >=0) 
#@variable(Dispatch_Model, vROR[T, r in ROR] >=0)
#@variable(Dispatch_Model, vFLOW[T,L])

```

```julia
#Constraints

#@constraint(Dispatch_Model, cMaxPower[t in T, c in C], 
# vGEN[t,c] <= generators[generators.R_ID .== c,:g_max][1])

@constraint(Dispatch_Model, cRORInjection[t in T, r in ROR], vROR[t,r] <= generators[r,:g_max] * variability[t,r])

#@expression(Dispatch_Model, eCosts, sum(vGEN[t,c]*generators.cost[c] for t in T, c in C))
#@objective(Dispatch_Model, Min, eCosts)

@constraint(Dispatch_Model, cEnergyBalance[t in T, z in Z],
            sum(vGEN[t,c] for c in intersect(generators[generators.node.==z, :R_ID],G))+
            sum(vRUNOFRIVER[t,r] for r in intersect(generators[generators.node.==z, :R_ID],G))
            -demand[t,z] -
            sum(lines[l,Symbol(string("z",z))] * vFLOW[t,l] for l in L) == 0
            )

```

```julia
@time optimize!(Dispatch_Model)

Model output message:

Coin0507I Presolve determined that the problem was infeasible with tolerance of 1e-08
Clp3003W Analysis indicates model infeasible or unbounded
Clp0006I 0 Obj 0 Primal inf 815 (5)
Clp0001I Primal infeasible - objective value 0
Clp0032I PrimalInfeasible objective 0 - 1 iterations time 0.012
  0.011694 seconds (2 allocations: 32 bytes)

Status
  Termination status : INFEASIBLE
  Primal status : UNKNOWN_RESULT_STATUS
  Dual status : INFEASIBILITY_CERTIFICATE
  Message from the solver:
  "1 - primal infeasible"

* Candidate solution
  Objective value : 0.0

* Work counters
  Solve time (sec) : 0.01200

```

Thanks

---

<div class="post-metadata">

**Author:** ![odow](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/odow/32/28685_2.png) [@odow](https://discourse.julialang.org/u/odow)\
**Post date:** [June 12, 2022, 7:11am UTC](https://discourse.julialang.org/t/infeasibilities-with-economic-dispatch-model/82571/2 "2022-06-12T07:11:31Z")

</div>

Do you expect there to always be a feasible solution? I can’t run your problem, because I don’t know which lines are meant to be commented or not.

Some tips to work on this:

- Simplify the problem. Remove constraints until it gets feasible
- Once you know which constraint causes the infeasibility, check the data. Is there a problem with the units? Did you get the sign of each term correct? It’s easy to have a `-` instead of a `+`.
- Make the problem smaller. Simplify the data, use less generators and less time periods
- Try changing the constraint. If it’s an equality, relax into `>=`, or add a slack variable that you penalize in the objective. What is the optimal value of the slack?

---

<div class="post-metadata">

**Author:** ![mbah](https://avatars.discourse-cdn.com/v4/letter/m/bcef8e/32.png) [@mbah](https://discourse.julialang.org/u/mbah)\
**Post date:** [June 14, 2022, 7:34am UTC](https://discourse.julialang.org/t/infeasibilities-with-economic-dispatch-model/82571/3 "2022-06-14T07:34:38Z")

</div>

Thanks. I re-posted the uncommented codes below. Only one of the constraints is causing the infeasibility ( **cRORInjection** ). The model is able to generate a feasible solution with the generator capacity constraint (cMaxPower) and the energy balance (cEnergyBalance) activated.

```julia
#using DataFrame, JuMP

generators = DataFrame(index = ["conv_gas", "conv_nuc", "conv_oil", "ror_plant"], 
                      node = ["zone_1", "zone_1", "zone_2", "zone_2"],
                      R_ID = 1:4, node_id = [1,1,2,2], 
                                g_max = [100,250,150,60], cost = [50,10,20,0], conventional = [1,1,1,0], 
                                ROR = [0,0,0,1])

demand = DataFrame(load_z1 = [100,150,80], load_z2 = [150,165,170])

variability = DataFrame(conv_gas = [1,1,1], conv_nuc = [1,1,1], conv_oil = [1,1,1], ror_plant = [0.6,0.6,0.6])

network_zones = ["z1","z2"]

lines = DataFrame(network_lines = 1, z1 = 1, z2 = -1)

```

```julia
# Sets
G = generators.R_ID #all generators
T = convert(Array{Int64}, 1:3)
Z = convert(Array{Int64}, 1:2) #zones
L = lines.network_lines 

#SUBSETS

C = generators.R_ID[generators.conventional.==1] #conventional generators
ROR = intersect(generators.R_ID[generators.ROR.==1], G) #Run-Of-River

```

```julia
Dispatch_Model = Model()

@variable(Dispatch_Model, vGEN[T, c in C] >=0) 
@variable(Dispatch_Model, vROR[T, r in ROR] >=0)
@variable(Dispatch_Model, vFLOW[T,L])

```

```julia
#Constraints

@constraint(Dispatch_Model, cMaxPower[t in T, c in C], 
       vGEN[t,c] <= generators[generators.R_ID .== c,:g_max][1])

@constraint(Dispatch_Model, cRORInjection[t in T, r in ROR], vROR[t,r] <= generators[r,:g_max] * variability[t,r])

@expression(Dispatch_Model, eCosts, sum(vGEN[t,c]*generators.cost[c] for t in T, c in C))
@objective(Dispatch_Model, Min, eCosts)

@constraint(Dispatch_Model, cEnergyBalance[t in T, z in Z],
            sum(vGEN[t,c] for c in intersect(generators[generators.node.==z, :R_ID],G))+
            sum(vRUNOFRIVER[t,r] for r in intersect(generators[generators.node.==z, :R_ID],G))
            -demand[t,z] -
            sum(lines[l,Symbol(string("z",z))] * vFLOW[t,l] for l in L) == 0
            )

@time optimize!(Dispatch_Model)

```

```julia
#Infeasibility results
Coin0507I Presolve determined that the problem was infeasible with tolerance of 1e-08
Clp3003W Analysis indicates model infeasible or unbounded
Clp0006I 0 Obj 0 Primal inf 815 (6)
Clp0001I Primal infeasible - objective value 0
Clp0032I PrimalInfeasible objective 0 - 1 iterations time 0.012
  7.858292 seconds (18.65 M allocations: 1.102 GiB, 4.93% gc time, 1.36% compilation time)

```

With regards to the slack variables, I have a couple of quick questions.

1. When adding the slack variable into the constraint, do i convert the constraint to an equality, like this?

```julia
@variable(Dispatch_Model, a)

@constraint(Dispatch_Model, cRORInjection[t in T, r in ROR], vROR[t,r] + a == generators[r,:g_max] * variability[t,r])

```

1. Should i use a unique slack variable for each constraint?

2. Could you direct me to a MWE for penalizing the slack variables in my objective function. I’ve tried doing it, but getting error messages.

Thanks!

---

<div class="post-metadata">

**Author:** ![odow](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/odow/32/28685_2.png) [@odow](https://discourse.julialang.org/u/odow)\
**Post date:** [June 14, 2022, 9:28am UTC](https://discourse.julialang.org/t/infeasibilities-with-economic-dispatch-model/82571/4 "2022-06-14T09:28:46Z")

</div>

If I print your energy balance constraint:

```julia
julia> cEnergyBalance
2-dimensional DenseAxisArray{ConstraintRef{Model, MathOptInterface.ConstraintIndex{MathOptInterface.ScalarAffineFunction{Float64}, MathOptInterface.EqualTo{Float64}}, ScalarShape},2,...} with index sets:
    Dimension 1, [1, 2, 3]
    Dimension 2, [1, 2]
And data, a 3×2 Matrix{ConstraintRef{Model, MathOptInterface.ConstraintIndex{MathOptInterface.ScalarAffineFunction{Float64}, MathOptInterface.EqualTo{Float64}}, ScalarShape}}:
 cEnergyBalance[1,1] : -vFLOW[1,1] = 100.0 cEnergyBalance[1,2] : vFLOW[1,1] = 150.0
 cEnergyBalance[2,1] : -vFLOW[2,1] = 150.0 cEnergyBalance[2,2] : vFLOW[2,1] = 165.0
 cEnergyBalance[3,1] : -vFLOW[3,1] = 80.0 cEnergyBalance[3,2] : vFLOW[3,1] = 170.0

```

`cEnergyBalance[1,1]` says that `vFLOW[1,1]` has to be `-100`, but `cEnergyBalance[1,2]` says that `vFLOW[1,1]` has to be `150.0`. They can’t both be true, so that’s the source of your infeasibility.

I suggest you take another look at your formulation.

---

<div class="post-metadata">

**Author:** ![mbah](https://avatars.discourse-cdn.com/v4/letter/m/bcef8e/32.png) [@mbah](https://discourse.julialang.org/u/mbah)\
**Post date:** [June 21, 2022, 7:25am UTC](https://discourse.julialang.org/t/infeasibilities-with-economic-dispatch-model/82571/5 "2022-06-21T07:25:19Z")

</div>

Thanks @odow for the insight. The main problem was correctly mapping plants to nodes and then slightly reformulating the energy balance. Model now generates a feasible solution without the need for slack variables.