使用钢退火算法解决斯坦纳问题。

引言

首先，我们来讨论什么是斯坦纳问题。我们有一个包含 n 个点的平面。我们可以随意添加 m 个点。之后，我们应该构建最小的生成树。该任务的问题在于如何放置 m 个点，以便获得最小的生成树。

那么，我们为什么不使用简单的暴力算法呢？因为它太慢了。 

爬山算法又如何呢？不幸的是，我们试图优化的函数有太多的局部最小值，所以爬山算法会停留在其中一个最小值上。

钢退火算法是一种用于寻找全局最小值的启发式算法。退火算法的主要思想与爬山算法相同：如果新解优于现有解，则接受它。但是，退火算法有时会接受不是最佳的解。这种“愚蠢”的决策可以避免陷入局部最小值陷阱。

例如，让我们看一下以下函数

 

如果我们从红色点开始： 

我们可能会停留在 (0,0) 点。为了达到真正的最小值，我们需要移动到大于局部最小值的位置。

因此，在本文中，我将描述如何使用退火算法实现斯坦纳问题的解决方案。我使用的编程语言是 C#。 

生成树长度计算的实现。

首先，让我们定义要最小化的函数。我们要最小化的函数是生成树中所有边的长度之和。为了构建生成树，我们将使用Prim算法。

Prim 算法的思想非常简单

1) 从树中取任意一个点，并将其添加到选定点的集合中。

2)  添加连接已选点和未选点的最轻边。将未选点添加到已选点集。

3) 如果未选点集不为空，则转到 2。

4) 完成。

以下代码说明了我对该算法的实现。 

Pair FindTheSmallestEdgeAndInitializeSets(HashSet<int> unUsed, HashSet<int> used, int problemNumber)
<pre>{
    double bestValue = double.MaxValue;
    Pair result = new Pair();
    for (int i = 0; i < data.points[problemNumber].Length; ++i)
    {
        for (int j = 0; j < data.points[problemNumber].Length; ++j)
        {
            if (bestValue > distances_[i, j] && (i != j))
            {
                bestValue = distances_[i, j];
                result.first = i;
                result.second = j;
            }
        }
    }
    used.Add(result.first);
    used.Add(result.second);
    unUsed.Remove(result.first);
    unUsed.Remove(result.second);
    return result;
}
private Pair FindBestEdgeToAdd(HashSet<int> unUsed, HashSet<int> used, int problemNumber)
{
    double bestValue = double.MaxValue;
    Pair result = new Pair();
    foreach (int i in used)
    {
        foreach (int j in unUsed)
        {
            if (bestValue > distances_[i, j])
            {
                bestValue = distances_[i, j];
                result.first = i;
                result.second = j;
            }
        }
    }
    unUsed.Remove(result.second);
    used.Add(result.second);
    return result;
}
private void PrimeAlgorithm(int problemNumber)
{
    HashSet<int> usedPoints = new HashSet<int>();
    HashSet<int> unUsedPoints = new HashSet<int>();
    initializeSet(unUsedPoints, data.points[problemNumber].Length);
    int edgesCount = 0;
    edges = new Pair[data.points[problemNumber].Length - 1];
    edges[edgesCount++] = FindTheSmallestEdgeAndInitializeSets(unUsedPoints, usedPoints, problemNumber);
    while (unUsedPoints.Count != 0)
    {
        edges[edgesCount++] = FindBestEdgeToAdd(unUsedPoints, usedPoints, problemNumber); 
    }
} 
		A brief description of how to use the article or code. The
		class names, the methods and properties, any tricks or tips.

退火算法实现。

正如前面提到的，我们并不总是选择最佳结果。让我们了解一下什么时候我们会选择一个坏的解决方案。“退火算法”这个名字的由来是因为该算法使用了钢冷却的思想。当金属冷却时，它的粒子也在冷却，但有时单个粒子的温度会升高。金属温度越高，粒子变热的可能性就越高。

因此，我们将有一些初始温度。每次迭代后，温度会越来越低。因此，接受坏结果的可能性会降低。

以下是我对退火算法的实现

private void Annealing(int problemNumber)
{
    double alpha = 0.99;
    double temperature = 100.0;
    double initialTemperature = temperature;
    double epsilon = 0.2;
    double deltha;
    double accepted = 0;
    double rejected = 0;
    double multiplicator = 0.05;
    best = GetTreeLength(problemNumber);
    Point[] bestPoints = new Point[data.points[problemNumber].Length];
    Pair[] bestEdges = new Pair[edges.Length];
    Array.Copy(data.points[problemNumber], bestPoints, bestPoints.Length);
    Array.Copy(edges, bestEdges, edges.Length);
    Random rand = new Random();
    
    while (temperature > epsilon)
    {
        double prevLength = GetTreeLength(problemNumber);
        double[,] offsets = new double[2, mutablePointsNumber]; 
        for (int i = 0; i < mutablePointsNumber; ++i)
        {
            offsets[0,i] = rand.NextDouble() * multiplicator - 0.5 * multiplicator;
            offsets[1,i] = rand.NextDouble() * multiplicator - 0.5 * multiplicator;
            int currentIndex = data.points[problemNumber].Length - i - 1;
            ModifyPoint(currentIndex, problemNumber, offsets[0, i], offsets[1, i]);
        }
            
            
        double currentLength = GetTreeLength(problemNumber);
        double proba = rand.NextDouble();
        deltha = currentLength - prevLength;
        double e = Math.Exp((-deltha * 13 * initialTemperature) / temperature);
        if (proba < e || deltha < 0)
        {
            accepted++;
            if (currentLength < best)
            {
                best = currentLength;
                Array.Copy(data.points[problemNumber], bestPoints, bestPoints.Length);
                Array.Copy(edges, bestEdges, edges.Length);
            }
        }
        else
        {
            rejected++;
            for (int i = 0; i < mutablePointsNumber; ++i)
            {
                int currentIndex = data.points[problemNumber].Length - i - 1;
                ModifyPoint(currentIndex, problemNumber, -offsets[0, i], -offsets[1, i]);
            }
        }
        
        chart1.Series[0].Points.AddY(currentLength);
        if (accepted > mutablePointsNumber * 2  || rejected > 4 * mutablePointsNumber )
        {
            accepted = 0;
            rejected = 0;
            temperature *= alpha;
        }
    }
    Array.Copy(bestPoints, data.points[problemNumber], bestPoints.Length);
    Array.Copy(bestEdges, edges, edges.Length);
} 

 

结果。

为了说明结果，我使用了以下图形 

这张图说明了算法迭代过程中长度的变化。在右下角，我们可以看到初始长度和结果长度之间的差异。

为了说明这些点的位置，我做了如下可视化。

 

该程序的代码和测试数据可以从下载源代码。