Gekko sum of intermediate variables condition

Viewed 271

I have a simple optimization problem, where I am trying to optimize the number of leads, each of which have a certain cost (which varies based on the number of leads, e.g. 1 lead may cost 1$ while 2 leads may cost 1.5$ each). The cost of each lead is calculated using the following function. The comments are tries I made to fix different issues.

# Here x is leads and returns spend
def objective_inverse(x, a, b, c):
     if c == 0:
         return 0
     else:
         return np.real(a + (x/c)**(1./b))

Here is the variable reg_feats , where each tuple in the list is (a,b,c) in the previous function.

reg_feats = [(0.0, 1.0, 0.0008237661018136143), (0.9999999999999811, 0.8536915943397881, 0.004454688129841911), (0.0, 1.0, 0.0004855869694355375), (0.0, 1.0, 0.00038427101866404334), (8.308985001840891e-10, 1.097417548409198, 0.0009170956872015353), (0.999999999998488, 1.9999999897072056, 1.144122017534284e-07), (0.9999999999825329, 0.8097037523302283, 0.20606857088111075), (0.9999999999999963, 1.1649279015402045, 0.009316713936903972), (0.0, 1.0, 0.0034430519212229715), (0.0, 1.0, 0.0007980573950249244), (0.0, 1.0, 0.0009069008844368589)]

and also have a list for each variable's bound.

bounds = [(0, 4), (0, 20), (0, 3), (0, 1), (0, 114), (0, 14), (0, 208), (0, 529), (0, 1), (0, 4), (0, 4)]

Here is the code for the optimization:

m = GEKKO(remote=not bool(i))

m.options.SOLVER=1
# optional solver settings with APOPT
m.solver_options = ['minlp_maximum_iterations 50000', \
                # minlp iterations with integer solution
                'minlp_max_iter_with_int_sol 100', \
                # treat minlp as nlp
                'minlp_as_nlp 0', \
                # nlp sub-problem max iterations
                'nlp_maximum_iterations 500', \
                # 1 = depth first, 2 = breadth first
                'minlp_branch_method 1', \
                # maximum deviation from whole number
                'minlp_integer_tol 0.05', \
                # covergence tolerance
                'minlp_gap_tol 0.01']

Here x is the number of leads for each campaign.

x = m.Array(m.Var, (len(bounds)),integer=True)
# Here x is the number of leads and is the decision variable
# x variables are nb of leads, not $ to be allocated
for i, xi in enumerate(x):
    lb, ub = bounds[i]
    xi.lower = lb
    xi.upper = ub
    xi.value = ub

This intermediate variable is to get the cost of each lead. This is used in the next constraint.

sums = [m.Intermediate(objective_inverse(xi,*y)) for xi,y in zip(x,reg_feats)]

max_budget = 150000
m.Maximize(m.sum(x))
c = m.Const(max_budget, 'Budget')
m.Equation(m.sum(sums) < c)
m.solve(disp=True)

Here, if I run this code, I get an error saying no solution has been found and Warning: no more possible trial points and no integer solution

Also, if I remove the max_budget constraint, I get a solution where the total amount spent is lower than the max_budget, which means that adding the constraint shouldn't change the solution.

I know there is a problem as to how I defined the constraint, but don;t know how to fix it.

1 Answers

The problem doesn't seem to be with the constraint but rather the starting point of x. When I ran your code, I got the solution (I didn't get any errors). But started getting the same error when I changed max_budget to 100000. I suspect this is because APOPT requires a good initial guess. Since I didn't have a good guess, I used np.random to generate the starting point. I tried with a few different max_budget and it works (although it took a few tries to get a good initial point).

import numpy as np
from gekko import GEKKO

# Here x is leads and returns spend
def objective_inverse(x, a, b, c):
     if c == 0:
         return 0
     else:
         return np.real(a + (x/c)**(1./b))

reg_feats = [(0.0, 1.0, 0.0008237661018136143), (0.9999999999999811, 0.8536915943397881, 0.004454688129841911), (0.0, 1.0, 0.0004855869694355375), (0.0, 1.0, 0.00038427101866404334), (8.308985001840891e-10, 1.097417548409198, 0.0009170956872015353), (0.999999999998488, 1.9999999897072056, 1.144122017534284e-07), (0.9999999999825329, 0.8097037523302283, 0.20606857088111075), (0.9999999999999963, 1.1649279015402045, 0.009316713936903972), (0.0, 1.0, 0.0034430519212229715), (0.0, 1.0, 0.0007980573950249244), (0.0, 1.0, 0.0009069008844368589)]
bounds = [(0, 4), (0, 20), (0, 3), (0, 1), (0, 114), (0, 14), (0, 208), (0, 529), (0, 1), (0, 4), (0, 4)]
m = GEKKO(remote=False)

m.options.SOLVER=1
# optional solver settings with APOPT
m.solver_options = ['minlp_maximum_iterations 50000', \
                # minlp iterations with integer solution
                'minlp_max_iter_with_int_sol 100', \
                # treat minlp as nlp
                'minlp_as_nlp 0', \
                # nlp sub-problem max iterations
                'nlp_maximum_iterations 500', \
                # 1 = depth first, 2 = breadth first
                'minlp_branch_method 1', \
                # maximum deviation from whole number
                'minlp_integer_tol 0.05', \
                # covergence tolerance
                'minlp_gap_tol 0.01']

x = m.Array(m.Var, (len(bounds)),integer=True)
# Here x is the number of leads and is the decision variable
# x variables are nb of leads, not $ to be allocated
for i, xi in enumerate(x):
    lb, ub = bounds[i]
    xi.lower = lb
    xi.upper = ub
    xi.value = np.random.randint(lb,ub)

print("Initial Value: ",x)

sums = [m.Intermediate(objective_inverse(xi,*y)) for xi,y in zip(x,reg_feats)]

max_budget = 90000
# max_budget = 100000
# max_budget = 120000
# max_budget = 150000
m.Maximize(m.sum(x))
c = m.Const(max_budget, 'Budget')
m.Equation(m.sum(sums) < c)
m.solve(disp=True)


x_ans =[i[0] for i in x]
print('x : ', x_ans)
print('Constraint: ', np.sum([(objective_inverse(xi,*y)) for xi,y in zip(x_ans,reg_feats)]))

Output:

Initial Value:  [2 7 0 0 113 5 205 151 0 3 2]
 ----------------------------------------------------------------
 APMonitor, Version 0.9.2
 APMonitor Optimization Suite
 ----------------------------------------------------------------


 --------- APM Model Size ------------
 Each time step contains
   Objects      :  1
   Constants    :  1
   Variables    :  24
   Intermediates:  11
   Connections  :  12
   Equations    :  24
   Residuals    :  13

 Number of state variables:    24
 Number of total equations: -  13
 Number of slack variables: -  1
 ---------------------------------------
 Degrees of freedom       :    10

 ----------------------------------------------
 Steady State Optimization with APOPT Solver
 ----------------------------------------------
Iter:     1 I:  0 Tm:      0.00 NLPi:   12 Dpth:    0 Lvs:    2 Obj: -8.85E+02 Gap:       NaN
Iter:     2 I:  0 Tm:      0.00 NLPi:    3 Dpth:    1 Lvs:    3 Obj: -8.85E+02 Gap:       NaN
Iter:     3 I:  0 Tm:      0.00 NLPi:    3 Dpth:    2 Lvs:    4 Obj: -8.85E+02 Gap:       NaN
Iter:     4 I:  0 Tm:      0.00 NLPi:    3 Dpth:    3 Lvs:    5 Obj: -8.85E+02 Gap:       NaN
Iter:     5 I:  0 Tm:      0.00 NLPi:    3 Dpth:    4 Lvs:    6 Obj: -8.84E+02 Gap:       NaN
Iter:     6 I:  0 Tm:      0.00 NLPi:    3 Dpth:    5 Lvs:    7 Obj: -8.84E+02 Gap:       NaN
--Integer Solution:  -8.84E+02 Lowest Leaf:  -8.85E+02 Gap:   8.39E-04
Iter:     7 I:  0 Tm:      0.00 NLPi:    2 Dpth:    6 Lvs:    7 Obj: -8.84E+02 Gap:  8.39E-04
 Successful solution

 ---------------------------------------------------
 Solver         :  APOPT (v1.0)
 Solution time  :  0.0677 sec
 Objective      :  -884.
 Successful solution
 ---------------------------------------------------


x :  [0.0, 18.0, 0.0, 0.0, 114.0, 14.0, 208.0, 529.0, 1.0, 0.0, 0.0]
Constraint:  89187.87455034524

I understand randomly generating starting points might not be feasible for a bigger problem without a good initial guess. But unfortunately, I have no idea as to how APOPT could be used in such a case.

Related