Integer constrained optimization with objective function expressed to execute another program call

Question

I'm working on an integer optimization problem and would like to try the GEKKO.

Problem Description : x1,x2,x3,...,x9,x10,x11 (which are integers in the range 1-16) are the eleven integer parameters of an objective function. I would like to find a set of x values to minimize the output of the objective function. However, the objective value is obtained by running another C/C++ program since the problem cannot be expressed in math. formula like your examples.

I suspect if the return value j of the objective function doesn't match the spec. Is j python float OK in the case ? But from the log and double check, the objective value of the init state is 0.06 which is correct. I don't know why it only runs for one iteration and then stop.

def f_per_particle(x1,x2,x3,x4,x5,x6,x7,x8,x9,x10,x11):
    inputs_list = [x1.value,x2.value,x3.value,x4.value,x5.value,x6.value,x7.value,x8.value,x9.value,x10.value,x11.value]
    #Use the inputs_list as input to run another SW program 
    #The SW program will return an objective value j which is float value
    #The goal is to find min. objective value with corresponding x1-x11 value
    return j



m = GEKKO() # Initialize gekko
m.options.SOLVER=1  # APOPT is an MINLP solver

# optional solver settings with APOPT
m.solver_options = ['minlp_maximum_iterations 500', \
                    # minlp iterations with integer solution
                    'minlp_max_iter_with_int_sol 10', \
                    # treat minlp as nlp
                    'minlp_as_nlp 0', \
                    # nlp sub-problem max iterations
                    'nlp_maximum_iterations 50', \
                    # 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']

# Initialize variables
# Integer constraints for x1 to x11
x1 = m.Var(value=8,lb=1,ub=16,integer=True)
x2 = m.Var(value=8,lb=1,ub=16,integer=True)
x3 = m.Var(value=8,lb=1,ub=16,integer=True)
x4 = m.Var(value=8,lb=1,ub=16,integer=True)
x5 = m.Var(value=8,lb=1,ub=16,integer=True)
x6 = m.Var(value=8,lb=1,ub=16,integer=True)
x7 = m.Var(value=8,lb=1,ub=16,integer=True)
x8 = m.Var(value=8,lb=1,ub=16,integer=True)
x9 = m.Var(value=8,lb=1,ub=16,integer=True)
x10 = m.Var(value=8,lb=1,ub=16,integer=True)
x11 = m.Var(value=8,lb=1,ub=16,integer=True)

m.Obj(f_per_particle(x1,x2,x3,x4,x5,x6,x7,x8,x9,x10,x11)) # Objective
m.solve() # Solve

Execution result as follows :
apm 130.226.87.11_gk_model0 <br><pre> ----------------------------------------------------------------
 APMonitor, Version 0.9.1
 APMonitor Optimization Suite
 ----------------------------------------------------------------


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

 ________________________________________________
 WARNING: objective equation           1 has no variables
 ss.Eqn(1)
 0 = 0.06
 ________________________________________________
 Number of state variables:             11
 Number of total equations: -            0
 Number of slack variables: -            0
 ---------------------------------------
 Degrees of freedom       :             11

 ----------------------------------------------
 Steady State Optimization with APOPT Solver
 ----------------------------------------------
Iter:     1 I:  0 Tm:      0.00 NLPi:    1 Dpth:    0 Lvs:    0 Obj:  6.00E-02 Gap:  0.00E+00
 Successful solution

 ---------------------------------------------------
 Solver         :  APOPT (v1.0)
 Solution time  :   1.410000000032596E-002 sec
 Objective      :   6.000000000000000E-002
 Successful solution
 ---------------------------------------------------

Results
x1: [8.0]
x2: [8.0]
x3: [8.0]
x4: [8.0]
x5: [8.0]
x6: [8.0]
x7: [8.0]
x8: [8.0]
x9: [8.0]
x10: [8.0]
x11: [8.0]
Objective: 0.06
(virtual-env) nan@TEK-CB-NAN-01:~/Swarm_Sandbox/GEKKO/mlp_fi_fwl/run$ 

Answer 1

Here is a complete but simplified version of your problem:

from gekko import GEKKO
def f(x1,x2):
    return (x1.value-3.3)**2 + (x2.value-5.7)**2
m = GEKKO() # Initialize gekko
m.options.SOLVER=1  # APOPT is an MINLP solver
x1 = m.Var(value=8,lb=1,ub=16,integer=True)
x2 = m.Var(value=8,lb=1,ub=16,integer=True)
m.Obj(f(x1,x2)) # Objective
m.solve() # Solve
print('Objective: ' + str((8-3.3)**2 + (8-5.7)**2))
print(x1.value[0],x2.value[0])

It produces a solution, similar to yours, where the solver does not move the values from the initial guess. This is because your objective function is returning a scalar, not a Gekko expression. Gekko requires an expression so that it can perform operator overloading to provide automatic differentiation for the solvers and provide exact 1st and 2nd derivatives in sparse form. You would need to modify your code to provide a function approximation of your objective function using something like a cspline or bspline based on sampling your objective function over a range of variables or else use a Gekko expression. Here is a modified version that uses a Gekko expression:

from gekko import GEKKO
def f(x1,x2):
    return (x1-3.3)**2 + (x2-5.7)**2
m = GEKKO() # Initialize gekko
m.options.SOLVER=1  # APOPT is an MINLP solver
x1 = m.Var(value=8,lb=1,ub=16,integer=True)
x2 = m.Var(value=8,lb=1,ub=16,integer=True)
m.Obj(f(x1,x2)) # Objective
m.solve() # Solve
print(x1.value[0],x2.value[0])

It produces the correct solution: x1=3.0 and x2=6.0.

 ----------------------------------------------
 Steady State Optimization with APOPT Solver
 ----------------------------------------------
Iter:     1 I:  0 Tm:      0.00 NLPi:    4 Dpth:    0 Lvs:    3 Obj:  2.05E-27 Gap:       NaN
--Integer Solution:   1.80E-01 Lowest Leaf:   2.05E-27 Gap:   1.80E-01
Iter:     2 I:  0 Tm:      0.00 NLPi:    1 Dpth:    1 Lvs:    2 Obj:  1.80E-01 Gap:  1.80E-01
Iter:     3 I:  0 Tm:      0.00 NLPi:    1 Dpth:    1 Lvs:    2 Obj:  9.00E-02 Gap:  1.80E-01
Iter:     4 I:  0 Tm:      0.00 NLPi:    1 Dpth:    1 Lvs:    1 Obj:  4.90E-01 Gap:  1.80E-01
--Integer Solution:   1.80E-01 Lowest Leaf:   5.80E-01 Gap:  -4.00E-01
Iter:     5 I:  0 Tm:      0.00 NLPi:    1 Dpth:    2 Lvs:    1 Obj:  5.80E-01 Gap: -4.00E-01
 Successful solution

 ---------------------------------------------------
 Solver         :  APOPT (v1.0)
 Solution time  :   1.390000000174041E-002 sec
 Objective      :   0.180000000000000     
 Successful solution
 ---------------------------------------------------