[Python] Handling on/off behavior in MILP optimization problems in GEKKO

I'm trying to understand how to implement ON/OFF behavior in optimizaiton problems to be solved in GEKKO effectively.

Consider the following scenario:

• 2 power generators, with upper and lower bounds
• one generator has a power loss of 30%
• the total generated power must meet an external power demand
• the goal is to minimize the energy loss

The above problem is implemented in Python:

import numpy as np
from gekko import GEKKO
import matplotlib.pyplot as plt

timesim = 24*1 # hours
timesteps = 1 # steps/hour
n = np.int64(timesim*timesteps + 1) # this is the vector length in GEKKO

Edh_demand =np.ones(timesim+1)*80;
Edh_demand[int(7*n/24):int(10*n/24)] = 200
Edh_demand[int(3*n/4)-1:int(3*n/4)+2] = 300


m = GEKKO(remote=False) # initialize gekko, körs lokalt
m.time = np.linspace(0, 24, 25) # Gekko-tid

Egen1 = m.Var(value = 30, lb = 30, ub = 200)
Egen2 = m.Var(value = 30, lb = 30, ub = 200)
Eloss = m.Var(value = 30*.3)
Eprod = m.Var(value = 30)
Edemand = m.Param(value = Edh_demand)

cost = m.Param(value = 10)

m.Equation(Eloss == Egen1*.3)
m.Equation(Eprod == Egen2 + Egen1 - Eloss)

m.Minimize(Eloss)

m.Equation(Eprod >= Edemand)

m.options.SOLVER = 3
m.options.IMODE = 6
# m.options.COLDSTART=2
m.options.COLDSTART=0

m.solve(disp=True)

plt.figure(5)
plt.subplot(2, 1, 1)
plt.plot(m.time, Egen1, label='gen 1')
plt.plot(m.time, Egen2, label='gen 2')
plt.legend()

plt.subplot(2, 1, 2)
plt.plot(m.time, Edemand, '--r')
plt.plot(m.time, Eprod)

plt.show()


I have the following questions/problems:

1. As it is implemented now, gen1 goes to its lower bound as expected. I want to constraint gen1 so that it can operate either at the lower bound or at 0, but not between 0 an gen1_lb.
2. I want to penalize turning off and on gen1.
3. I want to constraint gen1 so that every time it is turned off, it takes 2 sampling intervals before it can be turned on again.
4. Are there better ways (mainly for algoritm efficiency but also readibility) to implement such problem?

I believe one way to attack this problem is by creating an Integer variable (gen1_onoff) and multiply all occurances of gen1 by gen1_onoff. One can also add the following constraint

m.Equation(gen1 * gen1_onoff <= gen1_lb)


One could penalize variations of gen1_onoff in the objective function

dt_gen1_onoff = Var(value = 0)
m.Equation(dt_gen1_onoff == gen1_onoff.dt())
m.Minimize(dt_gen1_onoff)


I have no idea how to implement the constraint for problem 3, but suspect that m.if3() should be useful.

Thanks in advance for your help.

Continue reading...