为什么这段代码 return 使用不同的程序会产生不同的结果？
Why does this code return different results when use different programs?
这段代码 objective 在 python gekko 中的值是 568，但是这段代码的 MATLAB 版本给出了 70。我不明白为什么。也许关于求解器选项？
from gekko import GEKKO
import numpy as np
m = GEKKO(remote=False)
capacity=np.array([70, 55, 51, 43, 41, 80])
demand=np.array([60, 57, 62, 38, 70])
cost=np.array([[5, 4, 5, 7, 2],
[2, 9, 2, 6, 3],
[6, 5, 1, 7, 9],
[7, 3, 9, 3, 5],
[4, 8, 7, 9, 7],
[2, 5, 4, 2, 1]])
x = m.Array(m.Var,(6,5),lb=0,integer=True)
for i in range(6):
for j in range(5):
m.Minimize(cost[i,j]*x[i,j])
for i in range(6):
m.Equation(m.sum(x[i,:])<=capacity[i])
for j in range(5):
m.Equation(m.sum(x[:,j])==demand[j])
m.options.solver = 1
m.solve()
print('Objective Function: ' + str(m.options.objfcnval))
print(x)
这里是给出 fval=70 的代码的 MATLAB 版本。也许我写错了代码，但看起来大部分是一样的，我不明白为什么。感谢帮助。
clear
clc
capacity_points=6;
demand_points=5;
capacity=[70 55 51 43 41 80];
demand=[60 57 62 38 70];
cost=[5 4 5 7 2;
2 9 2 6 3;
6 5 1 7 9;
7 3 9 3 5;
4 8 7 9 7;
2 5 4 2 1];
x=optimvar('x',capacity_points,demand_points,'Type','integer','LowerBound',0);
%expr=optimexpr;
expr=0;
for i=1:capacity_points
for j= 1:demand_points
expr=expr+cost(i,j)*x(i,j);
end
end
%const1=optimconstr(capacity_points);
for i=1:capacity_points
const1=sum(x(i,:)) <=capacity(i);
end
%const2=optimconstr(demand_points);
for j=1:demand_points
const2=sum(x(:,j)) ==demand(j);
end
prob=optimproblem;
prob.Objective=expr;
prob.Constraints.const1=const1;
prob.Constraints.const2=const2;
sol=solve(prob)
[sol,fval] = solve(prob)
看来 Gekko 中的解决方案是正确的。您可能还想 post MATLAB 代码，以便我们可以查看问题陈述之间是否存在差异。
当您想确定解决方案是否最优时，这里有一些一般的故障排除策略：
判断问题是否是非凸的。此问题是符合 minimize c x 形式的线性规划 (LP) 问题，服从 A x = b 和 A x < b。 LP 是凸的，因此局部解也是全局解。如果问题是非凸的，那么您可以尝试不同的初始猜测以查看解决方案是否发生变化，或者使用全局优化器为您自动化。
尝试不同的求解器。您可以设置一个循环并使用 m.options.SOLVER 更改求解器。如果您不想更改代码，请使用 m.options.SOLVER = 0 尝试所有求解器。这是 SOLVER=0.
Solver Objective Solution Time Status
-------------- ------------ ------------- ---------
APOPT (v1.0) 5.68000E+02 0.036 Success
BPOPT (v1.0) 5.68001E+02 0.016 Success
IPOPT (v3.12) 5.68000E+02 0.017 Success
IPOPT (v2.3) 0.00000E+00 0.000 Skip
SNOPT (v6.1) 0.00000E+00 0.000 Skip
MINOS (v5.5) 0.00000E+00 0.000 Skip
-------------- ------------------------------------
免费提供的求解器都得到了相同的解决方案，因此这表明多种求解器方法都达到了相同的共识。
- 对于除小规模问题（<100 个变量）之外的任何问题，都可能难以理解优化解决方案。另一种策略是降低问题的复杂性，直到您可以得出解析解或直到解易于验证。对于您的情况，列必须等于
demand 而行必须保持在 capacity 约束之下。
capacity=np.array([70, 55, 51, 43, 41, 80])
demand =np.array([60, 57, 62, 38, 70])
cost =np.array([[5, 4, 5, 7, 2],
[2, 9, 2, 6, 3],
[6, 5, 1, 7, 9],
[7, 3, 9, 3, 5],
[4, 8, 7, 9, 7],
[2, 5, 4, 2, 1]])
如果我用np.min(cost,axis=0)挑出每一列满足需求的最低成本项，就是[2 3 1 2 1]。如果我将其乘以 demand 和 np.dot(np.min(cost,axis=0),demand)，它会给出 objective 函数值 499。这是没有 capacity 约束的最小 objective 函数。它表明除非减少 demand 或 cost，否则 70 的 objective 不可能解决此问题。
虽然这些策略特定于您的问题，但它们可以应用于其他优化问题，以诊断和解决非直观的解决方案。如果求解器说它找到了一个解决方案，那么 Karush-Kuhn-Tucker conditions 满足最优性并且它至少是一个局部解决方案。
对编辑的回应
MATLAB 代码的问题是只有最后一个约束被强制执行，因为 const1 和 const2 在每个循环中都被重新定义。这给出了解决方案 sol.x.
0 0 0 0 0
0 0 0 0 0
0 0 0 0 0
0 0 0 0 0
0 0 0 0 0
0 0 0 0 70
您需要在 MATLAB 中包含所有等式和不等式约束。
clear
clc
capacity_points=6;
demand_points=5;
capacity=[70 55 51 43 41 80];
demand=[60 57 62 38 70];
cost=[5 4 5 7 2;
2 9 2 6 3;
6 5 1 7 9;
7 3 9 3 5;
4 8 7 9 7;
2 5 4 2 1];
x=optimvar('x',capacity_points,demand_points,'Type','integer','LowerBound',0);
%expr=optimexpr;
expr=0;
for i=1:capacity_points
for j= 1:demand_points
expr=expr+cost(i,j)*x(i,j);
end
end
%const1=optimconstr(capacity_points);
for i=1:capacity_points
const1(i)=sum(x(i,:)) <=capacity(i);
end
%const2=optimconstr(demand_points);
for j=1:demand_points
const2(j)=sum(x(:,j)) ==demand(j);
end
prob=optimproblem;
prob.Objective=expr;
prob.Constraints.const1 = const1;
prob.Constraints.const2 = const2;
sol=solve(prob)
[sol,fval] = solve(prob)
这给出了与 gekko 相同的 objective 函数值 568。
这段代码 objective 在 python gekko 中的值是 568，但是这段代码的 MATLAB 版本给出了 70。我不明白为什么。也许关于求解器选项？
from gekko import GEKKO
import numpy as np
m = GEKKO(remote=False)
capacity=np.array([70, 55, 51, 43, 41, 80])
demand=np.array([60, 57, 62, 38, 70])
cost=np.array([[5, 4, 5, 7, 2],
[2, 9, 2, 6, 3],
[6, 5, 1, 7, 9],
[7, 3, 9, 3, 5],
[4, 8, 7, 9, 7],
[2, 5, 4, 2, 1]])
x = m.Array(m.Var,(6,5),lb=0,integer=True)
for i in range(6):
for j in range(5):
m.Minimize(cost[i,j]*x[i,j])
for i in range(6):
m.Equation(m.sum(x[i,:])<=capacity[i])
for j in range(5):
m.Equation(m.sum(x[:,j])==demand[j])
m.options.solver = 1
m.solve()
print('Objective Function: ' + str(m.options.objfcnval))
print(x)
这里是给出 fval=70 的代码的 MATLAB 版本。也许我写错了代码，但看起来大部分是一样的，我不明白为什么。感谢帮助。
clear
clc
capacity_points=6;
demand_points=5;
capacity=[70 55 51 43 41 80];
demand=[60 57 62 38 70];
cost=[5 4 5 7 2;
2 9 2 6 3;
6 5 1 7 9;
7 3 9 3 5;
4 8 7 9 7;
2 5 4 2 1];
x=optimvar('x',capacity_points,demand_points,'Type','integer','LowerBound',0);
%expr=optimexpr;
expr=0;
for i=1:capacity_points
for j= 1:demand_points
expr=expr+cost(i,j)*x(i,j);
end
end
%const1=optimconstr(capacity_points);
for i=1:capacity_points
const1=sum(x(i,:)) <=capacity(i);
end
%const2=optimconstr(demand_points);
for j=1:demand_points
const2=sum(x(:,j)) ==demand(j);
end
prob=optimproblem;
prob.Objective=expr;
prob.Constraints.const1=const1;
prob.Constraints.const2=const2;
sol=solve(prob)
[sol,fval] = solve(prob)
看来 Gekko 中的解决方案是正确的。您可能还想 post MATLAB 代码，以便我们可以查看问题陈述之间是否存在差异。
当您想确定解决方案是否最优时，这里有一些一般的故障排除策略：
判断问题是否是非凸的。此问题是符合
minimize c x形式的线性规划 (LP) 问题，服从A x = b和A x < b。 LP 是凸的，因此局部解也是全局解。如果问题是非凸的，那么您可以尝试不同的初始猜测以查看解决方案是否发生变化，或者使用全局优化器为您自动化。尝试不同的求解器。您可以设置一个循环并使用
m.options.SOLVER更改求解器。如果您不想更改代码，请使用m.options.SOLVER = 0尝试所有求解器。这是SOLVER=0. 问题的输出
Solver Objective Solution Time Status
-------------- ------------ ------------- ---------
APOPT (v1.0) 5.68000E+02 0.036 Success
BPOPT (v1.0) 5.68001E+02 0.016 Success
IPOPT (v3.12) 5.68000E+02 0.017 Success
IPOPT (v2.3) 0.00000E+00 0.000 Skip
SNOPT (v6.1) 0.00000E+00 0.000 Skip
MINOS (v5.5) 0.00000E+00 0.000 Skip
-------------- ------------------------------------
免费提供的求解器都得到了相同的解决方案，因此这表明多种求解器方法都达到了相同的共识。
- 对于除小规模问题（<100 个变量）之外的任何问题，都可能难以理解优化解决方案。另一种策略是降低问题的复杂性，直到您可以得出解析解或直到解易于验证。对于您的情况，列必须等于
demand而行必须保持在capacity约束之下。
capacity=np.array([70, 55, 51, 43, 41, 80])
demand =np.array([60, 57, 62, 38, 70])
cost =np.array([[5, 4, 5, 7, 2],
[2, 9, 2, 6, 3],
[6, 5, 1, 7, 9],
[7, 3, 9, 3, 5],
[4, 8, 7, 9, 7],
[2, 5, 4, 2, 1]])
如果我用np.min(cost,axis=0)挑出每一列满足需求的最低成本项，就是[2 3 1 2 1]。如果我将其乘以 demand 和 np.dot(np.min(cost,axis=0),demand)，它会给出 objective 函数值 499。这是没有 capacity 约束的最小 objective 函数。它表明除非减少 demand 或 cost，否则 70 的 objective 不可能解决此问题。
虽然这些策略特定于您的问题，但它们可以应用于其他优化问题，以诊断和解决非直观的解决方案。如果求解器说它找到了一个解决方案，那么 Karush-Kuhn-Tucker conditions 满足最优性并且它至少是一个局部解决方案。
对编辑的回应
MATLAB 代码的问题是只有最后一个约束被强制执行，因为 const1 和 const2 在每个循环中都被重新定义。这给出了解决方案 sol.x.
0 0 0 0 0
0 0 0 0 0
0 0 0 0 0
0 0 0 0 0
0 0 0 0 0
0 0 0 0 70
您需要在 MATLAB 中包含所有等式和不等式约束。
clear
clc
capacity_points=6;
demand_points=5;
capacity=[70 55 51 43 41 80];
demand=[60 57 62 38 70];
cost=[5 4 5 7 2;
2 9 2 6 3;
6 5 1 7 9;
7 3 9 3 5;
4 8 7 9 7;
2 5 4 2 1];
x=optimvar('x',capacity_points,demand_points,'Type','integer','LowerBound',0);
%expr=optimexpr;
expr=0;
for i=1:capacity_points
for j= 1:demand_points
expr=expr+cost(i,j)*x(i,j);
end
end
%const1=optimconstr(capacity_points);
for i=1:capacity_points
const1(i)=sum(x(i,:)) <=capacity(i);
end
%const2=optimconstr(demand_points);
for j=1:demand_points
const2(j)=sum(x(:,j)) ==demand(j);
end
prob=optimproblem;
prob.Objective=expr;
prob.Constraints.const1 = const1;
prob.Constraints.const2 = const2;
sol=solve(prob)
[sol,fval] = solve(prob)
这给出了与 gekko 相同的 objective 函数值 568。