模擬退火算法

1、模擬退火算法匯報(bào)人： 許炯樓2014xxxxxx 李娜2014200909冒亞婷2014200922李園園20142009231 1 31.1 模擬退火算法的來(lái)源及基本原理w算法的提出算法的提出 模擬退火算法最早的思想由模擬退火算法最早的思想由Metropolis等（等（1953）提出，）提出，1983年年Kirkpatrick等將其應(yīng)用于組合優(yōu)化。等將其應(yīng)用于組合優(yōu)化。w算法的目的算法的目的 解決解決NP復(fù)雜性復(fù)雜性問(wèn)題；問(wèn)題； 克服優(yōu)化過(guò)程陷入局部極??；克服優(yōu)化過(guò)程陷入局部極??； 克服初值依賴性?？朔踔狄蕾囆?。 41.1 模擬退火算法的來(lái)源及基本原理n物理退火過(guò)程：物理退火過(guò)程：退火是指

2、將固體加熱到足夠高的溫度，使分子呈隨機(jī)排列狀態(tài)，退火是指將固體加熱到足夠高的溫度，使分子呈隨機(jī)排列狀態(tài)，然后逐步降溫使之冷卻，最后分子以低能狀態(tài)排列，固體達(dá)到某然后逐步降溫使之冷卻，最后分子以低能狀態(tài)排列，固體達(dá)到某種穩(wěn)定狀態(tài)。種穩(wěn)定狀態(tài)。 加溫過(guò)程加溫過(guò)程增強(qiáng)粒子的熱運(yùn)動(dòng)，消除系統(tǒng)原先可能存在的非均勻態(tài)；等溫過(guò)程等溫過(guò)程對(duì)于與環(huán)境換熱而溫度不變的封閉系統(tǒng)，系統(tǒng)狀態(tài)的自發(fā)變化總是 朝自由能減少的方向進(jìn)行，當(dāng)自由能達(dá)到最小時(shí)，系統(tǒng)達(dá)到平衡態(tài)；冷卻過(guò)程冷卻過(guò)程使粒子熱運(yùn)動(dòng)減弱并漸趨有序，系統(tǒng)能量逐漸下降，從而得到低能的 晶體結(jié)構(gòu)。 51.1 模擬退火算法的來(lái)源及基本原理模仿自然界退火現(xiàn)象而得，利用

3、了物理中固體物質(zhì)的退火過(guò)程與一般優(yōu)化問(wèn)題的相似性； 從某一初始溫度開始，伴隨溫度的不斷下降，結(jié)合概率突跳特性在解空間中隨機(jī)尋找全局最優(yōu)解。 61.1 模擬退火算法的來(lái)源及基本原理w數(shù)學(xué)表述數(shù)學(xué)表述 在溫度在溫度T，分子停留在狀態(tài)，分子停留在狀態(tài)r滿足滿足Boltzmann概率分布概率分布DsBBBTksETZTZkrrEETkrETZrEEP)(exp)()(Boltzmann0)()(exp)(1)(子：為概率分布的標(biāo)準(zhǔn)化因常數(shù)。為的能量，表示狀態(tài)機(jī)變量，表示分子能量的一個(gè)隨 71.1 模擬退火算法的來(lái)源及基本思想w數(shù)學(xué)表述數(shù)學(xué)表述 在在同一個(gè)溫度同一個(gè)溫度T，選定兩個(gè)能量，選定兩個(gè)能量E1

4、01 81.1 模擬退火算法的來(lái)源及基本思想BoltzmanBoltzman概率分布概率分布告訴我們：告訴我們： （1）在同一個(gè)溫度，分子停留在能量小狀態(tài)的概率大于停留在能量大狀態(tài)）在同一個(gè)溫度，分子停留在能量小狀態(tài)的概率大于停留在能量大狀態(tài)的概率的概率 （2）溫度越高，不同能量狀態(tài)對(duì)應(yīng)的概率相差越??；溫度足夠高時(shí)，各狀）溫度越高，不同能量狀態(tài)對(duì)應(yīng)的概率相差越小；溫度足夠高時(shí)，各狀態(tài)對(duì)應(yīng)概率基本相同。態(tài)對(duì)應(yīng)概率基本相同。 （3）隨著溫度的下降，能量最低狀態(tài)對(duì)應(yīng)概率越來(lái)越大；溫度趨于）隨著溫度的下降，能量最低狀態(tài)對(duì)應(yīng)概率越來(lái)越大；溫度趨于0時(shí)，時(shí)，其狀態(tài)趨于其狀態(tài)趨于1 91.1 模擬退火算法的

5、來(lái)源及基本原理w數(shù)學(xué)表述數(shù)學(xué)表述 若若|D|為狀態(tài)空間為狀態(tài)空間D中狀態(tài)的個(gè)數(shù)，中狀態(tài)的個(gè)數(shù)，D0是具有最低能量的狀態(tài)集合：是具有最低能量的狀態(tài)集合： 當(dāng)溫度很高時(shí)，每個(gè)狀態(tài)概率基本相同，接近平均值當(dāng)溫度很高時(shí)，每個(gè)狀態(tài)概率基本相同，接近平均值1/|D|； 狀態(tài)空間存在超過(guò)兩個(gè)不同能量時(shí)，具有最低能量狀態(tài)的概率超出平均值狀態(tài)空間存在超過(guò)兩個(gè)不同能量時(shí)，具有最低能量狀態(tài)的概率超出平均值1/|D| ； 當(dāng)溫度趨于當(dāng)溫度趨于0時(shí)，分子停留在最低能量狀態(tài)的概率趨于時(shí)，分子停留在最低能量狀態(tài)的概率趨于1。能量最低狀態(tài)能量最低狀態(tài) 非能量最低狀態(tài)非能量最低狀態(tài) 101.1 模擬退火算法的來(lái)源及基本原理wM

6、etropolis準(zhǔn)則（準(zhǔn)則（1953）以概率接受新狀態(tài)以概率接受新狀態(tài) 固體在恒定溫度下達(dá)到熱平衡的過(guò)程可以用固體在恒定溫度下達(dá)到熱平衡的過(guò)程可以用Monte Carlo方法（計(jì)算機(jī)方法（計(jì)算機(jī)隨機(jī)模擬方法）加以模擬，雖然該方法簡(jiǎn)單，但必須大量采樣才能得到比較精隨機(jī)模擬方法）加以模擬，雖然該方法簡(jiǎn)單，但必須大量采樣才能得到比較精確的結(jié)果，計(jì)算量很大。確的結(jié)果，計(jì)算量很大。 若在溫度若在溫度T，當(dāng)前狀態(tài)，當(dāng)前狀態(tài)i 新狀態(tài)新狀態(tài)j; 若若Ej=randrom0,1 s=sj; Until 抽樣穩(wěn)定準(zhǔn)則滿足；抽樣穩(wěn)定準(zhǔn)則滿足； 退溫退溫tk+1=update(tk)并令并令k=k+1； Unti

7、l 算法終止準(zhǔn)則滿足；算法終止準(zhǔn)則滿足； 輸出算法搜索結(jié)果。輸出算法搜索結(jié)果。u基本步驟 121.1模擬退火算法的基本思想和步驟模擬退火算法的基本思想和步驟 給定初溫給定初溫t=t0，隨機(jī)產(chǎn)生初始狀態(tài)，隨機(jī)產(chǎn)生初始狀態(tài)s=s0，令，令k=0； Repeat Repeat 產(chǎn)生新狀態(tài)產(chǎn)生新狀態(tài)sj=Genete(s)； if min1,exp-(C(sj)-C(s)/tk=randrom0,1 s=sj; Until 抽樣穩(wěn)定準(zhǔn)則滿足；抽樣穩(wěn)定準(zhǔn)則滿足； 退溫退溫tk+1=update(tk)并令并令k=k+1； Until 算法終止準(zhǔn)則滿足；算法終止準(zhǔn)則滿足； 輸出算法搜索結(jié)果。輸出算法搜索結(jié)

8、果。u影響優(yōu)化結(jié)果的主要因素三函數(shù)兩準(zhǔn)則三函數(shù)兩準(zhǔn)則初始溫度初始溫度 131.1模擬退火算法的基本思想和步驟模擬退火算法的基本思想和步驟u模擬退火算法的步驟Step1 設(shè)定初始溫度設(shè)定初始溫度t = tmax, 任選初始解任選初始解r = r0Step2 內(nèi)循環(huán)內(nèi)循環(huán) Step2.1 從從r的鄰域中隨機(jī)選一個(gè)解的鄰域中隨機(jī)選一個(gè)解rt, 計(jì)算計(jì)算r和和rt對(duì)應(yīng)目標(biāo)函對(duì)應(yīng)目標(biāo)函 數(shù)值數(shù)值, 如如rt對(duì)對(duì) 應(yīng)目標(biāo)函數(shù)值較小，則令應(yīng)目標(biāo)函數(shù)值較小，則令r = rt; 否則若否則若 exp(E(rt)E(r)/t)random(0,1), 則令則令r=rt. Step2.2 不滿足內(nèi)循環(huán)停止條件時(shí)，重

9、復(fù)不滿足內(nèi)循環(huán)停止條件時(shí)，重復(fù)Step2.1Step3 外循環(huán)外循環(huán) Step3.1 降溫降溫t = decrease(t) Step3.2 如不滿足外循環(huán)停止條件，則轉(zhuǎn)如不滿足外循環(huán)停止條件，則轉(zhuǎn)Step2；否則算法結(jié)束；否則算法結(jié)束1. 達(dá)到終止溫度2. 達(dá)到迭代次數(shù)3. 最優(yōu)值連續(xù)若干步 保持不變1. 目標(biāo)函數(shù)均值穩(wěn)定2. 連續(xù)若干步的目標(biāo)值變化較小3. 固定的抽樣步數(shù)模擬退火算法的馬氏鏈描述模擬退火算法的馬氏鏈描述2 2 152.1馬爾科夫鏈u定義) 1()(Pr) 1(,) 1 (,) 0()(Pr )( )(,1021inXjnXinXiXiXjnXZnkXkkXss，滿足稱為馬爾

10、可夫鏈，若隨機(jī)序列時(shí)刻狀態(tài)變量的取值。為間，為所有狀態(tài)構(gòu)成的解空令一步轉(zhuǎn)移概率：一步轉(zhuǎn)移概率：n步轉(zhuǎn)移概率：步轉(zhuǎn)移概率：) 1()(Pr) 1(,inXjnXnpji若解空間有限，稱馬爾可夫鏈為若解空間有限，稱馬爾可夫鏈為有限狀態(tài)有限狀態(tài)；)0()(Pr)(,iXjnXpnji若若 ，稱馬爾可夫鏈為，稱馬爾可夫鏈為時(shí)齊的時(shí)齊的。) 1()(,npnpZnjiji 162.2 模擬退火算法與馬爾科夫鏈模擬退火算法與馬爾科夫鏈w模擬退火算法對(duì)應(yīng)了一個(gè)馬爾可夫鏈模擬退火算法對(duì)應(yīng)了一個(gè)馬爾可夫鏈 模擬退火算法：新狀態(tài)接受概率僅依賴于新狀態(tài)和當(dāng)前狀態(tài)，并由溫度加以控制。 若固定每一溫度，算法均計(jì)算馬氏鏈

11、的變化直至平穩(wěn)分布，然后下降溫度，則稱為時(shí)齊算法； 若無(wú)需各溫度下算法均達(dá)到平穩(wěn)分布，但溫度需按一定速率下降，則稱為非時(shí)齊算法。w分析收斂性分析收斂性模擬退火算法關(guān)鍵參數(shù)和操作的設(shè)計(jì)3 3 * *3 模擬退火算法關(guān)鍵參數(shù)和操作的設(shè)計(jì)原則原則 產(chǎn)生的候選解應(yīng)遍布全部解空間方法方法 在當(dāng)前狀態(tài)的鄰域結(jié)構(gòu)內(nèi)以一定概率方式（均勻分布、正態(tài)分布、指數(shù)分布等）產(chǎn)生 狀態(tài)產(chǎn)生函數(shù)狀態(tài)產(chǎn)生函數(shù) 狀態(tài)接受函數(shù)狀態(tài)接受函數(shù)原則原則 (1)在固定溫度下，接受使目標(biāo)函數(shù)下降的候選解的概率要大于使目標(biāo)函數(shù)上升的候選解概率； (2)隨溫度的下降，接受使目標(biāo)函數(shù)上升的解的概率要逐漸減小； (3)當(dāng)溫度趨于零時(shí)，只能接受目標(biāo)

12、函數(shù)下降的解。方法方法 具體形式對(duì)算法影響不大 一般采用一般采用min1,exp(-C/t) * *3 模擬退火算法關(guān)鍵參數(shù)和操作的設(shè)計(jì)收斂性分析收斂性分析 通過(guò)理論分析可以得到初溫的解析式，但解決實(shí)際問(wèn)題時(shí)難以得到精確的參數(shù)； 初溫應(yīng)充分大；實(shí)驗(yàn)表明實(shí)驗(yàn)表明 初溫越大，獲得高質(zhì)量解的機(jī)率越大，但花費(fèi)較多的計(jì)算時(shí)間；方法方法 （1）均勻抽樣一組狀態(tài)，以各狀態(tài)目標(biāo)值得方差為初溫； （2）隨機(jī)產(chǎn)生一組狀態(tài)，確定兩兩狀態(tài)間的最大目標(biāo)值差，根據(jù)差值，利用一定的函數(shù)確定初溫； （3）利用經(jīng)驗(yàn)公式。 初溫初溫 * *3 模擬退火算法關(guān)鍵參數(shù)和操作的設(shè)計(jì)時(shí)齊算法的溫度下降函數(shù)時(shí)齊算法的溫度下降函數(shù) （1）

13、，越接近1溫度下降越慢，且其大小可以不斷變化； （2） ，其中t0為起始溫度，K為算法溫度下降的總次數(shù)。 溫度更新函數(shù)溫度更新函數(shù)10 , 0 ,1kttkk0tKkKtk * *3 模擬退火算法關(guān)鍵參數(shù)和操作的設(shè)計(jì)非時(shí)齊模擬退火算法非時(shí)齊模擬退火算法 每個(gè)溫度下只產(chǎn)生一個(gè)或少量候選解時(shí)齊算法時(shí)齊算法常用的常用的Metropolis抽樣穩(wěn)定準(zhǔn)則抽樣穩(wěn)定準(zhǔn)則 （1）檢驗(yàn)?zāi)繕?biāo)函數(shù)的均值是否穩(wěn)定； （2）連續(xù)若干步的目標(biāo)值變化較??； （3）按一定的步數(shù)抽樣。 內(nèi)循環(huán)終止準(zhǔn)則內(nèi)循環(huán)終止準(zhǔn)則 外循環(huán)終止準(zhǔn)則外循環(huán)終止準(zhǔn)則常用方法常用方法 （1）設(shè)置終止溫度的閾值； （2）設(shè)置外循環(huán)迭代次數(shù)； （3）算法

14、搜索到的最優(yōu)值連續(xù)若干步保持不變； （4）概率分析方法。 實(shí)例計(jì)算4 4成都市三環(huán)線-繞城高速西北區(qū)域公交路線的設(shè)計(jì) * *4.1 問(wèn)題描述發(fā)展成都市經(jīng)濟(jì)文化的發(fā)展導(dǎo)致三環(huán)路以外居民不斷普及，繞城高速區(qū)聚集了越來(lái)越多的企業(yè)、學(xué)校和商業(yè)網(wǎng)點(diǎn)。123經(jīng)濟(jì)的發(fā)展凸顯了公共交通網(wǎng)絡(luò)的不完善，成都市現(xiàn)有的公交線路覆蓋范圍主要涉及主城區(qū)，但郊區(qū)的公交線路和數(shù)量都很少，無(wú)法滿足人們的出行要求。鑒于此，對(duì)未開通公交車的地區(qū)進(jìn)行線路的優(yōu)化是十分有必要的。優(yōu)化目標(biāo)：明確選擇成都三環(huán)線繞城高速的西北區(qū)域范圍（必須包括大豐鎮(zhèn)、安靖鎮(zhèn)、犀浦鎮(zhèn)和高新西區(qū)等重點(diǎn)區(qū)域），利用數(shù)學(xué)模型，在考慮人口密度、交通路況等基礎(chǔ)上，對(duì)該區(qū)

15、域的公交路線進(jìn)行合理的規(guī)劃設(shè)計(jì)，以滿足居民在該區(qū)域內(nèi)或區(qū)域外工作和生活的需要； * *4.2 調(diào)查分析 在成都三環(huán)線繞城高速的西北區(qū)域范圍（包括大豐鎮(zhèn)、安靖鎮(zhèn)、犀浦鎮(zhèn)和高新西區(qū)等重點(diǎn)區(qū)域）內(nèi)： 一共收集到了通過(guò)該區(qū)域的公交車輛3838輛（線路7676條，往返線路可以不一致），在該區(qū)域范圍內(nèi)的公交站點(diǎn)，以及每條公交線路的總長(zhǎng)、在區(qū)域范圍內(nèi)的實(shí)際長(zhǎng)度和空間直線長(zhǎng)度，同時(shí)也刻畫了公交站點(diǎn)在區(qū)域范圍內(nèi)的相對(duì)位置。注：在進(jìn)行線路優(yōu)化時(shí)，總的公交線路條數(shù)是不變的（7676條），不添 加額外的公交線路。 時(shí)間代價(jià)矩陣 * *4.3 模型建立（1）數(shù)據(jù)準(zhǔn)備客流矩陣OD距離矩陣 公交網(wǎng)絡(luò)矩陣21路況矩陣345（

16、2）公交空間網(wǎng)絡(luò)解空間M的建立記原來(lái)區(qū)域內(nèi)的路線組成的交通網(wǎng)絡(luò)為S0，根據(jù)S0的每一條線路在區(qū)域上的起訖點(diǎn)，運(yùn)用 前n條最短路算法，確定S0的每一條路線的n條備用路線，建立區(qū)域內(nèi)公交網(wǎng)絡(luò)解空間M。記原來(lái)區(qū)域內(nèi)的路線組成的交通網(wǎng)絡(luò)為S0，顯然 S0 屬于M，所設(shè)計(jì)的新公交路線便從M中選取。（生成代碼見附錄2）Dijkstra其中： * *4.3 模型建立（3）目標(biāo)函數(shù)的確定對(duì)于每一個(gè)交通網(wǎng)絡(luò) 根據(jù)相應(yīng)的數(shù)據(jù)庫(kù)得出的數(shù)據(jù)所確定時(shí)間成本函數(shù) ，我們的目標(biāo)即在解空間內(nèi)找到最優(yōu)解S使得 最小，即求： 表示從站點(diǎn)i到j(luò)站點(diǎn)的第m條直達(dá)公交線路上的人流量； 表示從站點(diǎn)i到j(luò)站點(diǎn)的第m條需換乘公交線路上的人流

17、量； 表示第m條直達(dá)公交線路上從站點(diǎn)i到站點(diǎn)j所需行車時(shí)間； 表示第m條需換乘公交線路上從站點(diǎn)i到站點(diǎn)j所需行車時(shí)間。( )mmmmmijmijrrtrtrijijijiji N j N rDRi N j N trTRMinT Sdtdt mrijdmtrijdmrijtmtrijt,1,1( ,1)( , )mmmmrrrijk kk kk kri jtDv ,1,1( ,1)( , )r( , )mmlmtrtrrijk kk klmk ktri jtDvtr i j +5K （v ：第m條公交線路上相鄰站點(diǎn)k與k+1之間的平均速度；v ：第m條公交線路上相鄰站點(diǎn)k與k+1之間的距離；,1

18、mrk kv,1mrk kD * *4.3 模型建立v ：第m條需換乘公交線路上相鄰站點(diǎn)k與k+1之間的距離；v ：第m條需換乘公交線路上換乘的次數(shù)。,1mtrk kDK轉(zhuǎn)化為解空間的表達(dá)z( )( ,)( )( ,)( )mmmmmijmijrrtrtrijijijiji N j N rDRi N j N trTRMinSdS odtSdS odtS 客戶評(píng)價(jià)內(nèi)容S：解空間產(chǎn)生的一個(gè)公交路線網(wǎng)絡(luò)M：公交路線網(wǎng)絡(luò)解空間 Pih：第i條路線的原起始點(diǎn) li：解空間中第i條線路的一條備選路線 Pit：第i條路線的原始終點(diǎn)Di ：第i條路的原始距離Cover1：交通網(wǎng)絡(luò)的1級(jí)覆蓋率Cover2：交通

19、網(wǎng)絡(luò)的2級(jí)覆蓋率. .1( )0.92( )0.5SMstCover SCoverS511.51.4ihiitiiiiiPlPlll D . .st * *4 利用算法求解第一步 確定目標(biāo)函數(shù)第二步 利用模擬退火算法進(jìn)行最優(yōu)解計(jì)算根據(jù)路況矩陣，以及odod矩陣和相關(guān)參數(shù)確定目標(biāo)函數(shù)T(S)T(S) 設(shè)置初始溫度T，和截止溫度T （min），以及截止指數(shù)k（初始為0）退火指數(shù)r20時(shí)輸出S，否則返回（2）（4）)() (STSTt0t SSS TkzbetP )( S)( tPRSS 第三步 列出所得最優(yōu)解，進(jìn)行合理性測(cè)試。 * *4.4 算法程序function S ,L3new = tuih

20、uoNew(Sc,OD, condition , L,L3,L4 ) % L3 原始解 L4 L4 解空間S=Sc;S=Sc;S2=Sc;S2=Sc;k=0;T=T=100;00;k_ _max=21T_min=20;T_min=20;r=0.95;r=0.95;while Twhile TT T_minmin & k0if de0 S=S2; S=S2; k= =0;elseelse if ( exp( de/T ) rand ( 1 ) ) if ( exp( de/T ) rand ( 1 ) ) S=S2; S=S2; k= =0; end endendendT=r*T;k=k+1;e

21、ndendL3new=L3;L3new=L3;endend模擬退火算法function A = buildA( S )function A = buildA( S )%UNTITLED2 Summary of this function goes here%UNTITLED2 Summary of this function goes here% Detailed explanation goes here% Detailed explanation goes here%到了第三問(wèn)這里有buildA for count version two.buildA for count version

22、two.A=zeros(244,244);A=zeros(244,244);for i=1:244for i=1:244for j=1:244for j=1:244if length(Si,j)=1if length(Si,j)=1 temp=(Si,j); temp=(Si,j); aaa=mean(temp,1); aaa=mean(temp,1); A(i,j)=aaa(2); A(i,j)=aaa(2);endendendendendendfor i=1:244for i=1:244 for j=1:244 for j=1:244 if A(i,j)=0&i=j if A(i,j)=0

23、&i=j A(i,j)=inf; A(i,j)=inf; end end end endendend代價(jià)矩陣A代碼 * *4.4 算法程序function pathaim,dist=dijkstra2(D,s,aim)function pathaim,dist=dijkstra2(D,s,aim)%Dijkstra最短路算法MatlabMatlab程序用于求從起始點(diǎn)s s到其它各點(diǎn)的最短路%D為賦權(quán)鄰接矩陣%d為s s到其它各點(diǎn)最短路徑的長(zhǎng)度%DD記載了最短路徑生成樹tic;tic;m,n=size(D);m,n=size(D);d=inf.d=inf.* *ones(1,m);ones(1,

24、m);d(1,s)=0;d(1,s)=0;dd=zeros(1,m);dd=zeros(1,m);dd(1,s)=1;dd(1,s)=1;y=s;y=s;for iii=1:mfor iii=1:m pathiii,1(1)=iii; pathiii,1(1)=iii;endendwhile length(find(dd=1)mwhile length(find(dd=1)m for i=1:m for i=1:m if dd(i)=0 if dd(i)=0 d(i)=min(d(i),d(y)+D(y,i); d(i)=min(d(i),d(y)+D(y,i); if d(i)=d(y)+D(y,i) if d(i)=d(y)+D(y,i) pathi=pathy;pathi(length(pathy,1)+1)=i; pathi=pathy;pathi(length(pathy,1)+1)=i; end end end end end end ddd=inf; d