Pyomo does not optimize

Question

I am really new in the pyomo. But it is not looking so comlicated and still I am not sure where is the error. The following optimizer is just an example and it does not reach a solution. The iteration stops from where it starts.

from pyomo.environ import *
from pyomo.opt import SolverFactory
import numpy as np
from scipy.interpolate import interp1d


x=[1,2,3,4,5,6,7]
y=[1,2,3,4,5,6,7]

# Adding objective function to the Pyomo model
def pyomo_objective(model):
    global x, y
    
    # Interpolation fit
    y_interp = interp1d(x, [value(model.x[j]) for j in model.I])
        
    y_new = list(y_interp(x))

    # Calculate the error using vectorized operations
    return np.sum(np.abs(np.array(y_new) - np.array(y)))

init_guess=[1,1,1,1,1,1,1]
# Define the Pyomo model
model = ConcreteModel()

# Define the Set
model.I = Set(initialize=range(len(init_guess)))  # Example indices for variables

# Define the Variables with initial values
model.x = Var(model.I, initialize=init_guess, within=NonNegativeReals)


# Add objective function to the Pyomo model
model.objective = Objective(rule=pyomo_objective, sense=minimize)

# Create a solver
glpk_solver = SolverFactory('glpk', tee=True)

# Set solver options
glpk_solver.options['tmlim'] = 3600  
glpk_solver.options['mipgap'] = 0.00001 

# Solve the optimization problem using IPOPT solver
results = glpk_solver.solve(model, tee=True)

# Check the solver status
if results.solver.termination_condition == TerminationCondition.optimal:
    print("Optimization successful")
    # Access optimal solution
    optimal_values = [model.x[i].value for i in model.I]
    print("Optimal values:", optimal_values)
    print("Objective value:", model.objective())
else:
    print("Solver failed to converge")

My initial suggestion was that interpolation is not fitting the optimizer but I am not sure. Please do not remove interpolation. It is a must in my context.

Answer 1

It is not working because 'glpk' is a Mixed Integer Linear Programming (MILP) solver and use of interpolation in the objective function makes the problem non-linear. There are two options to make this work, either find a solver that can handle interpolation inside the objective function or transform the objective function to piece-wise linear interpolation (difficulty of doing this depends if the dependence of y on x is convex or concave)