Using Gekko in Python to run time series linear optimization

Viewed 588

I am trying to run an electricity arbitrage model in Gekko. I have an electricity price array for every hour of a year (8760 total hours), a battery of energy size E, and for every hour I want to decide whether to charge the battery or discharge it, based on minimizing the electricity cost, and keep track of the energy in the battery constrained to hold no less than 0 energy and no more than E.

I've tried many times, most recently getting the error that the optimize equation exceeds the limits

from gekko import Gekko

m = Gekko()

#variables
E_battery = m.Var(lb=0, ub=366.2, value=0) #energy in battery at time t, battery size 366 MWh
Pc = m.Var(lb=0, ub=50) #charge power, 50 MW max
Pd = m.Var(lb=0, ub=36.6)  #discharge power, max 36 MW
E_price = m.Param(electricity_price[:,1])
m.time = np.linspace(0,8759, 8760)

m.Equation(E_battery.dt() == (1-delta)*E_battery + roundtrip_eff*(Pc - Pd))

m.Obj(sum(E_price[i]*Pc for i in range(8760)))
m.options.IMODE = 7
m.solve()
1 Answers

If you are optimizing over the entire time horizon then you'll want to switch to IMODE=6. You may be having a problem with:

m.Obj(sum(E_price[i]*Pc for i in range(8760)))

because it creates an objective expression that is very long. For dynamic optimization problems with IMODE=6 you can use the following instead:

m.Obj(E_price*Pc)

Gekko automatically does the summation of all time points in that expression. The script isn't a complete example so it is hard to verify what isn't working. If you can post a minimal and complete example that shows the issue then it will be easier to provide helpful feedback. There are additional examples in the machine learning and dynamic optimization course, especially with economic optimization or some of the other benchmark problems.

Summation

You can sum a variable such as E_battery with something like the following:

from gekko import GEKKO
import numpy as np
import matplotlib.pyplot as plt  
m = GEKKO(remote=False)
m.time = [0,1,2,3,4,5]
x = m.Param(value=m.time)
E_battery = m.Var()
m.Equation(E_battery.dt()==x)
m.options.IMODE = 4
m.solve(disp=False)
print('x = ' + str(x.value))
print('E_battery = ' + str(E_battery.value))

This produces the solution:

x = [0.0, 1.0, 2.0, 3.0, 4.0, 5.0]
E_battery = [0.0, 1.0, 3.0, 6.0, 10.0, 15.0]
Related