How to handle optimizing across vectors in Gekko

Question

I'm trying to use Gekko to optimize across time (a vector) to generate schedules based on several constraints, namely price and the number of times an event can occur (say we limit it to 5 per x7 below).

I'm getting the following error on my sum call: TypeError: x must be a python list of GEKKO parameters, variables, or expressions.

Even if I remove the sum, I also get the following error: TypeError: @error: Model Expression *** Error in syntax of function string: Invalid element: ame

The following variables are vectors with length 42: sp, base, neg_ln,

The following variables are scaler values (length 1): min, max, co_ln, co1,

I can get this to work by looking at only 1 row of the 42 rows at a time to optimize x1, but I'm unable to figure out how to translate this doing it across multiple rows in 1 go while constraining frequency (x7). Any guidance in appreciated.

m = GEKKO(server='XXX', remote=False)
m.options.IMODE = 2
m.options.MAX_ITER = 1000

#initialize variables
x1, x2, x3, x4, x5, x6, x7 = [m.Var() for i in range(7)]

x1.value = 1
x2.value = 1
x3.value = 1
x4.value = 0
x5.value = 0
x6.value = 0
x7.value = 1

neg_ln  =  m.Intermediate(-m.log(x1/sp))

vol1  =  m.Intermediate(co1+ base + (neg_ln*co_ln))
vol2  =  m.Intermediate(co2+ base + (neg_ln*co_ln))
vol3  =  m.Intermediate(co3+ base + (neg_ln*co_ln))
vol4  =  m.Intermediate(base + (neg_ln*co_ln))

total_vol = m.Intermediate((
(m.max2(0,base*(m.exp(vol1)-1)) * x3 +
m.max2(0,base*(m.exp(vol2)-1)) * x4 +
m.max2(0,base*(m.exp(vol3)-1)) * x5 +
m.max2(0,base*(m.exp(vol4)-1)) * x6) + base) * x7)

m.Equation(x3+x4+x5+x6 == 1) 
m.Equation(x3+x4+x5+x6 >= 0)
m.Equation(m.sum(x7)<= 5)

m.Equation(x1 >= min)
m.Equation(x1 <= max)

z1 = m.Var(lb = 3, ub = 15)
m.Equation(z1 == x1)
z2 = m.Var(lb = 3, ub = 15)
m.Equation(z2 == x2)
z3 = m.Var(lb = 0,ub = 1, integer=True)
m.Equation(z3 == x3)
z4 = m.Var(lb = 0,ub = 1, integer=True)
m.Equation(z4 == x4)
z5 = m.Var(lb = 0,ub = 1, integer=True)
m.Equation(z5 == x5)
z6 = m.Var(lb = 0,ub = 1, integer=True)
m.Equation(z6 == x6)
z7 = m.Var(lb = 0,ub = 1, integer=True)
m.Equation(z7 == x7)

m.Maximize(m.sum(simu_total_volume))

try:
m.solve(disp = True)
except:
"error"

Answer 1

Try using list comprehensions to define the Intermediate values. I plugged in sample random values for constants that are not defined in the question:

sp = np.random.rand(42)
base = np.random.rand(42)
co_ln = 2; co1 = 5; co2=6; co3=7

It is unclear if x values should be vectors. This is solved with IMODE=3 where the time discretization is explicitly defined versus IMODE=6 where the time discretization is defined by m.time. Here is a sample problem that solves successfully.

import numpy as np
from gekko import GEKKO
m = GEKKO(remote=False)
sp = np.random.rand(42)
base = np.random.rand(42)

x1,x2 = m.Array(m.Var,2,lb=3,ub=15,value=1)
x3,x4,x5,x6,x7 = m.Array(m.Var,5,lb=0,ub=1,value=1,integer=True)
# change default values and upper bound
x3.value = 1; x7.value=1; x7.upper = 5
x = [x1,x2,x3,x4,x5,x6,x7]

n = 42
co_ln = 2; co1 = 5; co2=6; co3=7
neg_ln=[m.Intermediate(-m.log(x[0]/sp[i])) for i in range(n)]
vol1  =[m.Intermediate(co1+ base[i] + (neg_ln[i]*co_ln)) for i in range(n)]
vol2  =[m.Intermediate(co2+ base[i] + (neg_ln[i]*co_ln)) for i in range(n)]
vol3  =[m.Intermediate(co3+ base[i] + (neg_ln[i]*co_ln)) for i in range(n)]
vol4  =[m.Intermediate(base[i] + (neg_ln[i]*co_ln)) for i in range(n)]

total_vol = [m.Intermediate((
  (m.max2(0,base[i]*(m.exp(vol1[i])-1)) * x[2] +
   m.max2(0,base[i]*(m.exp(vol2[i])-1)) * x[3] +
   m.max2(0,base[i]*(m.exp(vol3[i])-1)) * x[4] +
   m.max2(0,base[i]*(m.exp(vol4[i])-1)) * x[5]) + base[i]) * x[6]) for i in range(n)]

m.Equation(sum(x[2:7])==1) 
m.Equation(sum(x[2:7])>=0)
m.Equation(x7<=5)

m.Maximize(m.sum(total_vol))

m.options.SOLVER=1
m.solve(disp = True)

for i in range(7):
    print(f'x[{i+1}]: {x[i].value[0]}')

This produces a solution:

 Number of state variables:    892
 Number of total equations: -  718
 Number of slack variables: -  2
 ---------------------------------------
 Degrees of freedom       :    172

 ----------------------------------------------
 Steady State Optimization with APOPT Solver
 ----------------------------------------------
Iter:     1 I:  0 Tm:      0.32 NLPi:   61 Dpth:    0 Lvs:    3 Obj: -4.87E+02 Gap:       NaN
--Integer Solution:  -2.10E+01 Lowest Leaf:  -4.87E+02 Gap:   1.83E+00
Iter:     2 I:  0 Tm:      0.01 NLPi:    2 Dpth:    1 Lvs:    2 Obj: -2.10E+01 Gap:  1.83E+00
--Integer Solution:  -2.10E+01 Lowest Leaf:  -4.87E+02 Gap:   1.83E+00
Iter:     3 I:  0 Tm:      0.01 NLPi:    2 Dpth:    1 Lvs:    1 Obj: -2.10E+01 Gap:  1.83E+00
Iter:     4 I:  0 Tm:      0.01 NLPi:    2 Dpth:    1 Lvs:    0 Obj:  0.00E+00 Gap:  1.83E+00
 No additional trial points, returning the best integer solution
 Successful solution

 ---------------------------------------------------
 Solver         :  APOPT (v1.0)
 Solution time  :  0.36019999999999985 sec
 Objective      :  -20.951499798000548
 Successful solution
 ---------------------------------------------------


x[1]: 3.0
x[2]: 3.0
x[3]: 0.0
x[4]: 0.0
x[5]: 0.0
x[6]: 0.0
x[7]: 1.0

A few suggestions:

• Do not define z variables to add on integer constraints. Define integer constraints when declaring x variables.
• Use m.Array() to efficiently create arrays of variables. The x values can be defined as a vector m.Array(m.Var,7) or matrix m.Array(m.Var,(7,42)).
• Use lb and ub to define the bounds at initialization or use x[6].upper=5 to define a bound after it is initialized. Using m.Equation(x[6]<=5) is less efficient because the inequality constraint is solved implicitly.
• For future questions, please include a minimal and complete example of the problem so that the error can be reproduced. Several values were missing from the question.