旅行商问题
旅行推销员问题（英语：Travelling salesman problem, TSP）是这样一个问题：给定一系列城市和每对城市之间的距离，求解访问每一座城市一次并回到起始城市的最短回路。它是组合优化中的一个NP难问题，在运筹学和理论计算机科学中非常重要。
——百度百科
模拟退火算法
模拟退火算法（Simulate Anneal Arithmetic，SAA）是一种通用概率演算法，用来在一个大的搜寻空间内找寻命题的最优解。模拟退火是S.Kirkpatrick, C.D.Gelatt和M.P.Vecchi在1983年所发明。而V.Černý在1985年也独立发明此演算法。模拟退火算法是解决TSP问题的有效方法之一。
——百度百科
效果展示
初始解
最终解
代码实现：
SAPoint
public class SAPoint
{
public double X { get;set; }
public double Y { get;set; }
public SAPoint()
{
}
public SAPoint(double _x,double _y)
{
X=_x;
Y=_y;
}
/// <summary>
/// 计算欧式距离
/// </summary>
/// <param name="point"></param>
/// <returns></returns>
public double DistanceTo(SAPoint point)
{
return Math.Sqrt((point.X - X) * (point.X - X) + (point.Y - Y) * (point.Y - Y));
}
}SimulatedAnnealing
public class SimulatedAnnealing
{
/// <summary>
/// 城市坐标
/// </summary>
public List<SAPoint> points;
/// <summary>
/// 温度下降参数
/// </summary>
public double a = 0.99;
/// <summary>
/// 初始温度
/// </summary>
public double t0 = 97;
/// <summary>
/// 最低温度
/// </summary>
public double tf = 3;
/// <summary>
/// Boltzmann常数
/// </summary>
public double Kb = 1;
/// <summary>
/// 马尔可夫链长度
/// </summary>
public double Markov_length = 1000;
/// <summary>
/// 当前解
/// </summary>
public List<int> sol_current = new List<int>();
/// <summary>
/// 最好解
/// </summary>
public List<int> sol_best = new List<int>();
/// <summary>
/// i到j的距离
/// </summary>
private double[][] distMap;
/// <summary>
/// 解决旅行商问题
/// </summary>
public void Run()
{
// 计算每个点的距离
InitDist();
// 产生初始温度
double t = t0;
// 产生初始解
List<int> sol_new = new List<int>();
for (int i = 0; i < points.Count; i++)
{
sol_new.Add(i);
}
int amount = sol_new.Count;
Draw(sol_new, "new");
// 当前解对应的回路距离
double E_current = double.MaxValue;
// 最优解对应的回路距离
double E_best =double.MaxValue;
double E_new;
sol_current.Clear();
sol_current.AddRange(sol_new);
sol_best.Clear();
sol_best.AddRange(sol_new);
Random rand = new Random();
while (t >= tf)
{
for (int r = 1; r < Markov_length; r++)
{
// 产生随机扰动（产生新解）
// 随机决定是进行两交换还是三交换
if (rand.NextDouble() < 0.5)
{// 两交换-u,v交换位置
int ind1 = 0;
int ind2 = 0;
while (ind1 == ind2)
{
ind1 = rand.Next(amount);
ind2 = rand.Next(amount);
}
int temp = sol_new[ind1];
sol_new[ind1] = sol_new[ind2];
sol_new[ind2] = temp;
}
else
{// 三交换 u,v,w；uv之间的路径放到w后面
int ind1 = 0;
int ind2 = 0;
int ind3 = 0;
while (ind1 == ind2|| ind1 == ind3|| ind2 == ind3 || Math.Abs(ind1-ind2)==1)
{
ind1 = rand.Next(amount);
ind2 = rand.Next(amount);
ind3 = rand.Next(amount);
}
// 确保ind1<ind2<ind3
if (ind1 > ind2) Swap(ref ind1, ref ind2);
if (ind1 > ind3) Swap(ref ind1, ref ind3);
if (ind2 > ind3) Swap(ref ind2, ref ind3);
// 交换
List<int> tempList = new List<int>();
for (int i = ind1 + 1; i <= ind2 - 1; i++)
{
tempList.Add(sol_new[i]);
}
for (int i = ind1 + 1, j = ind2; i <= ind1 + ind3 - ind2 + 1; i++, j++)
{
sol_new[i] = sol_new[j];
}
for (int i = ind1 + ind3 - ind2 + 2, j = 0; i <= ind3; i++, j++)
{
sol_new[i] = tempList[j];
}
}
// 计算内能（目标函数）
E_new = 0;
for (int i = 0; i < amount; i++)
{
E_new += distMap[sol_new[i]][sol_new[(i + 1) % amount]];
}
// 检查是否满足约束
if(E_new < E_current)
{
E_current = E_new;
sol_current.Clear();
sol_current.AddRange(sol_new);
if (E_new < E_best)
{
// 把冷却过程中最好的解保存下来
E_best = E_new;
sol_best.Clear();
sol_best.AddRange(sol_new);
}
}
else
{
// 若新解的目标函数小于当前解，则仅以一定概率接受新解
if (rand.NextDouble() < Math.Exp(-(E_new - E_current) / (Kb * t)))
{
E_current = E_new;
sol_current.Clear();
sol_current.AddRange(sol_new);
}
else
{
sol_new.Clear();
sol_new.AddRange(sol_current);
}
}
}
// 控制温度t减少为原来的a倍
t *= a;
}
Draw(sol_best, "best"+E_best);
Draw(sol_current, "current"+E_current);
}
/// <summary>
/// 交换ab的值
/// </summary>
/// <typeparam name="T"></typeparam>
/// <param name="a"></param>
/// <param name="b"></param>
private void Swap<T>(ref T a, ref T b)
{
T temp = a;
a = b;
b = temp;
}
/// <summary>
/// 初始化距离矩阵
/// </summary>
private void InitDist()
{
distMap = new double[points.Count][];
for (int i = 0;i<points.Count; i++)
{
distMap[i] = new double[points.Count];
}
for (int i = 0; i < points.Count - 1; i++)
{
for (int j = i+1; j < points.Count; j++)
{
double dist = i == j ? 0 : points[i].DistanceTo(points[j]);
distMap[i][j] = dist;
distMap[j][i] = dist;
}
}
}
/// <summary>
/// 可视化-输出图片
/// </summary>
/// <param name="sol"></param>
/// <param name="name"></param>
private void Draw(List<int> sol,string name = "sol")
{
MapData mapData = new MapData(1);
List<Line> lines = PointsToLines(sol);
mapData.CreateMap(mapData.GetMaxBound(lines));
mapData.AddLines(lines, (int)AstarType.Red);
mapData.MapToBitmap(mapData.mapData, "E://" + name + ".bmp");
}
/// <summary>
/// 点转线
/// </summary>
/// <param name="sol"></param>
/// <returns></returns>
private List<Line> PointsToLines(List<int> sol)
{
int k = 10;
List<Line> lines = new List<Line>();
for (int i = 0; i < sol.Count; i++)
{
SAPoint p1 = points[sol[i]];
SAPoint p2 = points[sol[(i+1)%sol.Count]];
lines.Add(new Line(new XYZ(p1.X * k, p1.Y * k, 0), new XYZ(p2.X * k, p2.Y * k, 0)));
}
return lines;
}
}其中画图的方法参考我之前的代码：
测试方法：
/// <summary>
/// 测试模拟退火算法
/// </summary>
public void Test()
{
List<SAPoint> points = new List<SAPoint>()
{
new SAPoint(37,52),
new SAPoint(49,49),
new SAPoint(52,64),
new SAPoint(20,26),
new SAPoint(40,30),
new SAPoint(21,47),
new SAPoint(17,63),
new SAPoint(31,62),
new SAPoint(52,33),
new SAPoint(51,21),
new SAPoint(42,41),
new SAPoint(31,32),
new SAPoint(5,25),
new SAPoint(12,42),
new SAPoint(36,16),
new SAPoint(52,41),
new SAPoint(27,23),
new SAPoint(17,33),
new SAPoint(13,13),
new SAPoint(57,58),
new SAPoint(62,42),
new SAPoint(42,57),
new SAPoint(16,57),
new SAPoint(8,52),
new SAPoint(7,38),
new SAPoint(27,68),
new SAPoint(30,48),
new SAPoint(43,67),
new SAPoint(58,48),
new SAPoint(58,27),
new SAPoint(37,69),
new SAPoint(38,46),
new SAPoint(46,10),
new SAPoint(61,33),
new SAPoint(62,63),
new SAPoint(63,69),
new SAPoint(32,22),
new SAPoint(45,35),
new SAPoint(59,15),
new SAPoint(5,7),
new SAPoint(10,17),
new SAPoint(21,10),
new SAPoint(5,64),
new SAPoint(30,15),
new SAPoint(39,10),
new SAPoint(32,39),
new SAPoint(25,32),
new SAPoint(25,55),
new SAPoint(48,28),
new SAPoint(56,37),
new SAPoint(30,40),
};
SimulatedAnnealing sa = new SimulatedAnnealing()
{
points = points,
Markov_length = 2000
};
sa.Run();
}我的代码仓库
Algorithm: 简单算法小仓库https://gitee.com/imaginationWdq/algorithm
