Question: I am using the GEKKO optimization package (https://gekko.readthedocs.io/en/latest/index.html) in python. Since my problem sometimes includes integer variables I run with the following options to get optimization for a steady state mixed integer problem (APOPT): m.options.IMODE = 3 m.options.SOLVER = 1  The optimizer finds a solution that is not optimal. Consider the following minimal working example: from gekko import GEKKO m = GEKKO(remote=False) # Define objective def obj(x1, x2): return x1*0.04824647331436083 + x2*0.1359023029124411 + x1*x2*(0.5659570743890336) # Define variables x1 = m.Var(lb=-0.7280376309420935, ub=0.2719623690579065) x2 = m.Var(lb=-0.8912733608888134, ub=0.10872663911118663) u = m.Intermediate(obj(x1,x2)) m.Maximize(u) # Solve the problem m.options.SOLVER = 'APOPT' m.solve(disp=False) # Print the results print(f'x1: {x1.value[0]}') print(f'x2: {x2.value[0]}') print(f'u: {u.value[0]}')  This gives the output: x1: 0.27196236906 x2: 0.10872663911 u: 0.044632524297  However, consider the following solution and its utility: x1 = -0.7280376309420935 x2 = -0.8912733608888134 obj(x1, x2)  This results in a utility of 0.2109871851434535

Answer: A contour plot verifies that there are two local maxima. Gekko solvers report local minima or maxima when the Karush Kuhn Tucker (KKT) conditions are satisfied.  from gekko import GEKKO m = GEKKO(remote=False) # Define objective def obj(x1, x2): return x1*0.04824647331436083 + x2*0.1359023029124411 + x1*x2*(0.5659570743890336) # Define variables x1 = m.Var(lb=-0.7280376309420935, ub=0.2719623690579065) x2 = m.Var(lb=-0.8912733608888134, ub=0.10872663911118663) u = m.Intermediate(obj(x1,x2)) m.Maximize(u) # Solve the problem m.options.SOLVER = 'APOPT' m.solve(disp=False) # Print the results print(f'x1: {x1.value[0]}') print(f'x2: {x2.value[0]}') print(f'u: {u.value[0]}') import numpy as np import matplotlib.pyplot as plt # Generate grid of values x1v = np.linspace(-0.7280376309420935, 0.2719623690579065, 400) x2v = np.linspace(-0.8912733608888134, 0.10872663911118663, 400) x1g, x2g = np.meshgrid(x1v, x2v) ug = obj(x1g, x2g) # Create contour plot plt.figure(figsize=(6, 3.5)) contour = plt.contourf(x1g, x2g, ug, levels=50, cmap='viridis') plt.colorbar(contour); plt.xlabel('x1'); plt.ylabel('x2') plt.tight_layout(); plt.show()  For problems with multiple local minima or maxima, there are multi-start methods to assist in finding the global optimum: Multi-start implementation with GEKKO Below is an implementation of hyperopt that is used to find the global optimum. hyperopt is a Python package for performing hyperparameter optimization with a variety of optimization algorithms including random search, Tree-structured Parzen Estimator (TPE), and adaptive TPE, as well as a simple and flexible way to define the search space for the hyperparameters. from gekko import GEKKO from hyperopt import fmin, tpe, hp from hyperopt import STATUS_OK, STATUS_FAIL # Define the search space for the multi-start parameters space = {'x1': hp.quniform('x1', -0.72, 0.27, 0.1), 'x2': hp.quniform('x2', -0.89, 0.10, 0.1)} def objective(params): m = GEKKO(remote=False) x1 = m.Var(lb=-0.7280376309420935, ub=0.2719623690579065) x2 = m.Var(lb=-0.8912733608888134, ub=0.10872663911118663) x1.value = params['x1'] x2.value = params['x2'] u = m.Intermediate(x1*0.04824647331436083 + x2*0.1359023029124411 + x1*x2*(0.5659570743890336)) m.Maximize(u) m.options.SOLVER = 1 m.solve(disp=False,debug=False) obj = m.options.objfcnval if m.options.APPSTATUS==1: s=STATUS_OK else: s=STATUS_FAIL m.cleanup() return {'loss':obj, 'status': s, 'x1':x1.value[0],'x2':x2.value[0],'u':u.value[0]} best = fmin(objective, space, algo=tpe.suggest, max_evals=50) sol = objective(best) print(f"Solution Status: {sol['status']}") print(f"Objective: {-sol['loss']:.2f}") print(f"x1: {sol['x1']}") print(f"x2: {sol['x2']}") print(f"u: {sol['u']}")  This produces the global maximum with the TPE algorithm to find the best optimal solution by changing the initial guess values using a Bayesian approach. Solution Status: ok Objective: 0.21 x1: -0.72803763094 x2: -0.89127336089 u: 0.21098718514