簡述

算法設計課這周的作業：
趕緊寫了先，不然搞不完了。

文章目錄

• • 簡述
  • 算法理論部分
  • • 變量簡單分析
    • 從狀態轉移概率到狀態概率
    • 推導
    • 理解當溫度收斂到接近0的時候，收斂到結果
    • 理論部分的后記
  • python實現
  • • python實現模擬退火解極值
  • C++實現解函數值問題
  • TSP問題求解
  • • 代碼詳細解釋
    • • 定義了兩個宏，主要是為了加速運算
      • 定義全局變量
      • 函數解釋
      • 給一組數據和結果

算法理論部分

• 用粒子的排列或相應的能量表示物體所處的狀態，在溫度T下，物體(系統)所處的狀態具有一定的隨機性。主流趨勢是系統向能量較低的狀態發展，但粒子的不規則熱運動妨礙系統準確落入低能狀態。

簡單來說就是，在溫度T下，

有兩種可能

• $f(x_i)\ge f(x_j)$， 那沒得說有更好的肯定選更好的。
• $f(x_i)<f(x_j)$，這個時候我們也有一定的概率選這個。

概率為

$$P(X_{i\to j})=e^{\frac{f(x_i)?f(x_j)}{KT}}$$

其中K為玻爾茲曼常數（在這個算法上，我們經常使用數值1代替這個）

• $K>0$，一般設置為$K=1$
• $T>0$，我們在遞減的過程中，終止的T，一般認為是$T=1$

變量簡單分析

不難看出，隨著T下降，這個狀態轉移的概率是下降的。

• 原因：$f(x_i)?f(x_j)<0$

物理上解釋：

• 因為之前的這個轉移概率，認為是，在高溫情況下，分子可能會產生較為不穩定的隨機運動。且溫度越高，這個分子不規則運動的可能是更大的。這是在模擬這個過程。

• 所以，我們稱這個算法為，模擬退火算法

從狀態轉移概率到狀態概率

這個分析過程，其實是類似于推理馬爾可夫鏈的過程。

• 之前給給出的$P(X_{i\to j})$其實是條件概率。

我們假設第n個狀態為$X_i$，那么現在就是考慮下一個狀態的概率。

$$P(s_{n+1}=X_j|s_n=X_i)=P(X_{i\to j})=e^{\frac{f(x_i)?f(x_j)}{KT}}$$

我們用$s_n$，表示第n個狀態。

用條件概率公式得到

$$P(s_{n+1}=X_j|s_n=X_i)=\frac{P(s_{n+1}=X_j,s_n=X_i)}{P(s_n=X_i)}$$

$$P(s_{n+1}=X_j,s_n=X_i)=P(s_n=X_i)?P(s_{n+1}=X_j|s_n=X_i)$$
$$P(s_{n+1}=X_j)=\sum _{i}P(s_n=X_i)?P(s_{n+1}=X_j|s_n=X_i)$$

而解這個過程，就用到了以前解馬爾科夫過程的方法，這里我需要假設一個均衡態，在這種狀態下，經過任意次狀態轉移之后，整個模型的概率保持一致。

得到結果是

$$P_i(T)=\frac{e^{\frac{?f(x_i)}{KT}}}{Z_T}$$

$Z_T$只是一個歸一化的因子而已，就是把所有的這樣的分子的數值求個和。使得概率和為一而已。

這個分布稱之為 Boltzmann分布。

推導

$$P_i(T)?P_j(T)=\frac{e^{\frac{?f(i)}{KT}}}{Z_T}(1?e^{?\frac{f(j)?f(i)}{KT}})$$

理解當溫度收斂到接近0的時候，收斂到結果

其實只需要理解下面這個函數的導數就可以了

$$e^{?\frac{E(i)}{KT}}$$
關于T求導。

$$e^{?\frac{E(i)}{KT}}?\frac{E(i)}{k{T}^2}$$

那么，這就是這個變量關于T的變化導數。再考慮關于函數值$E(i)$

這個函數是關于$E(i)$的單調減函數。$E(i)>0$的情況下。

所以說，函數值稍微小的數值所對應的狀態概率，受到T的減小而導致的減小的幅度，其實是較為小的。
那么當T趨于0的時候，就可以得到，狀態概率，會集中在函數值較為小的數值點上（數值小的點的概率大）

也就是說，當模擬退火，溫度越低，越有可能收斂到正確解。而且，這是收斂的。

理論部分的后記

這里，我們證明了，當溫度越小，越會收斂到正確解。這只是在理論上的證明。但是我們都知道，當T趨于0，但是特別小的時候，作為分母時，這個精度就會變得非常低了（計算機上是離散的）。

• 所以，為了簡單，我們一般設置T到1就截止了。

python實現

先嘗試解決下面的這個圖

import numpy as np import matplotlib.pyplot as pltf = lambda x: (x - 1) ** 2 + x + 100 * np.sin(x)**2 x = np.linspace(-10., 10, 101) plt.plot(x, f(x)) plt.show()

python實現模擬退火解極值

• 先在這個區間上隨機生成一個起始點
• 每次在給定點的一個鄰近區域（自己設置），找到一個新的點。
• 然后這兩個點做比較。通過概率模擬來確定是否發生位置遷移。
• 直到溫度下降到一定的數值。

import numpy as np import matplotlib.pyplot as pltf = lambda x: (x - 1) ** 2 + x + 100 * np.sin(x) ** 2 x = np.linspace(-10., 10, 101) plt.plot(x, f(x))T = 1000 # 起始溫度 alpha = 0.9 # T_{k+1} = alpha * T_k方式更新溫度 limitedT = 1. # 最小值的T iterTime = 1000 # 每個溫度下迭代的次數 sigmaX = 3 # 每次的搜索半徑 K = 1 # 系數K X = 20 * np.random.rand() - 10 Y = f(X) while T > limitedT:for i in range(iterTime):xnew = X + sigmaX * (2 * np.random.rand() - 1)if xnew <= 10. and xnew >= -10.:fnew = f(xnew)if fnew < Y:Y = fnewX = xnewelse:res = Y-fnewp = np.exp(res / (K * T))if np.random.rand() < p:X = xnewY = fnewT *= alpha print(X, Y) plt.plot(X, Y, '*r') plt.show()

基本上準確~

C++實現解函數值問題

對于上面一個問題的C++解法
但是C++的精度似乎沒有Python的好。這個可以適當調節參數來獲得更高的精度，反正C++的速度也快很多。

#include <iostream> #include <cmath> #include <ctime> using namespace std; #define FUN(x) ((x-1)*(x-1) + x + 100 * sin(x) * sin(x)) #define RAND() ( (double)rand() / double(RAND_MAX))int main() {srand((unsigned)time(NULL));double T = 1000; // 起始溫度double alpha = 0.99; // T_{ k + 1 } = alpha * T_k方式更新溫度double limitedT = 1; // 最小值的Tint iterTime = 1000; // 每個溫度下迭代的次數double sigmaX = 3; // 每次的搜索半徑double K = 1; // 系數Kdouble X = 20 * RAND() - 10;double Y = FUN(X);double xnew, ynew, rest, p;while (T > limitedT) {for (int i = 0; i < iterTime; ++i) {xnew = X + sigmaX * (2 * RAND() - 1);if (xnew >= -10. && xnew <= 10) {ynew = FUN(xnew);if (ynew <= Y) {X = xnew;Y = ynew;}else {rest = Y-ynew;p = exp(rest / (K*T));if (RAND() < p) {X = xnew;Y = ynew;}}}}T *= alpha;}cout << X << " " << Y << endl;system("pause"); }

TSP問題求解

---

代碼詳細解釋

---

定義了兩個宏，主要是為了加速運算

• 第一個是RAND(b,e) 生成在$[b,e)$這個區間上的整數
• 第二個是RANDFLOAT() 是為了生成(0,1)之間的數。

---

定義全局變量

double **Mat; int *Path, *tempPath; double Value, tempValue; int N = 0;

• Mat是來存儲距離矩陣的
• Path存儲路徑
• Value存儲路勁所對應的函數。這里考慮的是具體的數值
• N表示有多少個點
• 所有前面加了temp都表示臨時變量。

---

函數解釋

• double CalValue(int *p) 函數用于給輸入的路勁下，求對應的value值。
• void refresh() 將temp的數組和具體數值恢復。為了下一次的考慮
• void change() 覆蓋掉原來的路勁。進行存儲。
• void initialPath() 初始化路徑。并算出初始值。

#include <iostream> #include <cmath> #include <ctime> #include <fstream> using namespace std; #define RAND(b, e) (rand() % (e-b) + b) // 左閉右開 #define RANDFLOAT() ((double)rand() / double(RAND_MAX)) // 0-1浮點數double **Mat; int *Path, *tempPath; double Value, tempValue; int N = 0;double CalValue(int *p) {double t = 0;for (int i = 1; i < N; ++i) {t += Mat[p[i - 1]][p[i]];}t += Mat[p[N - 1]][0];return t; }// 將temp重置為path void refresh() {for (int i = 0; i < N; ++i) {tempPath[i] = Path[i];}tempValue = Value; }// 覆蓋修改Path void change() {for (int i = 0; i < N; ++i) {Path[i] = tempPath[i];}Value = tempValue; }void initialPath() {for (int i = 0; i < N; ++i) {Path[i] = i;}// Path[i]表示路上第i個點的標記為Path[i]// Path[0] = 0int tx,t;for (int i = 1; i < N-1; ++i) {tx = RAND(i, N);if (tx != i) { // swapt = Path[i];Path[i] = Path[tx];Path[tx] = t;}}Value = CalValue(Path); }int main() {srand((unsigned)time(NULL));ifstream cin("data.txt");cin >> N;// initializeMat = new double *[N];for (int i = 0; i < N; ++i)Mat[i] = new double[N];for (int i = 0; i < N; ++i) {for (int j = 0; j < N; ++j) {cin >> Mat[i][j];}}Path = new int[N];tempPath = new int[N];double T = 1000; // 起始溫度double alpha = 0.99; // T_{ k + 1 } = alpha * T_k方式更新溫度double limitedT = 1e-9; // 最小值的Tint iterTime = 1000; // 每個溫度下迭代的次數double K = 1; // 系數Kdouble p = 0;initialPath();int tx, ty, t;while (T >= limitedT) {for (int i = 0; i < iterTime; ++i) {// 任意交換兩點的順序refresh();tx = RAND(1, N);ty = RAND(1, N);if (tx != ty) {t = tempPath[tx];tempPath[tx] = tempPath[ty];tempPath[ty] = t;tempValue = CalValue(tempPath);if (tempValue <= Value) {change();}else {p = exp((Value-tempValue) / (K*T));if (RANDFLOAT() < p) { change(); }}}}T *= alpha;}cout << "Value:" << Value << endl;for (int i = 0; i < N; ++i) {cout << Path[i] << " -> ";}cout << 0 << endl;// releasedelete[]tempPath;delete[] Path;for (int i = 0; i < N; ++i)delete[] Mat[i];delete[] Mat;system("pause"); }

給一組數據和結果

10 0 58 82 89 17 50 26 48 70 19 58 0 74 46 70 2 70 49 87 60 82 74 0 58 76 98 37 97 34 67 89 46 58 0 15 17 28 69 46 79 17 70 76 15 0 98 60 69 97 89 50 2 98 17 98 0 81 14 43 47 26 70 37 28 60 81 0 43 73 56 48 49 97 69 69 14 43 0 39 0 70 87 34 46 97 43 73 39 0 53 19 60 67 79 89 47 56 0 53 0

輸出：

Value:244 0 -> 4 -> 3 -> 6 -> 2 -> 8 -> 5 -> 1 -> 7 -> 9 -> 0

總結

以上是生活随笔為你收集整理的模拟退火算法理论+Python解决函数极值+C++实现解决TSP问题的全部內容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網站內容還不錯，歡迎將生活随笔推薦給好友。