pythonoptimizationlinear-programminggekkomixed-integer-programming

How to handle optimizing across vectors in Gekko


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"

Solution

  • 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: