google or-tools 无法获得最佳 LP 结果，如 gurobi 示例

Posted 2023-02-16

【中文标题】google or-tools 无法获得最佳 LP 结果，如 gurobi 示例

【英文标题】：google or-tools can not get to optimal LP result as in gurobi example

【发布时间】：2018-10-17 15:05:31

【问题描述】：

我正在尝试使用 Google 的 OR-Tools 解决 P.Williams 的经典数学编程书中的一个问题。 Gurobi 演示中使用了相同的示例：http://www.gurobi.com/resources/examples/food-manufacture-I

我的解决方案接近最优，但与 Gurobi 示例的答案不匹配（实际上这也是书中所说的正确答案）。

我对照约束检查了 OR-tools 解决方案，一切看起来都正确，但与最佳答案不匹配。

我做错了什么，还是因为 OR-tools GLOP 算法与 Gurobi 相比存在任何限制。

我的GitHub代码链接：https://github.com/APA092/optimum_global/blob/master/food_produce.py

from ortools.linear_solver import pywraplp

def main():
    data = [[110, 120, 130, 110, 115],
            [130, 130, 110, 90, 115],
            [110, 140, 130, 100, 95],
            [120, 110, 120, 120, 125],
            [100, 120, 150, 110, 105],
            [90, 100, 140, 80, 135]
        ];

    char = [8.8, 6.1, 2, 4.2, 5];

    solver = pywraplp.Solver('Linear_test', pywraplp.Solver.GLOP_LINEAR_PROGRAMMING)

#create variables
    buy = [[0 for x in range(len(data[0]))] for y in range(len(data))] 
    produce = [[0 for x in range(len(data[0]))] for y in range(len(data))] 
    store = [[0 for x in range(len(data[0]))] for y in range(len(data))] 

    for i in range(0, len(data)):
        for j in range(0, len(data[0])):
            buy[i][j] = solver.NumVar(0, solver.infinity(), 'buy')
            produce[i][j] = solver.NumVar(0, solver.infinity(), 'produce')
            store[i][j] = solver.NumVar(0, solver.infinity(), 'store')

#create objective
    objective = solver.Objective()
    for i in range(0, len(buy)):
        for j in range(0, len(buy[0])):
            objective.SetCoefficient(buy[i][j], data[i][j]*(-1))
            objective.SetCoefficient(produce[i][j], 150)
            objective.SetCoefficient(store[i][j], -5)
    objective.SetMaximization()

#create constraints
    #production not higher than capacity of machine 1
    constraint1 = [0]*len(produce)
    for i in range(0, 6):
        constraint1[i] = solver.Constraint(0, 200)
        for j in range(0, 2):
            constraint1[i].SetCoefficient(produce[i][j],1)

    #production not higher than capacity of machine 2
    constraint2 = [0]*len(produce)
    for i in range(0, 6):
        constraint2[i] = solver.Constraint(0, 250)
        for j in range(2, 5):
            constraint2[i].SetCoefficient(produce[i][j],1)

    #production not higher than resources available
    constraint3 = [[0 for x in range(len(data[0]))] for y in range(len(data))]
    for i in range(0, len(produce)):
        for j in range(0, len(produce[0])):
            constraint3[i][j] = solver.Constraint(0, solver.infinity())
            constraint3[i][j].SetCoefficient(produce[i][j], -1)
            constraint3[i][j].SetCoefficient(store[i][j], 1)
            constraint3[i][j].SetCoefficient(buy[i][j], 1)

    #storage limited to 1000 units
    constraint4 = [[0 for x in range(len(data[0]))] for y in range(len(data))]        
    for i in range(0, len(produce)):
        for j in range(0, len(produce[0])):    
            constraint4[i][j] = solver.Constraint(0, 1000)
            constraint4[i][j].SetCoefficient(store[i][j], 1)

    #initial storage
    constraint5 = [0]*len(store[0])
    for i in range(0, len(store[0])):
        constraint5[i] = solver.Constraint(500, 500)
        constraint5[i].SetCoefficient(store[0][i],1)
        constraint5[i].SetCoefficient(buy[0][i],-1)
        constraint5[i].SetCoefficient(produce[0][i],1)

    #final storage
    constraint6 = [0]*len(store[0])
    for i in range(0, len(store[0])):
        constraint6[i] = solver.Constraint(500, 500)
        constraint6[i].SetCoefficient(store[4][i],1)
        constraint6[i].SetCoefficient(buy[5][i],1)
        constraint6[i].SetCoefficient(produce[5][i],-1)

    #linking storage and production
    constraint7 = [[0 for x in range(len(data[0]))] for y in range(len(data))]
    for i in range(1,6):
        for j in range(0,len(data[0])):
            constraint7[i][j] = solver.Constraint(0, 0)
            constraint7[i][j].SetCoefficient(store[i-1][j],1)
            constraint7[i][j].SetCoefficient(store[i][j],-1)
            constraint7[i][j].SetCoefficient(buy[i][j],1)
            constraint7[i][j].SetCoefficient(produce[i][j],-1)

    #products characteristics HIGH
    constraint7 = [0]*len(produce)
    for i in range(0, len(produce)):
        constraint7[i] = solver.Constraint(-solver.infinity(), 0)
        for j in range(0, len(produce[0])):
            constraint7[i].SetCoefficient(produce[i][j], char[j]-6)

    #products characteristics LOW
    constraint8 = [0]*len(produce)
    for i in range(0, len(produce)):
        constraint8[i] = solver.Constraint(0, solver.infinity())
        for j in range(0, len(produce[0])):
            constraint8[i].SetCoefficient(produce[i][j], char[j]-3)




    solver.Solve()
    storage_cost = 0
    revenue = 0
    purchase_cost =0

    for i in range(0, len(produce)):
        for j in range(0, len(produce[0])):
            purchase_cost += data[i][j]*buy[i][j].solution_value()
            revenue += 150*produce[i][j].solution_value()
            storage_cost += 5*store[i][j].solution_value()


    profit = revenue - storage_cost - purchase_cost

    print "Profit - " + str(profit)
    print "Revenue - " + str(revenue)
    print "Storage Cost - " + str(storage_cost)
    print "Purchase Cost - " + str(purchase_cost)

【问题讨论】：

最好将实际代码粘贴到问题中，而不是链接。

感谢您的建议，已添加

我会说核心求解器不太可能被破坏，因为有很多数学保障（例如来自对偶理论）。在 API 方面更有可能，但这可能是一个纯粹的建模问题。这个问题可能已经复杂到让手动推理变得很麻烦。所以我建议做常见的事情：减少问题。产品;那些日子。然后推理就容易多了。当然，这可以使一些问题不可见，但也许其他一些减少可以提供帮助。全模型调试很麻烦。

另一种方法是在您将看到的单个/或多个约束块的注释之外，2、3 甚至更多的完整约束/块删除将导致解决方案，这已经低于您的目标。那么很容易说，剩下的块包含一些问题。

绝对是模型问题，我在上面的代码上运行了 GUROBI，它返回的值与 Glop 相同。

【参考方案1】：

似乎你发现你重新定义了constraint7。

作者自己打的补丁->https://github.com/APA092/optimum_global/commit/54b15836f5860ab56984ed6d139541e961088159