研究线性不可行性

这个例子展示了如何研究导致问题不可行的线性约束。有关这些技术的更多细节，请参见Chinneck［1］和［2］．

如果线性约束导致问题不可行，则可能希望找到不可行的约束子集，但删除子集的任何成员会使子集的其余部分可行。这样一个子集的名称是不可约约束的不可行子集，缩写为IIS。相反，您可能希望找到可行的约束的最大基数子集。这个子集叫做a最大可行的子集,缩写maxf。这两个概念是相关的，但并不完全相同。一个问题可以有许多不同的iis，其中一些具有不同的基数。

这个示例展示了查找IIS的两种方法，以及获取可行约束集的两种方法。这个例子假设所有给定的边界都是正确的，这意味着磅和乌兰巴托争论没有错误。

不可行的例子

创建一个随机矩阵一个表示大小为150 × 15的线性不等式。设置相应的向量b将这些值的5%更改为–10。

N=15；rng违约一个= randn ((10 * N, N));b = 10 *的(大小(A, 1), 1);Aeq = [];说真的= [];B (rand(size(B)) <= 0.05) = -10;f =的(N - 1);磅= - f;乌兰巴托= f;

检查问题这是不可行的。

[x，fval，exitflag，output，lambda]=linprog（f，A，b，Aeq，beq，lb，ub）；

没有找到可行的解决办法。Linprog停止了，因为没有点满足约束。

删除过滤器

要标识IIS，请执行以下步骤。给定一组编号为1到N的线性约束，其中所有问题约束都是不可行的:

为每一个我从1到N：

暂时删除约束我从这个问题。
测试由此产生的问题的可行性。
如果问题不受约束是可行的我，返回约束我问题的办法。
如果问题没有约束是不可行的我，不返回约束我这个问题

继续下一个我（最高价值N)．

在此过程结束时，问题中保留的约束将形成IIS。

实现此程序的MATLAB®代码，请参阅deletionfilter辅助函数在这个例子到此结束．

请注意：如果在本例中使用live脚本文件，则deletionfilter函数已经包含在文件的末尾。否则，您需要在.m文件的末尾创建这个函数，或者将其作为文件添加到MATLAB路径中。本例后面使用的其他helper函数也是如此。

看到…的效果deletionfilter在示例数据上。

[ineqs,方程式,ncall] = deletionfilter (A, b, Aeq,说真的,磅,乌兰巴托);

这个问题没有等式约束。求不等式约束的指标和值b (iis)．

iis =找到(ineqs)

iis = 114

b (iis)

ans = -10

只有一个不等式约束和有界约束导致问题不可行。约束条件是

A（iis，：）*x<=b（iis）．

为什么这个约束和边界都是不可行的?求这一行的绝对值的和一个．

disp(总和(abs ((iis,:))))

8.4864

由于边界的限制x向量的值在-1和1之间，等等:一个(iis) * x不能少于b (iis)= –10.

有多少linprog调用了deletionfilter执行?

disp (ncall)

这个问题有150个线性约束，所以这个函数叫做linprog150次。

弹性过滤器

作为删除过滤器(它检查每个约束)的替代方法，可以尝试弹性过滤器。这个过滤器的工作原理如下。

首先，允许每个约束我被正数违反e(我)，其中等式约束具有加和减的正弹性值。

$\begin{array}{l} {一个}_{我 n e 问} x \leq b_{我 n e 问} + e \\ {一个}_{e 问} x ＝ b_{e 问} + e_{1} - e_{2} \end{array}$

接下来，解决相关的线性规划问题（LP）

$\underset{x ， e}{闵} \sum e_{我}$

受列出的约束和with的约束 $e_{我} \geq 0$ ．

如果关联的LP有一个解，删除所有严格为正的关联约束 $e_{我}$ ，并在索引(潜在的IIS成员)列表中记录这些约束。返回到前面的步骤来解决新的、减少的关联LP。
如果相关的LP无解(不可行)或没有严格正相关 $e_{我}$ ，退出过滤器。

与删除过滤器相比，弹性过滤器可以以更少的迭代次数退出，因为它可以一次将多个索引引入IIS，并且可以在不查看整个索引列表的情况下停止。但是，该问题比原始问题有更多的变量，其结果索引列表可能比IIS更大。要在运行弹性筛选器后查找IIS，请对结果运行删除筛选器。

实现此滤波器的MATLAB®代码，请参阅elasticfilter辅助函数在这个例子到此结束．

看到…的效果elasticfilter在示例数据上。

[ineqselastic、eqselastic、ncall]=...Aeq elasticfilter (A, b,说真的,磅,乌兰巴托);

这个问题没有等式约束。求不等式约束的指标。

iiselastic =找到(ineqselastic)

iiselastic =5×12 60 82 97 114

弹性IIS列出了5个约束，而删除过滤器只找到一个。对返回的集合运行删除筛选器以找到真正的IIS。

A_1=A（ineqseastic>0，：）；b_1=b（ineqseastic>0）；[dineq_-iis，deq_-iis，ncall2]=deletionfilter（A_1，b_-1，Aeq，beq，lb，ub）；iiselasticdeletion=find（dineq_-iis）

iiselasticdeletion = 5

弹性滤波结果中的第五个约束，不等式114，为IIS。这个结果与删除过滤器的答案一致。这两种方法之间的区别在于，弹性过滤器和删除过滤器的组合使用更少linprog调用。显示总数linprog在删除过滤器之后，由弹性过滤器使用的调用。

disp (ncall + ncall2)

在循环中删除IIS

通常，获取单个IIS并不能使您找到优化问题失败的所有原因。要纠正不可行的问题，可以重复地找到IIS并将其从问题中删除，直到问题变得可行为止。

下面的代码显示了如何一次从问题中删除一个IIS，直到问题变得可行。在算法删除任何约束之前，该代码使用索引技术跟踪约束在原始问题中的位置。

代码使用布尔向量跟踪问题中的原始变量活跃的表示当前约束（行）的一个矩阵和布尔向量activeAeq来表示当前的约束Aeq矩阵。当添加或删除约束时，代码索引到一个或Aeq这样行数就不会改变，即使约束的数量改变了。

运行此代码将返回idx2中非零元素下标的向量活跃的：

idx2 =找到(activeA)

假设变量布尔向量的长度和idx2. 然后

idx2(找到(var))

表达变量作为原始问题变量的索引。通过这种方式，索引可以获取约束子集的子集，只处理较小的子集，并且仍然明确地引用原始问题变量。

选择= optimoptions (“linprog”，“显示”，“没有”）;activeA = true(大小(b));activeAeq = true(大小(说真的));[~, ~, exitflag] = linprog (f, A、b Aeq,说真的,磅,乌兰巴托,选择);ncl = 1;而Exitflag < 0 [ineqselastic,eqselastic,ncall] =...elasticfilter ((activeA:), b (activeA) Aeq (activeAeq:),说真的(activeAeq)磅,乌兰巴托);NCL = NCL + ncall;idxaa =找到(activeA);idxae =找到(activeAeq);tmpa = idxaa(找到(ineqselastic));tmpae = idxae(找到(eqselastic));AA = (tmpa:);bb = b (tmpa);AE = Aeq (tmpae:);是= beq (tmpae); [ineqs,eqs,ncall] =...deletionfilter (AA、bb AE,磅,乌兰巴托);NCL = NCL + ncall;activeA (tmpa (ineqs)) = false;activeAeq (tmpae(方程式))= false;disp ([消除不平等的, int2str ((tmpa (ineqs))”),"平等"int2str ((tmpae(方程式 ))')]) [~,~, exitflag] =...linprog (f (activeA:)、b (activeA) Aeq (activeAeq:),说真的(activeAeq)磅,乌兰巴托,选择);NCL = NCL + 1;结束

删除不等式114和等式删除不等式97和等式删除不等式64 82和等式删除不等式60和等式

流(' linprog调用的数量:%d\n'ncl)

linprog调用次数:28次

注意，该循环同时删除不等式64和82，这表明这两个约束形成了IIS。

找到maxf

获得可行约束集的另一种方法是直接找到MaxFS。正如Chinneck[1]所解释的，寻找MaxFS是一个NP-complete问题，这意味着这个问题不一定有有效的算法来寻找MaxFS。然而，Chinneck提出了一些有效的算法。

使用Chinneck's算法7.3求a覆盖集约束的集合，当移除时，给出一个可行集。算法实现在发电机换向辅助函数在这个例子到此结束．

(coversetineq coverseteq, nlp) = generatecover (A, b, Aeq,说真的,磅,乌兰巴托)

coversetineq =5×1114 97 60 82 2

coverseteq = []

nlp=40

消除这些约束并解决LP。

usemeineq = true(大小(b));usemeineq (coversetineq) = false;去掉不等式约束usemeeq=真（尺寸（beq））；usemeeq（coverseteq）=假；删除等式约束[x, fvals exitflags] =...linprog (f (usemeineq:)、b (usemeineq) Aeq (usemeeq),说真的(usemeeq)磅,乌兰巴托);

找到了最优解。

请注意，封面和iiselastic设置从弹性过滤器。通常，弹性过滤器会发现太大的覆盖集。Chinneck的算法7.3从弹性过滤器结果开始，然后只保留必要的约束。

Chinneck的算法7.3需要40次调用linprog来完成MaxFS的计算。这个数字比之前在循环删除IIS过程中使用的28个调用多一点。

另外，注意循环中删除的不等式与算法7.3中删除的不等式并不完全相同。循环删除不等式114、97、82、60和64，而算法7.3消除了不等式114、97、82、60和2.检查不等式82和64是否构成IIS（如中所示在循环中删除IIS)，不等式82和不等式2也构成IIS。

usemeineq=false（size（b））；usemeineq（[82,64]）=true；ineqs=deletionfilter（A（usemeineq，：）、b（usemeineq）、Aeq、beq、lb、ub）；disp（ineqs）

１ １

usemeineq = false(大小(b));usemeineq([2] 82年)= true;ineqs = deletionfilter (A (usemeineq:)、b (usemeineq) Aeq,说真的,磅,乌兰巴托);disp (ineqs)

１ １

参考文献

[1] 金内克，J.W。优化的可行性和不可行性:算法和计算方法。施普林格,2008年。

j·W·Chinneck《优化的可行性和不可行性》CP-AI-OR-07教程，布鲁塞尔，比利时。可以在https://www.sce.carleton.ca/faculty/chinneck/docs/CPAIOR07InfeasibilityTutorial.pdf．

辅助函数

此代码创建deletionfilterhelper函数。

函数[ineq_iis,eq_iis,ncalls] = deletionfilter(Aineq,bineq,Aeq,beq,lb,ub) ncalls = 0; / /删除服务器[mi, n] =大小(Aineq);%变量数与线性不等式约束f=零（1，n）；me=尺寸（Aeq，1）；%线性等式约束数选择= optimoptions (“linprog”，“算法”，“对偶单纯形”，“显示”，“没有”); ineq_iis=真（mi，1）；%从问题中的所有不等式开始eq_iis=真（me，1）；%从问题中的所有等式开始为I =1:mi inq_iis (I) = 0;去掉不等式i[~, ~, exitflag] = linprog (f, Aineq (ineq_iis:), bineq (ineq_iis),...Aeq,说真的,磅,乌兰巴托,[],选择);Ncalls = Ncalls + 1;如果exitflag = = 1%如果现在可行ineq_iis (i) = 1;返回i到问题结束结束为I =1: eq_iis(I) = 0;删除等号i[~, ~, exitflag] = linprog (f Aineq bineq,...Aeq (eq_iis:),说真的(eq_iis)磅,乌兰巴托,[],选择);Ncalls = Ncalls + 1;如果exitflag = = 1%如果现在可行等式（i）=1；返回i到问题结束结束结束

此代码创建elasticfilterhelper函数。

函数[ineq_iis,eq_iis,ncalls,fval0] = elasticfilter(Aineq,bineq,Aeq,beq,lb,ub) ncalls = 0;[mi, n] =大小(Aineq);%变量数与线性不等式约束我=大小(Aeq, 1);Aineq_r = [Aineq -1.0*eye(mi) zeros(mi,2*me)];Aeq_r = [Aeq zeros(me,mi) eye(me) -1.0*eye(me)];%每个等式约束有两个松弛项lb_r=[lb（：）；零（mi+2*me，1）]；ub_r=[ub（：）；inf（mi+2*me，1）]；ineq_松弛_偏移=n；等式位置松弛偏移=n+mi；等式负松弛偏移=n+mi+me；f=[零（1，n）一（1，mi+2*me）]；opts=options(“linprog”，“算法”，“对偶单纯形”，“显示”，“没有”）;托尔= 1平台以及;ineq_iis = false (mi, 1);eq_iis = false(我,1);[x, fval exitflag] = linprog (f Aineq_r bineq Aeq_r,说真的,lb_r, ub_r,[],选择);fval0 = fval;Ncalls = Ncalls + 1;而exitflag==1&&fval>tol%可行，有些宽松是非零的c = 0;为I = 1:mi j = ineq_slack_offset+ I;如果X (j) > tol ub_r(j) = 0.0;ineq_iis (i) = true;c = c + 1;结束结束为i=1:me j=eq\u pos\u slack\u offset+i；如果X (j) > tol ub_r(j) = 0.0;eq_iis (i) = true;c = c + 1;结束结束为I = 1:me j = eq_neg_slack_offset+ I如果X (j) > tol ub_r(j) = 0.0;eq_iis (i) = true;c = c + 1;结束结束[x，fval，exitflag]=linprog（f，Aineq_r，bineq，Aeq_r，beq，lb_r，ub_r，[]，opts）；如果fval>0 fval0=fval；结束Ncalls = Ncalls + 1;结束结束

此代码创建发电机换向辅助函数。该代码使用了与代码相同的索引技术来跟踪约束在循环中删除IIS代码。

函数(coversetineq coverseteq, nlp) = generatecover (Aineq、bineq Aeq,说真的,磅,乌兰巴托)%返回线性不等式的覆盖集，线性的覆盖集% =，以及调用linprog的总次数。%改编自Chinneck[1]算法7.3。步骤编号来自这本书。coversetineq = [];coverseteq = [];activeA = true(大小(bineq));activeAeq = true(大小(说真的));%算法7.3的步骤1[ineq_iis、eq_iis、NCALL]=弹性过滤器（Aineq、bineq、Aeq、beq、lb、ub）；nlp=NCALL；求和(iis)和(iis)和(nineq):；如果Ninf == 1 coversetineq = inq_iis;coverseteq = eq_iis;返回结束holdsetineq=find（ineq_-iis）；holdseteq=find（eq_-iis）；candidateineq=holdsetineq；候选资格Q=holdseteq；%算法7.3的步骤2而Sum (candidateineq(:)) + Sum (candidateq (:)) > 0 minsinf = inf;ineqflag = 0;为i = 1:length(candidateineq(:)) activeA(candidateineq(i)) = false;idx2 =找到(activeA);idx2eq =找到(activeAeq);[ineq_iis, eq_iis ncalls fval] = elasticfilter (Aineq (activeA:), bineq (activeA) Aeq (activeAeq),说真的(activeAeq)磅,乌兰巴托);NLP = NLP + ncalls;ineq_iis = idx2(找到(ineq_iis));eq_iis = idx2eq(找到(eq_iis));如果Fval == 0 coversetineq = [coversetineq;candidateineq(i)];返回结束如果Fval < minsinf inqflag = 1;赢家= candidateineq(我);minsinf = fval;holdsetineq = ineq_iis;如果numel（ineq_iis（：）+numel（eq_iis（：）==1 nextwinner=ineq_iis；nextwinner2=eq_iis；nextwinner=[nextwinner，nextwinner2]；其他的nextwinner = [];结束结束activeA (candidateineq (i)) = true;结束为i = 1:length(candidateeq(:)) activeAeq(candidateeq(i)) = false;idx2 =找到(activeA);idx2eq =找到(activeAeq);[ineq_iis, eq_iis ncalls fval] = elasticfilter (Aineq (activeA) bineq (activeA) Aeq (activeAeq:),说真的(activeAeq)磅,乌兰巴托);NLP = NLP + ncalls;ineq_iis = idx2(找到(ineq_iis));eq_iis = idx2eq(找到(eq_iis));如果fval==0 coverseteq=[coverseteq；候选资格（i）]；返回结束如果Fval < minsinf inqflag = -1;赢家= candidateeq(我);minsinf = fval;holdseteq = eq_iis;如果numel（ineq_iis（：）+numel（eq_iis（：）==1 nextwinner=ineq_iis；nextwinner2=eq_iis；nextwinner=[nextwinner，nextwinner2]；其他的nextwinner = [];结束结束activeAeq (candidateeq (i)) = true;结束%算法7.3的步骤3如果== 1 coversetineq = [coversetineq;赢家];activeA(冠军)= false;如果Nextwinner coversetineq = [coversetineq; Nextwinner];返回结束结束如果Ineqflag == -1 coverseteq = [coverseteq;赢家];activeAeq(冠军)= false;如果Nextwinner coverseteq = [coverseteq; Nextwinner];返回结束结束candidateineq = holdsetineq;candidateeq = holdseteq;结束结束