使用命名索引变量创建优化的初始点
此示例说明如何为具有命名索引变量的优化问题创建初始点。对于命名索引变量，指定初始点的最简单方法通常是使用 findindex 函数。
该问题是一个涉及原料油和精炼油混合的多周期库存问题。目标是在生产和库存能力以及混合油的“硬度”等各种约束下实现利润最大化。该问题源自 Williams [1]。
问题描述
该问题涉及制造商可用来精炼成食用油的两种植物油原料和三种非植物油原料。该制造商每月最多可精炼 200 吨植物油和 250 吨非植物油。制造商每种原料油可以储存 1000 吨，这是很有利的，因为购买原料油的成本取决于月份以及原料油的品种。每种油都有一个称为“硬度”的属性。混合油的硬度为成份油的线性加权硬度。
由于加工限制，制造商限制一个月内精炼油的品种不超过三种。另外，如果一个月内要精炼一种油，那么至少要精炼 20 吨这种油。最后，如果一个月内精炼一种植物油，那么也必须同时精炼非植物油 3。
每售出一吨油的收入是一个常数。而成本则是采购各种原料油的成本，因油种和月份而异，并且每个月每吨油的储存成本是固定的。精炼油的过程没有成本，但制造商不能储存精炼油（必须售出）。
输入问题数据
为规划周期和不同油种创建命名索引变量。
months = {'January','February','March','April','May','June'};
oils = {'veg1','veg2','non1','non2','non3'};
vegoils = {'veg1','veg2'};
nonveg = {'non1','non2','non3'};创建一些变量，包含与储存及使用相关的参数。
maxstore = 1000; % Maximum storage of each type of oil maxuseveg = 200; % Maximum vegetable oil use maxusenon = 250; % Maximum nonvegetable oil use minuseraw = 20; % Minimum raw oil use maxnraw = 3; % Maximum number of raw oils in a blend saleprice = 150; % Sale price of refined and blended oil storecost = 5; % Storage cost per period per oil quantity stockend = 500; % Stock at beginning and end of problem hmin = 3; % Minimum hardness of refined oil hmax = 6; % Maximum hardness of refined oil
将原料油的硬度指定为此向量。
h = [8.8,6.1,2,4.2,5.0];
将原料油的成本指定为此数组。数组的每一行代表一个月的原料油成本。第一行代表一月份的成本，最后一行代表六月份的成本。
costdata = [...
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];创建变量
创建下面这些问题变量：
sell，每月销售的每种油的数量
store，每月末每种油的库存量buy，每月每种油的采购量
此外，为了考虑每月精炼和售出的油量以及最低产量的约束，创建一个辅助二元变量 induse，当油品在一个月内售出时，该变量恰好为 1。
sell = optimvar('sell', months, oils, 'LowerBound', 0); buy = optimvar('buy', months, oils, 'LowerBound', 0); store = optimvar('store', months, oils, 'LowerBound', 0, 'UpperBound', maxstore); induse = optimvar('induse', months, oils, 'Type', 'integer', ... 'LowerBound', 0, 'UpperBound', 1);
将每月售出的油量指定为 produce。
produce = sum(sell,2);
创建目标
创建该问题的目标函数，计算收入，并减去采购成本和储存油的成本。
创建一个最大化优化问题，并将目标函数作为 Objective 属性包含在内。
prob = optimproblem('ObjectiveSense', 'maximize'); prob.Objective = sum(saleprice*produce) - sum(sum(costdata.*buy)) - sum(sum(storecost*store));
目标表达式很长。如果您愿意，可以使用 showexpr(prob.Objective) 命令查看它。
创建约束
该问题有几个您需要设置的约束。
6 月份每种油的库存量为 500。使用下界和上界设置此约束。
store('June', :).LowerBound = 500; store('June', :).UpperBound = 500;
制造商任何一个月都不能精炼超过 maxuseveg 种植物油。使用 约束和方程的表达式 设置此约束和所有后续约束。
vegoiluse = sell(:, vegoils); vegused = sum(vegoiluse, 2) <= maxuseveg;
制造商每月不能精炼超过 maxusenon 种非植物油。
nonvegoiluse = sell(:,nonveg); nonvegused = sum(nonvegoiluse,2) <= maxusenon;
混合油的硬度必须在 hmin through hmax 之间。
hardmin = sum(repmat(h, 6, 1).*sell, 2) >= hmin*produce; hardmax = sum(repmat(h, 6, 1).*sell, 2) <= hmax*produce;
每种油月底库存量等于月初量加上购入量减去售出量。
initstockbal = 500 + buy(1, :) == sell(1, :) + store(1, :); stockbal = store(1:5, :) + buy(2:6, :) == sell(2:6, :) + store(2:6, :);
如果在一个月内要精炼油，则至少要精炼 minuseraw 并售出。
minuse = sell >= minuseraw*induse;
确保在油品被精炼时，其对应的 induse 变量正好为 1。
maxusev = sell(:, vegoils) <= maxuseveg*induse(:, vegoils); maxusenv = sell(:, nonveg) <= maxusenon*induse(:, nonveg);
制造商每月最多可以售出 maxnraw 种油。
maxnuse = sum(induse, 2) <= maxnraw;
如果精炼了植物油，则油 non3 也必须精炼后售出。
deflogic1 = sum(induse(:,vegoils), 2) <= induse(:,'non3')*numel(vegoils);在问题中包含约束表达式。
prob.Constraints.vegused = vegused; prob.Constraints.nonvegused = nonvegused; prob.Constraints.hardmin = hardmin; prob.Constraints.hardmax = hardmax; prob.Constraints.initstockbal = initstockbal; prob.Constraints.stockbal = stockbal; prob.Constraints.minuse = minuse; prob.Constraints.maxusev = maxusev; prob.Constraints.maxusenv = maxusenv; prob.Constraints.maxnuse = maxnuse; prob.Constraints.deflogic1 = deflogic1;
创建并使用可行的初始点
创建正确维度的初始点。
x0.buy = zeros(size(buy)); x0.induse = zeros(size(induse)); x0.store = zeros(size(store)); x0.sell = zeros(size(sell));
设定起始点，只卖植物油 veg2 和非植物油 non3。要适当地设置这个初始点，请使用 findindex 函数。
numMonths = size(induse,1);
[idxMonths,idxOils] = findindex(induse,1:numMonths,{'veg2','non3'});
x0.induse(idxMonths,idxOils) = 1;满足植物油和非植物油的最大用量约束。
x0.sell(:,idxOils) = repmat([200,250],numMonths,1)
x0 = struct with fields:
buy: [6×5 double]
induse: [6×5 double]
store: [6×5 double]
sell: [6×5 double]
设定起始点，第一个月不采购原料油。
x0.buy(1,:) = 0;
基于每种油的初始储存量为 500，并且第一个月不采购，且持续使用 initstockbal 和 veg2，满足第一个月的 non3 约束。
x0.store(1,:) = [500 300 500 500 250];
使用 findindex 函数满足剩余库存余额约束 stockbal。
[idxMonths,idxOils] = findindex(store,2:6,{'veg2'});
x0.store(idxMonths,idxOils) = [100;0;0;0;500];
[idxMonths,idxOils] = findindex(store,2:6,{'veg1','non1','non2'});
x0.store(idxMonths,idxOils) = 500;
[idxMonths,idxOils] = findindex(store,2:6,{'non3'});
x0.store(idxMonths,idxOils) = [0;0;0;0;500];
[idxMonths,idxOils] = findindex(buy,2:6,{'veg2'});
x0.buy(idxMonths,idxOils) = [0;100;200;200;700];
[idxMonths,idxOils] = findindex(buy,2:6,{'non3'});
x0.buy(idxMonths,idxOils) = [0;250;250;250;750];检查初始点是否可行
issatisfied(prob,x0)
ans = logical
1
所有问题约束均在 x0 处得到满足，这意味着初始点是可行的。
使用初始点解决问题。
[sol,fval,exitstatus,output] = solve(prob,x0)
Solving problem using intlinprog.
Running HiGHS 1.7.1: Copyright (c) 2024 HiGHS under MIT licence terms
Coefficient ranges:
Matrix [1e-01, 2e+02]
Cost [5e+00, 2e+02]
Bound [1e+00, 1e+03]
RHS [3e+00, 5e+02]
Assessing feasibility of MIP using primal feasibility and integrality tolerance of 1e-06
Solution has num max sum
Col infeasibilities 0 0 0
Integer infeasibilities 0 0 0
Row infeasibilities 0 0 0
Row residuals 0 0 0
Presolving model
126 rows, 115 cols, 368 nonzeros 0s
96 rows, 85 cols, 508 nonzeros 0s
96 rows, 85 cols, 508 nonzeros 0s
MIP start solution is feasible, objective value is 39250
Solving MIP model with:
96 rows
85 cols (30 binary, 0 integer, 0 implied int., 55 continuous)
508 nonzeros
Nodes | B&B Tree | Objective Bounds | Dynamic Constraints | Work
Proc. InQueue | Leaves Expl. | BestBound BestSol Gap | Cuts InLp Confl. | LpIters Time
0 0 0 0.00% 1074750 39250 2638.22% 0 0 0 0 0.0s
R 0 0 0 0.00% 107512.962963 57078.571429 88.36% 0 0 0 80 0.0s
C 0 0 0 0.00% 105440.917048 94225 11.90% 142 27 0 150 0.1s
L 0 0 0 0.00% 104903.646216 99741.666667 5.18% 303 41 0 178 0.2s
10.0% inactive integer columns, restarting
Model after restart has 90 rows, 82 cols (27 bin., 0 int., 0 impl., 55 cont.), and 481 nonzeros
0 0 0 0.00% 104744.514099 99741.666667 5.02% 20 0 0 633 0.2s
0 0 0 0.00% 104744.514099 99741.666667 5.02% 20 20 0 685 0.2s
L 0 0 0 0.00% 104643.726006 100278.703704 4.35% 113 22 0 986 0.7s
B 18 0 6 3.52% 104643.726006 100278.703704 4.35% 132 22 27 3282 0.7s
Solving report
Status Optimal
Primal bound 100278.703704
Dual bound 100279.004373
Gap 0.0003% (tolerance: 0.01%)
Solution status feasible
100278.703704 (objective)
0 (bound viol.)
2.16004991671e-14 (int. viol.)
0 (row viol.)
Timing 1.16 (total)
0.01 (presolve)
0.00 (postsolve)
Nodes 405
LP iterations 7677 (total)
1240 (strong br.)
445 (separation)
2491 (heuristics)
Optimal solution found.
Intlinprog stopped because the objective value is within a gap tolerance of the optimal value, options.RelativeGapTolerance = 0.0001. The intcon variables are integer within tolerance, options.ConstraintTolerance = 1e-06.
sol = struct with fields:
buy: [6×5 double]
induse: [6×5 double]
sell: [6×5 double]
store: [6×5 double]
fval = 1.0028e+05
exitstatus =
OptimalSolution
output = struct with fields:
relativegap: 2.9983e-06
absolutegap: 0.3007
numfeaspoints: 6
numnodes: 405
constrviolation: 1.8929e-11
algorithm: 'highs'
message: 'Optimal solution found.↵↵Intlinprog stopped because the objective value is within a gap tolerance of the optimal value, options.RelativeGapTolerance = 0.0001. The intcon variables are integer within tolerance, options.ConstraintTolerance = 1e-06.'
solver: 'intlinprog'
参考资料
[1] Williams, H. Paul.Model Building in Mathematical Programming.Fourth edition.J. Wiley, Chichester, England.Problem 12.1, "Food Manufacture1."1999.
版权所有 2012–2024 The MathWorks, Inc.
