黄金分割法求解局部最小值——python实现

黄金分割法搜索局部最小值的原理是基于单峰函数的特性。

1、单峰函数

定义：

设f是定义在闭区间[a,b]上的一元函数，x是f在[a,b]上的极小点，并且对于任意的x1,x2属于[a,b],x1<x2;
若当x2<=x 时，f(x1)>f(x2), 当x*<=x1时，f(x2)>f(x1) 则称f是在闭区间[a,b]上的单峰函数。

形象地描述如图：

图中f(x)和g(x)均为单峰函数，即在一个区间上存在一个极值点可以将函数明确地分成两半具有单调性的曲线。可以看到，单峰只是一个形象的称呼，也有单谷函数（区间）的称谓

2、黄金分割法

黄金分割法就是在单峰（谷）区间内搜索出最优解的算法，具体算法步骤如下图：

3、对黄金分割法的理论阐述

黄金分割法，顾名思义就是借用了黄金分割的办法。在具体说明黄金分割法之前，需要先说明试探法。
对于单峰（谷）区间，搜索出最优点的方法，可以从压缩区间考虑，如果搜索的区间足够小，那么在这个区间内的每一个点都可以当做是极小点。
基于这样的思路，开始进行试探。
在[a,b]，选取两个试探点，x1, x2。计算出两个点的函数值，f(x1),f(x2)
比较f(x1)和f(x2):如果f(x1)>f(x2)存在两种情况
情况（1）

情况（2）

当f(x1)>f(x2)时存在上述两种情况，即x1和x2同在极小点的左侧，还是x1在左侧而x2在右侧。而可以肯定的一点就是:x1肯定不会在极小点的右侧。于是这时候就可以将极小点的范围从[a,b]缩小到[x1,b]。
同样的道理，当f(x1)<f(x2)时，可以将极小点的范围从[a,b]缩小到[a,x2].
以上就是试探法的基本原理，在我们选取试探点的时候，我们当然可以随意地选取，可以取两个试探点是[a,b]范围的随机数，这样也可以达到不断地缩小搜索区间的目的。但是作为一种算法，这存在了大量的随机，其算法的效率就是随机性的。因此，我们需要一种稳定的算法，最好就是每次都可以将区间按已知的比率压缩。
这时候我们可以在上一步的搜索区间中线性地设置试探点，即，两个试探点别是左侧的L分位点（L<0.5)，和区间右侧的L分位点（L<0.5)，这样，每次搜索都可以将区间压缩到原来区间的L倍。这样最终就可以得到足够小的搜索区间。
回到黄金分割数0.618.
于是这时候，我们的疑问转化为为什么L需要设置成1-0.618, 即为何每次要将区间压缩0.618倍。因为毕竟每次的L可以小到非常接近0.5。使得在迭代中压缩空间地更快。
这里面主要有节省算法中资源的考虑，在1-L为黄金分割数时，在进行下一次试探点的选取时候，可以将其中一个试探点取自上一时刻的点，节省了资源。推导过程如下：

4、python实现及验证

初学python练手，先用黄金分割法求解了f=exp(-x)+x**2在单谷区间（0，1）的最小值，最后调用了scipy库中的fminbound验证了程序有效

import math
from scipy.optimize import fmin,fminbound
# one dimension golden search
def f ( x ):
    L = math.exp(-x)+x**2
    return L
def golden(f ,*args):
    if not('a' in dir()):
        a=[]
        b=[]
    a.append(args[0])
    b.append(args[1])
    L=1e-16  # tolrence for convergence
    n=80   # max steps for iteration
    #a,b is the region containing opt point
    lambda1=a[0]+0.382*(b[0]-a[0])
    miu1=a[0]+0.618*(b[0]-a[0])
    #lambda1 miu1 is the test point of golden search method
    for k in range(0,n):
        if abs(b[k]-a[k])<=L:
            solve=(a[k]+b[k])/2
            break
        f_lambda1=f(lambda1)
        f_miu1=f(miu1)
        if f_lambda1>f_miu1:
            a.append(lambda1)
            b.append(b[k])
            lambda2=miu1
            miu2=a[k+1]+0.618*(b[k+1]-a[k+1])
        else:
            a.append(a[k])
            b.append(miu1)
            miu2=lambda1
            lambda2=a[k+1]+0.382*(b[k+1]-a[k+1])
        lambda1=lambda2
        miu1=miu2
    print('optimum point  is ：',solve)
    return solve
print(golden(f,0,1))

验证程序结果：
这里调用的scipy库中的fminboud 函数。

# using embed fuction
min_global=fminbound(f,0,1)
print('min_global= ',min_global)

程序结果：

optimum point  is ： 0.35173284738859334
0.35173284738859334
min_global=  0.35173538036861207

求解正确