广东哪家网站建设哪家公司好,会员网站建设,wordpress 模块开发教程,长沙flash网站制作1. 大规模tsp问题的挑战
数学模型和精确解法见《运筹系列65#xff1a;TSP问题的精确求解法概述》和《运筹系列80:使用Julia精确求解tsp问题》#xff1a;
variable(m, x[1:n,1:n], Bin,Symmetric) # 0-1约束
objective(model, Min, sum(x.*distmat)/2)
constraint(model, …1. 大规模tsp问题的挑战
数学模型和精确解法见《运筹系列65TSP问题的精确求解法概述》和《运筹系列80:使用Julia精确求解tsp问题》
variable(m, x[1:n,1:n], Bin,Symmetric) # 0-1约束
objective(model, Min, sum(x.*distmat)/2)
constraint(model, [i 1:n], sum(x[i, :]) 2)
constraint(model, [j 1:n], sum(x[:, j]) 2)
for subset in subsets # 防subtour约束需要遍历所有的子集合constraint(model,sum(x[subset[i],subset[i]])2*length(subset[i])-2)
end
optimize!(m)
objective_value(model)主要的挑战包括
求解整数规划比较耗时并且遍历子集难度很大大规模问题时需要自定义高效的cut和price策略变量数量过多需要使用列生成方法来减少求解变量
2. 使用列生成减少求解变量
关于列生成的列子可以参考《运筹系列8Set partitioning问题的列生成方法》。对于tsp问题我们使用列生成的步骤是
首先将问题变量限制在每个点最邻近的k条边上求解限制主问题迭代求解子问题sub problem。如果目标函数最优值0就将新生成的列yk1转入基变量生成新的限制主问题进行求解。如此往复直至子问题的目标函数≥0。
我们以u159问题为例原始模型一共有25281个变量
using TSPLIB, LinearAlgebra, JuMP,HiGHS,PyPlot
#读取数据
tsp readTSPLIB(:u159)
n tsp.dimension
distmat [round(Int64,norm(tsp.nodes[i,:] - tsp.nodes[j,:])) for i in 1:n, j in 1:n];
for i in 1:n;distmat[i,i] 10^9;end#完整模型
model Model(HiGHS.Optimizer)
set_silent(model)
variable(model, 0x[i in 1:n,j in js[i]]1)
info variable size:,length(x)
objective(model, Min, sum([x[i,j]*distmat[i,j]/2 for i in 1:n for j in js[i]]))
c constraint(model, [i 1:n], sum([x[i, j] for j in js[i]]) 2)
e constraint(model, [i 1:n], sum([x[j, i] for j in js[i]]) 2)
optimize!(model)
info objective:,objective_value(model)主问题限制在每个点最近的3条边中此时只有604个变量第一轮列生成后为686个变量此时已经得到最优解了最终经过7轮列生成迭代变量个数为716
# js为初始下标合集我们将变量下标限制在这个集合内。注意需要保持对称性
idx Index(2)
add(idx, tsp.nodes)
c search(idx,tsp.nodes,4)[2];
js Vector{Vector{Int}}()
for i in 1:npush!(js,c[i,2:end])
end
for i in 1:nfor j in c[i,2:end]union!(js[j],i)end
end求解主问题然后找到检验数0的列。令对偶变量 p C B B − 1 pC_B B^{-1} pCBB−1, 则检验数 σ C N − p ∗ N \sigma C_N-p*N σCN−p∗N 。我们将所有检验数小于0的列进基迭代直至an 0此部分完整代码如下
an 1
while an 0model Model(HiGHS.Optimizer)set_silent(model)variable(model, 0x[i in 1:n,j in js[i]]1)info variable size:,length(x)objective(model, Min, sum([x[i,j]*distmat[i,j]/2 for i in 1:n for j in js[i]]))c constraint(model, [i 1:n], sum([x[i, j] for j in js[i]]) 2)e constraint(model, [i 1:n], sum([x[j, i] for j in js[i]]) 2)optimize!(model)info objective_value(model)p1 [dual(c[i]) for i in 1:n]p2 [dual(e[i]) for i in 1:n];an 0for i in 1:nfor j in setdiff(1:n,js[i])sigma distmat[i,j]/2 - p1[i]-p2[j]if sigma 0;push!(js[i],j);push!(js[j],i);an2;endendendinfo add size:,an
end3. 使用行生成添加约束
我们观察att48使用lp求解后的结果
有很多子环需要添加去除环的约束。 首先添加subtour约束去除所有的子环约束
using Graphs
paths strongly_connected_components(SimpleDiGraph(round.(x_val,digits 2)))
if length(paths)1for path in pathsconstraint(model,sum(x[path,path])2*length(path)-2)endoptimize!(model)x_val round.(JuMP.value.(x),digits 2)paths strongly_connected_components(SimpleDiGraph(x_val))
end迭代求解并绘制结果如下
然后添加comb约束约束介绍参见《运筹系列65TSP问题的精确求解法概述》julia代码参见《运筹系列80: 使用Julia精确求解tsp问题》的4.2节求解后得到两个comb 然后添加进约束迭代求解 迭代得到如下图
4. 使用分枝定界求解整数规划问题
我们这里使用priority queue存储分枝节点按照最简单的下标顺序对所有非整数变量进行分枝。