最优化-黄金分割法原理及代码实现

1.黄金分割法介绍

黄金分割法，也叫0.618法，主要用于单峰函数，通过不断地分割迭代，从而找到近似于最小值的函数值。

1.1单峰函数：

如图所示，

设f(x)是定义在D上的实值函数，如果存在 $x_0\epsilon D$ ，使得对D中当任何 x_1<x_2<x_0 时， f(x_1)<f(x_2) ，当 x_0<x_3<x_4 时有 f(x_3)<f(x_4) ,那么就说是定义在D上的单峰函数

1.2算法原理

为了找出最小点，我们通过不断试探的方法，首先我们需要给定区间 [a,b] ，然受在给定区间内取两点

x_1,x_2 ,

$x_1 = a + (1-\tau )(b-a)$

$x_2 = b - (1-\tau)(b-a)$

$\tau = \frac{\sqrt{5}-1}{2}\approx 0.618$

黄金分割的核心原理就是在于对点的选取：

对于x_1和x_2的选取分为以下两种情况：

$f(x_1) < f(x_2) : b =x_2 ,x_2 = x_1,x_1= a+(b-a)*\tau$

$f(x_1) \geq f(x_2): a= x_1,x_1 =x_2,x_2 = b - (b-a) * \tau$

通过两式的不断迭代，a，b在逐渐靠近，我们设置一个精度epi，当

a-b<epi 停止迭代。

2.代码实现

2.1输入函数

我们自己定义区间a，b的值，并且输入所要求的方程形式如（x^2 + 5*x +2）

a = float(input("初始左端点"))
b = float(input("初始右端点"))
def function(a):
     fx = str_fx.replace("x", "a")  # 所有的"x"换为方法所输入的参数a
     return eval(fx)  #eval类型用于对字符串公式进行转换
init_str = input("请输入一个函数，默认变量为x：\n")  # 输入的最初字符串
str_fx = init_str.replace("^", "**")  # 将所有的“^"替换为python的幂形式"**"

2.2画出函数并判断是否为单峰函数

通过定义drawf函数，并且传入参数a，b和间隔interp，来画出图像，通过图像判断是否可以使用黄金分割法进行一维线性搜索。

from pylab import *
def drawf(a, b, interp):
    x = [a+ele*interp for ele in range(0, int((b-a)/interp))]
    y = [function(ele) for ele in x]
    # y = [function(x)]
    plt.figure(1)
    plt.plot(x, y)
    xlim(a, b)
    #title(color="b")
    plt.show()
draw(a,b,0.05)

如图像所示，我们输入函数（x^2 + 5*x +2）在区间a，b上为凸函数

2.3不断迭代求出最小值和最小值点

def golden_section_search(a, b, eps=1e-8):
    phi = (1 + 5 ** 0.5) / 2  # 黄金分割常数phi
    x1 = b - (b - a) / phi
    x2 = a + (b - a) / phi
    while abs(b - a) > eps:  # 终止条件为区间长度小于给定的值
        if function(x1) < function(x2):
            b = x2
            x2 = x1
            x1 = a + (b - a) / phi
            #plt.plot(x1, function(x1), 'y*')
        else:
            a = x1
            x1 = x2
            x2 = b - (b - a) / phi
            #plt.plot(x2, function(x2), 'y*')
    return (a + b) / 2  # 返回最小值所在的位置
#这里我们设置eps为1e-8
x_min = golden_section_search(a,b)
y_min = function(x_min)
print('输入函数最小值点为：{:.5f}\n最小值为{:.5f}'.format(x_min,y_min))

输出结果为