模擬退火演算法求函式最大、小值——python實現

阿新 • • 發佈：2019-01-28

模擬退火演算法（Simulate Anneal，SA）是一種通用概率演演算法，用來在一個大的搜尋空間內找尋命題的最優解。模擬退火是由S.Kirkpatrick, C.D.Gelatt和M.P.Vecchi在1983年所發明的。V.Černý在1985年也獨立發明此演演算法。模擬退火演算法是解決TSP問題的有效方法之一。

模擬退火的出發點是基於物理中固體物質的退火過程與一般組合優化問題之間的相似性。模擬退火演算法是一種通用的優化演算法，其物理退火過程由加溫過程、等溫過程、冷卻過程這三部分組成。（by 百度百科）

本文討論用python實現簡單模擬退火法，來求函式在某個區間範圍內的最大/最小值。

1. 模擬退火演算法

基本流程圖：
這裡寫圖片描述

下面是模擬退火演算法的主程式simAnneal.py，實現了算一維，二維及多維函式在給定區間內的最大/最小值。

# -*- coding: utf-8 -*-

'''
=========================================
|                kiterun                |
|               2017/08/11              |
|             [email protected]           |
=========================================
''' 

from random import random
import math
import sys

class SimAnneal(object):
    '''
    Simulated annealing algorithm 
    '''
    def __init__(self, target_text='min'):
        self.target_text = target_text

    def newVar(self, oldList, T):
        '''
        :old : list
        :return : list, new solutions based on old solutions
        :T   : current temperature
        ''' 

        newList = [i + (random()*2-1) for i in oldList]
        return newList

    def juge(self, func, new, old, T):
        '''
        matropolise conditions: to get the maximun or minmun
        :new : new solution data from self.newX
        :old : old solution data
        :T   : current temperature

        '''
        dE = func(new) - func(old) if self.target_text == 'max' else func(old) - func(new)
        if dE >= 0:
            x, ans = new, func(new)
        else:
            if math.exp(dE/T) > random():
                x, ans = new,func(new)
            else:
                x, ans = old, func(old)
        return [x, ans]

class OptSolution(object):
    '''
    find the optimal solution.

    '''
    def __init__(self, temperature0=100, temDelta=0.98,
                 temFinal=1e-8, Markov_chain=2000, 
                 result=0, val_nd=[0]):
        # initial temperature
        self.temperature0 = temperature0
        # step factor for decreasing temperature
        self.temDelta = temDelta
        # the final temperature
        self.temFinal = temFinal
        # the Markov_chain length (inner loops numbers)
        self.Markov_chain = Markov_chain
        # the final result
        self.result = result
        # the initial coordidate values: 1D [0], 2D [0,0] ...
        self.val_nd = val_nd

    def mapRange(self, oneDrange):
        return (oneDrange[1]-oneDrange[0])*random() + oneDrange[0]


    def soulution(self, SA_newV, SA_juge, juge_text, ValueRange, func):
        '''
        calculate the extreme value: max or min value
        :SA_newV : function from class SimAnneal().newVar
        :SA_juge : function from class SimAnneal().juge_*
        :ValueRange : [], range of variables, 1D or 2D or 3D...
        :func : target function obtained from user

        '''
        Ti = self.temperature0
        ndim = len(ValueRange)
        f = max if juge_text=='max' else min
        loops =0

        while Ti > self.temFinal:
            res_temp, val_temp = [], []
            preV = [[self.mapRange(ValueRange[j]) for i in range(self.Markov_chain)] for j in range(ndim)]
            newV = [ SA_newV(preV[j], T=Ti) for j in range(ndim)]

            for i in range(self.Markov_chain):
                boolV = True
                for j in range(ndim):
                    boolV &= (ValueRange[j][0]<= newV[j][i] <= ValueRange[j][1])
                if boolV == True:
                    res_temp.append(SA_juge(new=[newV[k][i] for k in range(ndim)], func=func, old=[preV[k][i] for k in range(ndim)], T=Ti)[-1])
                    val_temp.append(SA_juge(new=[newV[k][i] for k in range(ndim)], func=func, old=[preV[k][i] for k in range(ndim)], T=Ti)[0])
                else:
                    continue
                loops += 1

            # get the index of extreme value
            idex = res_temp.index(f(res_temp))
            result_temp = f(self.result, f(res_temp))
            # update the cooordidate of current extrema value 
            self.val_nd = self.val_nd if result_temp == self.result else val_temp[idex]
            # update the extreme value
            self.result = result_temp
            # update the current temperature
            Ti *= self.temDelta
            print(self.val_nd, self.result)

        print(loops)

        #print('The extreme value = %f' %self.result[-1])

2. 利用這個程式我們來找一下一個二維函式的最大/最小值問題。

二維函式：

f(x,y)=ysin(2πx)+xcos(2πy)
example.py


# -*- coding: utf-8 -*-

'''
======================================
|           kiterun                  |
|          2017/08/11                |
|        [email protected]             |
======================================
'''

from random import random
import math
import sys
from time import time
from simAnneal_dev import SimAnneal
from simAnneal_dev import OptSolution
import matplotlib.pyplot as plt
import matplotlib.animation as animation
from mpl_toolkits.mplot3d import Axes3D
import numpy as np

def func(w):
    x, = w
    fx = x + 10*math.sin(5*x) + 7*math.cos(4*x)
    return fx

def func2(w):
    x, y = w
    fxy = y*np.sin(2*np.pi*x) + x*np.cos(2*np.pi*y)
    return fxy

if __name__ == '__main__':

    targ = SimAnneal(target_text='max')
    init = -sys.maxsize # for maximun case
    #init = sys.maxsize # for minimun case
    xyRange = [[-2, 2], [-2, 2]]
    xRange = [[0, 10]]
    t_start = time()

    calculate = OptSolution(Markov_chain=1000, result=init, val_nd=[0,0])
    output = calculate.soulution(SA_newV=targ.newVar, SA_juge=targ.juge, 
                        juge_text='max',ValueRange=xyRange, func=func2)

    '''
    with open('out.dat', 'w') as f:
        for i in range(len(output)):
            f.write(str(output[i]) + '\n')
    '''
    t_end = time()
    print('Running %.4f seconds' %(t_end-t_start))

下圖紅點是找到的-2 <= x,y <= 2區間內最大值位置，程式執行結果在x,y=[1.7573972092698966, -1.9991309314219978]找到最大值3.7543430598423946：
這裡寫圖片描述

下圖紅點是這該區間的最小值位置，執行結果在x,y = [-1.7612279505916202, -1.9998457808015955]找到最小值：-3.7560984312466141
這裡寫圖片描述

3. 程式優化

上面的程式直接採用python的list做資料處理，眾所周知numpy處理資料的速度快很多，所以對上面的程式進行小修改,主程式更改後為simAnneal_dev.py

# -*- coding: utf-8 -*-

'''
=========================================
|                kiterun                |
|               2017/08/11              |
|             [email protected]           |
=========================================
'''
from random import random
import math
import sys
import numpy as np

class SimAnneal(object):
    '''
    Simulated annealing algorithm 
    '''
    def __init__(self, target_text='min'):
        self.target_text = target_text

    def newVar(self, oldList, T):
        '''
        :old : list
        :return : list, new solutions based on old solutions
        :T   : current temperature
        '''
        newList = [i + (random()*2-1) for i in oldList]
        return newList

    def juge(self, func, new, old, T):
        '''
        matropolise conditions: to get the maximun or minmun
        :new : new solution data from self.newX
        :old : old solution data
        :T   : current temperature

        '''
        dE = func(new) - func(old) if self.target_text == 'max' else func(old) - func(new)
        if dE >= 0:
            x, ans = new, func(new)
        else:
            if math.exp(dE/T) > random():
                x, ans = new, func(new)
            else:
                x, ans = old, func(old)
        return [x, ans]

class OptSolution(object):
    '''
    find the optimal solution.

    '''
    def __init__(self, temperature0=100, temDelta=0.98,
                 temFinal=1e-8, Markov_chain=2000, 
                 result=0, val_nd=[0]):
        # initial temperature
        self.temperature0 = temperature0
        # step factor for decreasing temperature
        self.temDelta = temDelta
        # the final temperature
        self.temFinal = temFinal
        # the Markov_chain length (inner loops numbers)
        self.Markov_chain = Markov_chain
        # the final result
        self.result = result
        # the initial coordidate values: 1D [0], 2D [0,0] ...
        self.val_nd = val_nd

    # create unifrom distributed x,y ..., depend on value range
    def mapRange(self, oneDrange):
        return (oneDrange[1]-oneDrange[0])*random() + oneDrange[0]

    def soulution(self, SA_newV, SA_juge, juge_text, ValueRange, func):
        '''
        calculate the extreme value: max or min value
        :SA_newV : function from class SimAnneal().newVar
        :SA_juge : function from class SimAnneal().juge_*
        :ValueRange : [[],], range of variables, 1D or 2D or 3D...
        :func : target function obtained from user

        '''
        Ti = self.temperature0
        ndim = len(ValueRange)
        f = max if juge_text=='max' else min
        nf = np.amax if juge_text=='max' else np.amin
        loops = 0

        while Ti > self.temFinal:
            res_temp = []
            preV = [[self.mapRange(ValueRange[j]) for i in range(self.Markov_chain)] for j in range(ndim)]
            newV = [ SA_newV(preV[j], T=Ti) for j in range(ndim)]

            for i in range(self.Markov_chain):
                boolV = True
                for j in range(ndim):
                    boolV &= (ValueRange[j][0]<= newV[j][i] <= ValueRange[j][1])
                if boolV == True:
                    res_temp.append(SA_juge(new=[newV[k][i] for k in range(ndim)], 
                                    func=func, old=[preV[k][i] for k in range(ndim)], 
                                    T=Ti))
                else:
                    continue
                loops += 1

            # change list to numpy.array
            sol_temp = np.array(res_temp)
            # find the extreme value
            extreme_temp = nf(sol_temp[:, 1])
            # find the row No. and column No. of the extreme value
            re = np.where(sol_temp == extreme_temp)

            result_temp = f(self.result, extreme_temp)
            # update the cooordidate of current extrema value 
            self.val_nd = self.val_nd if result_temp == self.result else sol_temp[re[0][0], 0]
            # update the extreme value
            self.result = result_temp
            # update the current temperature
            Ti *= self.temDelta

        output = [self.val_nd, self.result]
        print(output)
        print('Total loops = %d' %loops)
        return output

執行的速度明顯提升，快了將近一倍（沒有做統計平均，只跑了單次結果）：
這裡寫圖片描述

模擬退火演算法求函式最大、小值——python實現

模擬退火演算法（Simulate Anneal，SA）是一種通用概率演演算法，用來在一個大的搜尋空間內找尋命題的最優解。模擬退火是由S.Kirkpatrick, C.D.Gelatt和M.P.Vecchi在1983年所發明的。V.Černý在1985年也獨

Matlab遺傳演算法求函式最大值

主函式函式main.m global Bitlength%定義3個全域性變數 global boundsbegin global boundsend boundsbegin=-2; boundsend=2; precision=0.0001;%運算精確度

遺傳演算法上機系列之用遺傳演算法求函式最值問題（附自己寫的程式碼）

本文基於下面的最值問題進行求解： maxf(x1,x2)=21.5+x1sin(4πx1)+x2sin(20πx2)\ max f(x_1,x_2)=21.5+x_1sin(4\pi x_1)+x_2sin(20\pi x_2)maxf(x1,x2)=21.

python 利用模擬退火演算法求解TSP最短路徑

而迪傑斯特拉演算法演算法強調的是從全域性的解當中，每次選擇最短路徑的節點作為最優解，因而它的路徑無疑是全域性最優的。但是，當資料量很大的時候，它的效能消耗也是非常嚇人的。對於本例來說，如果採用迪傑斯特拉演算法，它的時間複雜度是O(N!).

應用遺傳演算法求函式最小值

1、遺傳演算法概論遺傳演算法（GA）可能是最早開發出來的模擬生物遺傳系統的演算法模型。它首先由Fraser提出，後來有Bremermann和Reed等人再次提出。最後，Holland對遺傳演算法做了大量工作並使之推廣，因此被認為是遺傳演算法的奠基人。遺傳

Mathematica求函式最大值

具體如圖：可以看出：Maximize和MaxValue求得是精確解，並且Maximize給出最大值時變數的取值，而FindMaxValue給的數值解（非精確解）。如果想求一個定義域上的最大值，如下：

[01字典樹]求序列完美度(求區間最大異或值)

函數表字典 style targe efi cnblogs main code blank https://nanti.jisuanke.com/t/15531 解題關鍵：01字典樹模板，用字典樹保存每個數的二進制表示，從而動態維護區間上的最大異或值，註意添加和刪除都可

Tingq 模糊查詢共多少條數據最大、小、平均、求和值升、降序

string sys post nat sender type asp idt acl 頁面代碼 <%@ Page Language="C#" AutoEventWireup="true" CodeFile="Default3.aspx.cs" Inherits="

最大熵模型及其python實現

剛開始學習最大熵模型的時候，自以為書中的推導都看明白了。等到自己實現時才發現問題多多。因此，這篇部落格將把重點放在python程式的解讀上，為什麼說是解讀呢，因為這個程式不是我寫的（輕點噴~~），這個程式參考了網上的一篇部落格，地址：http://blog.cs

從暴力求解到動態規劃—— 7 種方法求解連續子陣列的最大和問題（python實現）

問題描述已知一個數組 a[n]，裡面存放著浮點數，可能是正數、負數或0。求它的所有連續子陣列中的最大和。連續子陣列：指的是陣列的一個連續切片，即可以表示為 a[i:j],0≤i≤j<n。連續子陣列的和：比如連續子陣列為 a[i:j] ，則和為

類動態規劃求解較小規模的最大團問題（Python實現）

1.圖：由點、邊（點與點之間連線），組成的集合，如點集V=[0,1,2,3,4]，邊集E=[[1,3,4],[2,3,4],[4],[4],[]]，則（V，E）就是一個圖，其表達的意思如下：該圖中含有5個端點，分別為0,1,2,3,4，這些點存在V中，如端點1對應V

[leetcode]53. Maximum Subarray 最大連續子串python實現【medium】

http://blog.csdn.net/zl87758539/article/details/51676108 題目： Maximum Subarray Find the contiguous subarray within an array (contain

模擬退火演算法實現尋找函式最值

模擬退火的演算法思想：模擬退火演算法從某一較高初溫出發，伴隨溫度引數的不斷下降,結合概率突跳特性在解空間中隨機尋找目標函式的全域性最優解，即在區域性最優解能概率性地跳出並最終趨於全域性最優。模擬退火演算法模板：初始溫度 T=100 冷卻速率 rate=0.99 w

模擬退火演算法理論+Python解決函式極值+C++實現解決TSP問題

簡述演算法設計課這周的作業：趕緊寫了先，不然搞不完了。文章目錄簡述演算法理論部分變數簡單分析從狀態轉移概率到狀態概率推導理解當溫度收斂到接近0的時候，收斂到結果理論

動態規劃硬幣問題，求最大、最小可以組成一個錢數的數

核心思想：在我們從 1 元開始依次找零時，可以嘗試一下當前要找零的面值（這裡指 1 元）是否能夠被分解成另一個已求解的面值的找零需要的硬幣個數再加上這一堆硬幣中的某個面值之和，如果這樣分解之後最終的硬幣數是最少的，

陣列求矩陣最大元素值（打擂臺演算法）

有一個3*4的矩陣，要求程式設計序求出其中值最大的那個元素的值，以及其所在的行號和列號。打擂臺，首先上去一個一個比較厲害的boxer，接下來和剩餘的boxer對打，贏著留下，輸者淘汰。 #include<stdio.h> int main() { int i,j;

java語言求陣列最大值、最小值、總和，列印，翻轉、擷取等操作

//Java陣列章節練習題 public class ArrayUtils{ //1.計算陣列中最大值 public static int arrayMaxElement(int[] data){ int max=data[0];

模擬退火法求橢圓離原點最短距離

有一個橢圓方程 ax^2+by^2+cz^2+dyz+exz+fxy ,其中輸入為，a,b,c,d,e,求結果精確到小數點後8位百度了一下，學習到了模擬退火法，寫這篇部落格進行記錄 #include <cstdio> #include <algorithm> #in

Python實現遺傳演算法（二進位制編碼）求函式最優值

目標函式 maxf(x1,x2)=21.5+x1sin(4πx1)+x2sin(20πx2) max f({x_1},{x_2}) = 21.5 + {x_1}\sin (4\pi {x_1}) + {x_2}\sin (20\pi {x_2})

資料結構面試題總結1——陣列：求最大、次大值

一般大家一開始想到的辦法就是一次迴圈，記錄下最大值和最小值。或者就是用兩次冒泡，找到最大值和次大值。這兩種方法實踐複雜度差不多都是O(2n)，如果陣列很長，效率還是不夠高的。注意：直接排序，再選擇最大的兩個值，這並不是一個好辦法，因為我們只需要前兩個數有序，不需要後N-

模擬退火演算法求函式最大、小值——python實現

1. 模擬退火演算法

2. 利用這個程式我們來找一下一個二維函式的最大/最小值問題。

3. 程式優化

相關推薦