黄金分割法(0.618法)python实现

黄金分割法的python实现

在准备学优化方法期末考试的时候拿到习题没有答案,就想写一个程序然后和自己算出来的答案对一下就写了个程序,考完试了留着也没用索性写一篇出来造福学弟学妹,黄金分割法通俗的讲就是给一个区间给一个函数在这个区间里求函数的最小值

1.一维搜索黄金分割法

这里搬一下书上的概念书上的概念严谨但确实不好懂做两道例题就好了其实
考察一维搜索问题:

{ m i n φ ( λ ) s . t . a 1 ≤ λ ≤ b 1 \begin{cases} min&\varphi(\lambda)\\ s.t.&a_1\le\lambda\le b_1 \end{cases} {mins.t.φ(λ)a1λb1
其中 φ ( λ ) \varphi(\lambda) φ(λ)为凸函数

1 。 1^。 1 ε > 0 \varepsilon>0 ε>0为允许的搜索区间最后的长度令 α \alpha α=0.618则黄金分割法计算步骤如下

λ 1 = a 1 + ( 1 − α ) ( b 1 − a 1 ) φ ( λ 1 ) μ 1 = a 1 + α ( b 1 − a 1 ) φ ( μ 1 ) \lambda_1=a_1+(1-\alpha)(b_1-a_1) \quad \varphi(\lambda_1)\\ \mu_1=a_1+\alpha(b_1-a_1) \qquad \varphi(\mu_1) λ1=a1+(1α)(b1a1)φ(λ1)μ1=a1+α(b1a1)φ(μ1)
k = 1 k=1 k=1

2 。 2^。 2 b k − a k < ε b_k-a_k<\varepsilon bkak<ε,则计算结束,最优解 λ ∗ ∈ [ a k , b k ] \lambda^*\in[a_k,b_k] λ[ak,bk],可取 λ ∗ ≈ ( 1 / 2 ) ( a k + b k ) \lambda^*\approx(1/2)(a_k+b_k) λ(1/2)(ak+bk);否则 φ ( λ k ) > φ ( μ k ) \varphi(\lambda_k)>\varphi(\mu_k) φ(λk)>φ(μk), 则转 3 。 3^。 3;若 φ ( λ k ) ≤ φ ( μ k ) \varphi(\lambda_k)\le\varphi(\mu_k) φ(λk)φ(μk),则转 4 。 4^。 4

3 。 3^。 3 a k + 1 = λ k , b k + 1 = b k a_{k+1}=\lambda_k,b_{k+1}=b_k ak+1=λk,bk+1=bk,再令 λ k + 1 = μ k , μ k + 1 = a k + 1 + α ( b k + 1 − a k + 1 ) \lambda_{k+1}=\mu_k, \mu_{k+1}=a_{k+1}+\alpha(b_{k+1}-a_{k+1}) λk+1=μk,μk+1=ak+1+α(bk+1ak+1),计算 φ ( μ k + 1 ) \varphi(\mu_{k+1}) φ(μk+1) 5 。 5^。 5

4 。 4^。 4 a k + 1 = a k , b k + 1 = μ k a_{k+1}=a_k,b_{k+1}=\mu_k ak+1=ak,bk+1=μk,再令 μ k + 1 = λ k , λ k + 1 = a k + 1 + ( 1 − α ) ( b k + 1 − a k + 1 ) \mu_{k+1}=\lambda_k, \lambda_{k+1}=a_{k+1}+(1-\alpha)(b_{k+1}-a_{k+1}) μk+1=λk,λk+1=ak+1+(1α)(bk+1ak+1),计算 φ ( μ k + 1 ) \varphi(\mu_{k+1}) φ(μk+1) 5 。 5^。 5

5 。 5^。 5 k = k + 1 k=k+1 k=k+1,返回 2 。 2^。 2

2.python代码实现

例:求解 min ⁡ φ ( λ ) = λ 2 + 2 λ , s . t . − 3 ≤ λ ≤ 5 \min\varphi(\lambda)=\lambda^2+2\lambda,s.t. -3\le\lambda\le 5 minφ(λ)=λ2+2λ,s.t.3λ5
ε = 0.2 \varepsilon=0.2 ε=0.2取精度为0.001
取三位小数用python的round(x,3)实现

# 只需改如下四个参数即可
f = lambda x: x**2 + 2*x
# 定义变量
e = 0.2       # 终止线
a = -3           # 左边
b = 5          # 右边


# 黄金分割法计算法则
lanb = lambda a, b : a + 0.382* (b -a)
mui = lambda a, b : a + 0.618* (b -a)

Mui = round(mui(a, b),3)
Mui_value = round(f(Mui), 3)

Lamb = round(lanb(a, b) ,3)
Lamb_value= round(f(Lamb), 3)

print(f'a={a} b={b}')
print(f'lanb:{Lamb}')
print(f'mui:{Mui}')
print(f'f(lanb)= {Lamb_value}')
print(f'f(mui)= {Mui_value}')
print() #换行最后看着舒服

while b - a > e: #判断是否到达终止线
    if Lamb_value < Mui_value:
        print('f(lanb) < f(mui)')
        print()
        b = Mui
        Mui = Lamb
        Lamb = round(lanb(a, b), 3)
    else:
        print('f(mui) < f(lanb) ')
        print()
        a = Lamb
        Lamb = Mui
        Mui = round(mui(a, b), 3)

    Mui_value = round(f(Mui), 3)
    Lamb_value = round(f(Lamb), 3)
    print(f'a={a} b={b}')
    print(f'lanb:{Lamb}')
    print(f'mui:{Mui}')
    print(f'f(lanb)= {Lamb_value}')
    print(f'f(mui)= {Mui_value}')


print(f'b-a绝对值为{abs(b -a)} 小于终止线 算法停止')
print(f'X⭐包含于[{a}, {b}]')
print(f'X⭐ = {round((a + b)/ 2, 3)}')

3.运行结果

a=-3 b=5
lanb:0.056
mui:1.944
f(lanb)= 0.115
f(mui)= 7.667
f(lanb) < f(mui)

a=-3 b=1.944
lanb:-1.111
mui:0.056
f(lanb)= -0.988
f(mui)= 0.115
f(lanb) < f(mui)

a=-3 b=0.056
lanb:-1.833
mui:-1.111
f(lanb)= -0.306
f(mui)= -0.988
f(mui) < f(lanb) 

a=-1.833 b=0.056
lanb:-1.111
mui:-0.666
f(lanb)= -0.988
f(mui)= -0.888
f(lanb) < f(mui)

a=-1.833 b=-0.666
lanb:-1.387
mui:-1.111
f(lanb)= -0.85
f(mui)= -0.988
f(mui) < f(lanb) 

a=-1.387 b=-0.666
lanb:-1.111
mui:-0.941
f(lanb)= -0.988
f(mui)= -0.997
f(mui) < f(lanb) 

a=-1.111 b=-0.666
lanb:-0.941
mui:-0.836
f(lanb)= -0.997
f(mui)= -0.973
f(lanb) < f(mui)

a=-1.111 b=-0.836
lanb:-1.006
mui:-0.941
f(lanb)= -1.0
f(mui)= -0.997
f(lanb) < f(mui)

a=-1.111 b=-0.941
lanb:-1.046
mui:-1.006
f(lanb)= -0.998
f(mui)= -1.0
b-a绝对值为0.17 小于终止线 算法停止
X⭐包含于[-1.111, -0.941]
X⭐ = -1.026

结束

用的书上的例子,书上的 λ 2 \lambda_2 λ2算成-1.112算错了应该是-1.111所以后面的结果不一样
在这里插入图片描述

评论 2
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包

打赏作者

浩浩的科研笔记

这我为您答疑发送资源的最大动力

¥1 ¥2 ¥4 ¥6 ¥10 ¥20
扫码支付:¥1
获取中
扫码支付

您的余额不足,请更换扫码支付或充值

打赏作者

实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值