coin-or / python-mip

Python-MIP: collection of Python tools for the modeling and solution of Mixed-Integer Linear programs
Eclipse Public License 2.0
516 stars 90 forks source link

MIP clipping function #194

Closed FlorianStraubIng closed 2 years ago

FlorianStraubIng commented 3 years ago

Dear Python-MIP team,

Thank you for providing your tool, which has been very helpful to me. Currently, I’m trying to integrate a clipping function into a mip. I do this by formulating constraints that clip input values to the interval [clip_lower, clip_upper] using the BigM notation:

"""
Demonstration of mip implementation of a clipping function. Constrains clip input values to the interval
[clip_lower, clip_upper]
"""
from typing import List
import mip
import numpy as np
import matplotlib.pyplot as plt

def solve_mip(x, upper_clip=1, lower_clip=0):
    """
    Solves a mip that clips the input value between 0 and the defined upper border.

    Args:
        x (float): input value
        upper_clip (float): Max value where the input value is to be clipped
        lower_clip (float): Min value where the input value is to be clipped

    Returns:
        upper_bound and lower_bound of the input_value. (Note that we want this mip to have a defined output, so
        lower_bound and upper_bound are expected to coincide.)
    """
    m = mip.Model(solver_name='CBC', sense=mip.MAXIMIZE)
    y = m.add_var()
    d = m.add_var(var_type=mip.BINARY)
    d2 = m.add_var(var_type=mip.BINARY)
    M = 1000

    # d == 1 -> y = upper_clip; d == 0 -> y = x
    m.add_constr(y <= upper_clip)
    m.add_constr(y <= x + (1 - d2) * M)  # d2 == 1 -> y = x
    m.add_constr(y >= x - M * d)
    m.add_constr(y >= upper_clip - M * (1 - d))

    # lower boundary
    m.add_constr(y >= lower_clip)
    m.add_constr(y <= lower_clip + d2 * M)  # d2 == 0 -> y = lower_clip

    # To make sure that the solution that is found is the only solution, the mip is solved minimizing and maximizing the
    # output values -> if lower == upper we can guarantee that the mip's solution is well defined.
    m.objective = y
    m.optimize()
    upper = y.x

    m.objective = -y
    m.optimize()
    lower = y.x
    return lower, upper

def plot_values(plotdata: List[tuple]):
    """
    Plots multiple datasets in the same figure.

    Args:
        plotdata (List[tuple]): Each tuple contains (x_data, y_data, color_definition). color_definition can be 'r',
            'b', ...
    """
    fig = plt.figure()

    ax = fig.add_subplot(111)

    for x, y, color in plotdata:
        ax.plot(x, y, color)
        ax.spines['left'].set_position('zero')
        ax.spines['right'].set_color('none')
        ax.spines['bottom'].set_position('zero')
        ax.spines['top'].set_color('none')

        # remove the ticks from the top and right edges
        ax.xaxis.set_ticks_position('bottom')
        ax.yaxis.set_ticks_position('left')
    plt.show()

if __name__ == '__main__':
    x_in = np.linspace(-2, 2, 500)
    y_lower = []
    y_upper = []

    # run the mip for all values on x-axis
    for x in x_in:
        lower, upper = solve_mip(x)
        y_lower.append(lower)
        y_upper.append(upper)

    # If the constraints are strict enough, we expect the data from the lower and upper boundary to be identical
    plot_data = [
        (x_in, y_lower, 'blue'),
        (x_in, y_upper, 'red')
    ]

    print('lower and upper boundaries are identical: ', y_lower == y_upper)
    plot_values(plot_data)

clip_values

The problem is solved using CBC. As you can see in the plot, this mip can be solved for most of the input values. However there are some values that return None/Infeasible. Since these values are only a small subset of points in the linear domain of the problem, I assume that the mip is actually feasible in these regions. Might this be a numerical issue? Are there any settings in Python-MIP that I could adjust to solve this issue? Any hints on what might cause this issue would be much appreciated.

Thank you in advance!

Best,

Florian

tuliotoffolo commented 2 years ago

Since this has been fixed in CBC (coin-or/Cbc#412), I'll close this issue.