Getting unoptimal solution for ILP problem

[timvandam]

Describe the bug A model is not returning an optimal (i.e. maximum objective value) solution, even though it is feasible

To Reproduce I've created the following reproduction code. It runs my model twice. The first time it lets the model do whatever it wants to find a solution. The second time I manually set some values in order to make it find a better solution. Running this script shows that the second time it indeed gets a better objective value, despite all the constraints being the exact same (except for it now forcing some variable values)

from typing import *
from mip import *

def solve_ilp_fair(trader_counts: List[int], total_card_counts: List[int], force_values: bool) -> Tuple[int, List[List[int]]]:
    """
    Distributes cards among traders fairly
    :param trader_counts: How many total cards each trader has
    :param total_card_counts: How many cards there are of each unique card
    :return: The distribution of cards. First axis is traders, second axis is unique cards
    """
    num_traders = len(trader_counts)
    num_unique_cards = len(total_card_counts)

    model = Model()
    model.verbose = 0

    # create variables
    card_counts = [[model.add_var(var_type=INTEGER) for _ in range(num_unique_cards)] for _ in range(num_traders)]

    # create constraints
    # ensure each trader has the correct number of cards
    for trader_idx in range(num_traders):
        model += xsum(card_counts[trader_idx]) == trader_counts[trader_idx]

    # ensure the total number of cards of each type is correct
    for card_idx in range(num_unique_cards):
        model += xsum(card_counts[trader_idx][card_idx] for trader_idx in range(num_traders)) == total_card_counts[card_idx]

    # set fairness objective
    minimum_complete_collections = model.add_var(var_type=INTEGER)
    for trader_idx in range(num_traders):
        for card_idx in range(num_unique_cards):
            model += minimum_complete_collections <= card_counts[trader_idx][card_idx]
    model.objective = maximize(minimum_complete_collections)

    # optimize fairness
    status = model.optimize()

    if status != OptimizationStatus.OPTIMAL:
        raise Exception(f"No fair solution found. Status: {status}")

    # set maximum complete collections objective, ensuring fairness
    model += minimum_complete_collections == int(minimum_complete_collections.x)
    total_complete_collections = model.add_var(var_type=INTEGER)
    trader_complete_collections = [model.add_var(var_type=INTEGER) for _ in range(num_traders)]
    for trader_idx in range(num_traders):
        for card_idx in range(num_unique_cards):
            model += trader_complete_collections[trader_idx] <= card_counts[trader_idx][card_idx]
    model += total_complete_collections == xsum(trader_complete_collections)
    model.objective = maximize(total_complete_collections)

    if force_values:
        model += card_counts[0][0] == 4579
        model += card_counts[0][1] == 4579
        model += card_counts[1][0] == 2743
        model += card_counts[1][1] == 2742
        model += card_counts[2][0] == 5411
        model += card_counts[2][1] == 3822

    # optimize maximum complete collections
    status = model.optimize()

    if status != OptimizationStatus.OPTIMAL:
        raise Exception(f"No fair + maximum solution found. Status: {status}")

    return int(model.objective_value), [[int(card_counts[trader_idx][card_idx].x) for card_idx in range(num_unique_cards)] for trader_idx in range(num_traders)]

def main():
    # a test case
    trader_counts = [9158, 5485, 9233]
    total_card_counts = [12733, 11143]

    # solve
    print("Solving with ILP...")
    objective_value, solution = solve_ilp_fair(trader_counts, total_card_counts, False)

    print(f"Solution with objective value {objective_value}:")
    for row in solution:
        print(row)

    print("Solving with ILP, forcing custom values...")
    objective_value, solution = solve_ilp_fair(trader_counts, total_card_counts, True)

    print(f"Solution with objective value {objective_value}:")
    for row in solution:
        print(row)

if __name__ == "__main__":
    main()

Expected behavior I would expect it to find an optimal solution. It does not appear to be doing this currently

Desktop (please complete the following information):

• Operating System, version: Ubuntu 20.04.5 LTS x86_64
• Python version: Python 3.8
• Python-MIP version (we recommend you to test with the latest version): 1.14.1, 1.14.2

Additional context Add any other context about the problem here.

[timvandam]

It appears that replacing

    model += total_complete_collections == xsum(trader_complete_collections)
    model.objective = maximize(total_complete_collections)

with

    model.objective = maximize(xsum(trader_complete_collections))

solves the issue. Strange

[timvandam]

I found another case where even this change does not fix it. I am using MIP to verify the correctness of an algorithm, but it seems that I'm using the algorithm to verify MIP now 😁

[tuliotoffolo]

Hi Tim,

Here’s the log of your run:

Cbc0011I Exiting as integer gap of 1 less than 1e-10 or 0.01%

Indeed, the solver (CBC in this case) considered the solution to be optimal since it has optimality gap smaller than 0.01%. Concerning the different result when changing the objective, this probably happened due to some level of "randomness" or high sensitivity to initial condition by the solver. For more info, I recommend reading this very nice paper by M. Fischetti and M. Monaci (http://www.dei.unipd.it/~fisch/papers/exploiting_erraticism_in_search.pdf).

To get the result you want and verify the correctness of your algorithm using MIP, consider setting parameter max_mip_gap to a smaller value, such as 0.00001. ;)

Happy holidays!

[timvandam]

Excellent, this worked. Thank you for the quick response :smile: Happy holidays to you too