我已经编写了以下代码enter image description here
如何在以下公式中解释辅助约束和子行程消除约束？enter image description here
[非常有效的模型是https://www.ibm.com/support/knowledgecenter/SSSA5P_12.8.0/ilog.odms.ide.help/examples/html/opl/models/TravelingSalesmanProblem/tsp.mod.html处OPL示例的一部分
但是一些用户出于许多不同的原因需要更多：
using CP;
int n = ...;
range Cities = 1..n;
int realCity[i in 1..n+1]=(i<=n)?i:1;
// Edges -- sparse set
tuple edge {int i; int j;}
setof(edge) Edges = {<i,j> | ordered i,j in 1..n};
setof(edge) Edges2 = {<i,j> | i,j in 1..n+1}; // node n+1 is node 1
int dist[Edges] = ...;
int dist2[<i,j> in Edges2]=(realCity[i]==realCity[j])?0:
((realCity[i]<realCity[j])?dist[<realCity[i],realCity[j]>]:dist[<realCity[j],realCity[i]>]);
dvar interval itvs[1..n+1] size 1;
dvar sequence seq in all(i in 1..n+1) itvs[i];
execute
{
cp.param.TimeLimit=60;
var f = cp.factory;
cp.setSearchPhases(f.searchPhase(seq));
}
tuple triplet { int c1; int c2; int d; };
{triplet} Dist = {
<i-1,j-1,dist2[<i ,j >]>
| i,j in 1..n+1};
minimize endOf(itvs[n+1]) - (n+1);
subject to
{
startOf(itvs[1])==0; // break sym
noOverlap(seq,Dist,true); // nooverlap with a distance matrix
last(seq, itvs[n+1]); // last node
}
int x[<i,j> in Edges]=prev(seq,itvs[i],itvs[j])+prev(seq,itvs[j],itvs[i]);
int isPrevFromNPlus1[i in 1..n]=prev(seq,itvs[i],itvs[n+1]);
int l=first({i | i in 1..n : isPrevFromNPlus1[i]==1});
edge el=<1,l>;
execute
{
isPrevFromNPlus1;
x;
x[el]=1;
}
// Let us check here that the constraints of the IP model are ok
assert forall (j in Cities)
as:sum (<i,j> in Edges) x[<i,j>] + sum (<j,k> in Edges) x[<j,k>] == 2;
// Let us compute here the objective the IP way
int cost=sum (<i,j> in Edges) dist[<i,j>]*x[<i,j>];
execute
{
writeln(cost);
}
让我在这里提供更多选项，我将依赖与tsp示例相同的.dat格式1）我们可以直接移除所有电路，但这不是很有效，可以检查：
// Cities int n = ...; range Cities = 1..n; // Edges -- sparse set tuple edge {int i; int j;} setof(edge) Edges = {<i,j> | ordered i,j in Cities}; int dist[Edges] = ...; // Decision variables dvar boolean x[Edges]; {int} nodes={i.i | i in Edges} union {i.j | i in Edges}; range r=1..-2+ftoi(pow(2,card(nodes))); {int} nodes2 [k in r] = {i | i in nodes: ((k div (ftoi(pow(2,(ord(nodes,i))))) mod 2) == 1)}; /***************************************************************************** * * MODEL * *****************************************************************************/ // Objective minimize sum (<i,j> in Edges) dist[<i,j>]*x[<i,j>]; subject to { // Each city is linked with two other cities forall (j in Cities) sum (<i,j> in Edges) x[<i,j>] + sum (<j,k> in Edges) x[<j,k>] == 2; // Subtour elimination constraints. forall(k in r) // all subsets but empty and all sum(e in Edges:(e.i in nodes2[k]) && (e.j in nodes2[k])) x[e]<=card(nodes2[k])-1; }
2）更好且依赖于CPLEX的是MTZ模型（Miller-Tucker-Zemlin公式）
// Cities int n = ...; range Cities = 1..n; // Edges -- sparse set tuple edge {int i; int j;} setof(edge) Edges = {<i,j> | ordered i,j in Cities}; int dist[Edges] = ...; setof(edge) Edges2 = {<i,j> | i,j in Cities : i!=j}; int dist2[<i,j> in Edges2] = (<i,j> in Edges)?dist[<i,j>]:dist[<j,i>]; // Decision variables dvar boolean x[Edges2]; dvar int u[1..n] in 1..n; /***************************************************************************** * * MODEL * *****************************************************************************/ // Objective minimize sum (<i,j> in Edges2) dist2[<i,j>]*x[<i,j>]; subject to { // Each city is linked with two other cities forall (j in Cities) { sum (<i,j> in Edges2) x[<i,j>]==1; sum (<j,k> in Edges2) x[<j,k>] == 1; } // MTZ u[1]==1; forall(i in 2..n) 2<=u[i]<=n; forall(e in Edges2:e.i!=1 && e.j!=1) (u[e.j]-u[e.i])+1<=(n-1)*(1-x[e]); }; {edge} solution={e | e in Edges2 : x[e]==1}; execute { writeln("path ",solution); }
3）没有计划的CPO
using CP; int n = ...; range Cities = 1..n; int realCity[i in 1..n+1]=(i<=n)?i:1; // Edges -- sparse set tuple edge {int i; int j;} setof(edge) Edges = {<i,j> | ordered i,j in 1..n}; setof(edge) Edges2 = {<i,j> | i,j in 1..n}; // node n+1 is node 1 int dist[Edges] = ...; int dist2[i in 1..n][j in 1..n] = (i==j)?0:((i<j)?dist[<i,j>]:dist[<j,i>]); execute { cp.param.TimeLimit=60; } dvar int x[1..n] in 1..n; dvar int obj; minimize obj; // x means who is on i th position subject to { x[1]==1; allDifferent(x); obj==sum(i in 1..n-1) dist2[x[i]][x[i+1]]+dist2[x[1]][x[n]]; }
u(i)分配给每个节点。约束u(i) - u(j) + n x(i,j) <= n - 1
可以解释为
u(j) >= u(i) + 1 - M*(1-x(i,j))
或
x(i,j) = 1 ==> u(j) >= u(i) + 1
[如果您有子路线，请说2-3-2，这将无法成立：
2-3 means u(3) >= u(2)+1 3-2 means u(2) >= u(3)+1
因为应该允许整个游览，所以我们只是从此检查中删除第一个节点。