最优化方法Python计算:一元函数搜索算法——黄金分割法

news2025/1/20 4:49:32

黄金分割法是包围策略的经典用例。算法思路为:假定目标函数 f ( x ) f(x) f(x)最优解 x 0 x_0 x0含于长度为 λ \lambda λ的区间 ( a 0 , b 0 ) (a_0,b_0) (a0,b0)内。在区间内插入两个备选点 a 1 ′ , b 1 ′ ∈ ( a 0 , b 0 ) a'_1,b'_1\in(a_0,b_0) a1,b1(a0,b0),使得 a 1 ′ < b 1 ′ a'_1<b'_1 a1<b1
a 1 ′ − a 0 = b 0 − b 1 ′ = ρ ( b 0 − a 0 ) = ρ λ a'_1-a_0=b_0-b'_1=\rho(b_0-a_0)=\rho\lambda a1a0=b0b1=ρ(b0a0)=ρλ
其中, ρ < 1 2 \rho<\frac{1}{2} ρ<21。比较 f ( a 1 ′ ) f(a'_1) f(a1) f ( b 1 ′ ) f(b'_1) f(b1),若 f ( a 1 ′ ) < f ( b 1 ′ ) f(a'_1)<f(b'_1) f(a1)<f(b1),则 x 0 ∈ [ a 0 , b 1 ′ ] x_0\in[a_0,b'_1] x0[a0,b1](见下图(a))。否则,即 f ( b 1 ′ ) ≤ f ( a 1 ′ ) ) f(b'_1)\leq f(a'_1)) f(b1)f(a1)),则 x 0 ∈ [ a 1 ′ , b 0 ] x_0\in[a'_1,b_0] x0[a1,b0](见下图(b)))。对于前者,令 [ a 1 , b 1 ] = [ a 0 , b 1 ′ ] [a_1,b_1]=[a_0,b'_1] [a1,b1]=[a0,b1]。相仿地,对后者令 [ a 1 , b 1 ] = [ a 1 ′ , b 0 ] [a_1,b_1]=[a'_1,b_0] [a1,b1]=[a1,b0]。无论那种情形,变换后的区间 [ a 1 , b 1 ] [a_1,b_1] [a1,b1]其长度缩短为 λ ( 1 − ρ ) \lambda(1-\rho) λ(1ρ),且 x 0 ∈ [ a 1 , b 1 ] x_0\in[a_1,b_1] x0[a1,b1]。继续用上述方法,在 [ a 1 , a 2 ] [a_1,a_2] [a1,a2]中插入备选点 a 2 ′ a'_2 a2 b 2 ′ b'_2 b2,可得长度为 λ ( 1 − ρ ) 2 \lambda(1-\rho)^2 λ(1ρ)2且包含 x 0 x_0 x0的压缩区间 [ a 2 , b 2 ] [a_2,b_2] [a2,b2]。按此方式迭代 k k k次,得到的含有 x 0 x_0 x0的压缩区间 [ a k , b k ] [a_k,b_k] [ak,bk]长度为 λ ( 1 − ρ ) k \lambda(1-\rho)^k λ(1ρ)k。对给定的容错误差 ε > 0 \varepsilon>0 ε>0,若 ( 1 − ρ ) k λ < ε (1-\rho)^k\lambda<\varepsilon (1ρ)kλ<ε,则停止迭代,当前区间 [ a k , b k ] [a_k,b_k] [ak,bk]内任一点均可充当最优解 x 0 x_0 x0的近似值。否则继续进行相同的迭代计算,直至满足精度要求。
请添加图片描述
事实上,我们在第一次迭代确定区间 [ a 1 , b 1 ] [a_1,b_1] [a1,b1]时,有一个备选点 a 1 ′ ∈ [ a 1 , b 1 ] a'_1\in[a_1,b_1] a1[a1,b1] = [ a 0 , b 1 ′ ] =[a_0,b'_1] =[a0,b1]见上图(a)所示情形)或 b 1 ′ ∈ [ a 1 , b 1 ] b'_1\in[a_1,b_1] b1[a1,b1] = [ a 1 ′ , b 0 ] =[a'_1,b_0] =[a1,b0]上图(b)),在该点处的函数值已计算过。我们可以利用这个点,作为第二次迭代时要插入的备选点 b 2 ′ b'_2 b2(或 a 2 ′ a'_2 a2)。这样可以减少一次函数值 f ( b 2 ′ ) f(b'_2) f(b2)(或 f ( a 2 ′ ) f(a'_2) f(a2))的计算。以上图(a)情形为例,在第一次迭代中所取的点 a 1 , b 1 a_1,b_1 a1,b1 a 0 , b 1 ′ a_0,b'_1 a0,b1,第二次迭代时,以 a 1 ′ a'_1 a1作为 b 2 ′ b'_2 b2。若选择 ρ = 3 − 5 2 ≈ 0.382 \rho=\frac{3-\sqrt{5}}{2}\approx0.382 ρ=235 0.382,即可使得 [ b 2 ′ , b 1 ′ ] [b'_2,b'_1] [b2,b1]的长度为 ρ ( 1 − ρ ) \rho(1-\rho) ρ(1ρ),如下图所示。
在这里插入图片描述
Python用于科学计算的工具包scipy中,有一个optimize模块,提供大量用于最优化问题解决方案。其中有一个用于计算指定一元函数 f ( x ) f(x) f(x)局部最优解的函数minimize_scalar。该函数常用的接口为
minimize_ scalar(fun, bracket, method) \text{minimize\_{}scalar(fun, bracket, method)} minimize_scalar(fun, bracket, method)
其中,参数fun表示目标函数 f ( x ) f(x) f(x),bracket表示 f ( x ) f(x) f(x)的单峰区间信息,method表示所要采用的搜索算法。我们只要将表示目标函数、单峰区间(可用myBracket函数计算,详见博文《连续函数的单峰区间计算》)及实现搜索算法的函数传递给minimize_scalar的参数fun、bracket和method即可望算得目标函数的局部最优解。用户自定义搜索算法函数需符合下列的接口规范
custmin(fun, bracket, args=(), ..., **options) \text{custmin(fun, bracket, args=(), ..., **options)} custmin(fun, bracket, args=(), ..., **options)
其中函数名custmin可任取,形式参数表fun, bracket, args=(), …中命名参数(如此处的args及其后续的参数)需排列在任意参数(此处的fun,bracket)之后。特殊的options参数是minimize_{}scalar在调用时向本函数传递所需的自定义实际参数的机制。
下列代码实现黄金分割搜索算法。

from scipy.optimize import OptimizeResult
def myGolden(fun,bracket,gtol=1e-6,**options):
   a0,b0=bracket									#初始化区间[a0,b0]
   rho=0.382										#ρ=1-0.618
   lam=b0-a0										#区间长度λ
   a1,b1=a0+lam*rho,b0-lam*rho						#首次插入点
   f1,f2=fun(a1),fun(b1)							#区间端点函数值
   k=1												#迭代次数
   while lam>gtol:									#重复迭代
      if f1<f2:										#情形(a)
         t,b0,f2=0,b1,f1							#保留左端点,更新右端点
      else:											#情形(b)
         t,a0,f1=1,a1,f2							#保留右端点,更新左端点
      lam=lam*(1-rho)								#更新区间长度
      k+=1											#迭代次数自增1
      if t==0:										#情形(a)
         b1=a1										#插入点b1更新为a1
         a1=a0+lam*rho								#更新插入点a1
         f1=fun(a1)									#计算a1处函数值
      else:											#情形(b)
         a1=b1										#插入点a1更新为b1
         b1=b0-lam*rho								#更新插入点b1
         f2=fun(b1)									#计算b1处函数值
   bestx=(a0+b0)/2									#计算最优解近似值   
   besty=fun(bestx)									#计算最优值近似值
   return OptimizeResult(fun=besty, x=bestx, nit=k)

程序的第2~26行定义实现黄金分割算法的Python函数myGolden。参数fun表示目标函数 f ( x ) f(x) f(x),参数bracket表示单峰区间 [ a 0 , b 0 ] [a_0,b_0] [a0,b0],gtol表示容错误差 ε \varepsilon ε,缺省值为 1 0 − 6 10^{-6} 106。minimize_scalar可利用参数**options将自命名参数gtol的具体值传递给由method接收的myGolden。
第3~8行执行初始化操作:第3行从参数bracket读取函数 f ( x ) f(x) f(x)的单峰区间端点 ( a 0 , b 0 ) (a_0,b_0) (a0,b0)赋予a0,b0。第4行设置缩放系数 ρ = 0.382 \rho=0.382 ρ=0.382赋予rho。第5行将区间长度 λ \lambda λ初始化为 b 0 − a 0 b_0-a_0 b0a0赋予lam。第6行用区间端点处的函数值 f ( a 0 ) , f ( b 0 ) f(a_0),f(b_0) f(a0),f(b0)初始化f1,f2。第8行将迭代次数k初始化为1。
第9~23行的while循环执行迭代操作:第10~13行的if-else分支根据条件
f ( a k ) < f ( b k ) f(a_k)<f(b_k) f(ak)<f(bk)
是否成立,是否成立,确定是前图中所示情形(a)还是(b)。若为前者置标志t为0,保留左端点a0(仅将右端点b0更新为插入点b1),并更新右端点处的函数值f2。否则置t为1,保留右端点b0(将左端点a0更新为插入点a1),更新左端点处的函数值f1。第14行计算新的压缩区间 [ a 0 , b 0 ] [a_0,b_0] [a0,b0]的长度 λ \lambda λ,第15行将迭代次数 k k k自增1。第16~{}23行的{\bf{if-else}}分支根据表示不同情形的t值(0或1),确定新的插入点的计算a1、b1的计算,并计算需更新的插入点处的函数值。
当压缩区间 [ a 0 , b 0 ] [a_0,b_0] [a0,b0]的长度 λ \lambda λ小于容错误差 ε \varepsilon ε时,迭代完成。第24行取 [ a 0 , b 0 ] [a_0,b_0] [a0,b0]的中点作为最优解 x 0 x_0 x0的近似值,赋予bestx。第25行计算目标函数在最优解处的近似值,赋予besty。第26行返回值为
OptimizeResult(fun=besty, x=bestx, nit=k) \text{OptimizeResult(fun=besty, x=bestx, nit=k)} OptimizeResult(fun=besty, x=bestx, nit=k)
是用前面算得的最优解 x 0 x_0 x0的近似值bestx,最优解处的函数值 f ( x 0 ) f(x_0) f(x0)的近似值besty以及迭代次数k创建的OptimizeResult类(第1行导入)对象。
例1 用myGolden方法计算函数 f ( x ) = x 4 − 14 x 3 + 60 x 2 − 70 x f(x)=x^4-14x^3+60x^2-70x f(x)=x414x3+60x270x x = 0 x=0 x=0近旁的局部最优解。
:下列代码计算本例

from scipy.optimize import minimize_scalar								#导入minimize_scalar
f=lambda x:x**4-14*x**3+60*x**2-70*x									#设置目标函数
bracket=myBracket(f,0)													#计算单峰区间
res=minimize_scalar(f,bracket,method=myGolden,options={'bracket': bracket, 'gtol':1.48e-8})	#计算最优解
print(res)

程序的第2行定义目标函数 f ( x ) f(x) f(x),第3行调用myBracket函数(详见博文《连续函数的单峰区间计算》)计算 f ( x ) f(x) f(x) x = 0 x=0 x=0附近的单峰区间 [ a 0 , b 0 ] [a_0,b_0] [a0,b0]赋予bracket。第4行调用minimize_scalar,传递 f ( x ) f(x) f(x)给参数fun,传递 [ a 0 , b 0 ] [a_0,b_0] [a0,b0]给bracket,传递myGolden给参数method,传递字典型数据{‘bracket’: bracket, ‘gtol’:1.48e-8}给参数options,籍此向myGolden传递初始区间 [ a 0 , b 0 ] [a_0,b_0] [a0,b0]和容错误差 ε = 1.48 \varepsilon=1.48 ε=1.48 1 0 − 8 10^{-8} 108。运行程序,输出

fun: -24.369601567349775
nit: 38
x: 0.7808836405154187

意味着myGolden以容错误差 ε = 1.48 × 1 0 − 8 \varepsilon=1.48\times10^{-8} ε=1.48×108,迭代38次,算得最优解近似值为0.7808836405154187,最优解处函数近似值-24.369601567349775。
Scipy.optimization模块为minimize_scalar提供了三个常用的搜索算法:
brent bounded golden \begin{array}{l} \text{brent}\\ \text{bounded}\\ \text{golden} \end{array} brentboundedgolden
供程序员选择使用,其中的golden方法就是实现的黄金分割算法。
例2 用Python提供的golden方法计算例1中函数 f ( x ) = x 4 − 14 x 3 + 60 x 2 − 70 x f(x)=x^4-14x^3+60x^2-70x f(x)=x414x3+60x270x x = 0 x=0 x=0近旁的局部最优解。
:下列代码计算本例

from scipy.optimize import minimize_scalar			#导入minimize_scalar
f=lambda x:x**4-14*x**3+60*x**2-70*x				#设置目标函数
bracket=myBracket(f,0)								#计算x=0近旁的单峰区间
res=minimize_scalar(f,bracket,method='golden')		#计算最优解
print(res)

注意第4行调用minimize_scalar时传递给method的参数为’golden’(要打引号),表示系统提供的golden方法。运行程序,输出

fun: -24.369601567355033
message: '\nOptimization terminated successfully;\nThe returned valuesatisfies the termination criteria\n(using xtol = 1.4901161193847656e-08 )'
nfev: 44
nit: 39
success: True
x: 0.7808840597145699

与例1的输出比较,我们的自己的myGolden无论是计算精度还是计算效率并不输给系统提供的golden。

本文来自互联网用户投稿,该文观点仅代表作者本人,不代表本站立场。本站仅提供信息存储空间服务,不拥有所有权,不承担相关法律责任。如若转载,请注明出处:http://www.coloradmin.cn/o/469328.html

如若内容造成侵权/违法违规/事实不符,请联系多彩编程网进行投诉反馈,一经查实,立即删除!

相关文章

【网络编程】网络基础

目录 一、协议分层 1、为什么要分层 2、OSI七层模型 3、TCP/IP四层协议&#xff08;五层协议&#xff09; 二、网络传输流程 1、同一个网段内的两台主机进行文件传输 2、跨网段的主机的文件传输 三、IP地址和MAC地址 1、IP地址 2、MAC地址 网络通信的本质就是跨主机的…

Python 输出日志 print 函数的应用(python专栏001)

在Python中&#xff0c;print()函数是一个用于输出内容到标准输出设备的函数&#xff0c;通常用于调试程序和显示程序运行结果 直接使用如下&#xff1a; print(5) print("早起的年轻人")print()函数可以接受多个参数&#xff0c;将它们用空格隔开&#xff0c;并输…

Docker CE介绍和CentOS7.5.1804离线安装Docker CE

文章目录 Docker CE介绍离线安装问题离线安装步骤Docker CE 的优点 Docker CE介绍 Docker CE&#xff08;CE&#xff0c;Community Edition是社区版本的意思&#xff09;是一个免费、开源的应用容器化平台&#xff0c;可帮助开发人员将应用程序打包成容器&#xff0c;并在任何…

【小白入门Linux】—目录结构介绍(一)

【小白入门Linux】—目录结构介绍&#xff08;一&#xff09; linux的方向&#xff1a; linux运维 linux嵌入式开发 在linux下做各种程序开发 应用领域 个人桌面领域 服务器领域&#xff08;java c程序部署到linux服务器&#xff09; 嵌入式&#xff08;运行稳定&#x…

GRE-MGRE综合实验

拓扑结构&#xff1a; 要求 1、R5为网络运营商&#xff08;ISP&#xff09;&#xff0c;接口IP地址均为公有地址&#xff1b; 2、R1与R5间使用PPP的PAP认证&#xff0c;R5为主认证方&#xff1b;R2与R5间使用PPP的chap认证&#xff0c;R5为主认证方&#xff1b;R3与R5间使用HD…

操作系统进程线程(一)—进程线程协程区别、多进程多线程、进程调度算法、进程线程通信

进程线程协程区别 定义上 进程&#xff1a;资源分配和拥有的基本单位&#xff0c;是调度的基本单位。运行一个可执行程序会创建一个或者多个进程&#xff1b;进程就是运行起来的程序线程&#xff1a;程序执行基本单位&#xff0c;轻量级进程。每个进程中都有唯一的主线程&…

React之Redux的使用

文章目录 Redux 介绍概述Redux 是什么&#xff1f;为什么要使用 Redux&#xff1f;我什么时候应该使用 Redux&#xff1f;Redux 库和工具React-ReduxRedux ToolkitRedux DevTools 扩展 Redux 术语和概念State 管理不可变性 Immutability术语ActionAction CreatorReducerReducer…

ChatGPT会一直火热下去吗?他会是下一个AR,区块链吗?

目录 前言 VR 热潮 区块链热潮 元宇宙热潮 ChatGPT 热潮 理智看待 ChatGPT 前言 如果在今年年底评选 2023 年的年度科技热词&#xff0c;以 ChatGPT 目前的热度&#xff0c;毫无疑问会是今年排名第一的科技热词。 即使今年才过去四个月&#xff0c;但我很难想象还有什么科…

JUC多并发编程 ThreadLocal

ThreadLocal 提供线程局部变量。这些变量与正常的变量不同, 因为每一个线程在访问 ThreadLocal 实例的时候(通过其 get 或 set 方法) 都有自己的&#xff0c;独立初始化的变量副本。ThreadLocal 实例通常是类中的私有静态字段, 使用它的目的是希望将状态(例如, 用户ID或事务ID)…

【OJ比赛日历】快周末了,不来一场比赛吗? #04.29-05.05 #16场

CompHub 实时聚合多平台的数据类(Kaggle、天池…)和OJ类(Leetcode、牛客…&#xff09;比赛。本账号同时会推送最新的比赛消息&#xff0c;欢迎关注&#xff01; 更多比赛信息见 CompHub主页 或 点击文末阅读原文 以下信息仅供参考&#xff0c;以比赛官网为准 目录 2023-04-29&…

SQL学习日记

目录 一、数据定义&#xff08;create&#xff0c;alter&#xff0c;drop&#xff09; 1.1数据类型 1.2定义基本表&#xff08;create&#xff0c;alter&#xff0c;drop&#xff09; 1.3约束 1.3.1主键约束 1.3.2外码约束 ​编辑 补充CASCADE 关键字 1.3.3Check约束 …

自定义类型:结构体

ok&#xff0c;兄弟们&#xff0c;今天来写关于自定义类型的博客&#xff0c;先来看结构体。 结构体 结构是一些值的集合&#xff0c;这些值称为成员变量。结构的每个成员可以是不同类型的变量。 1.结构体类型的声明 struct tag {member-list; }variable-list; 以上就是结构…

学习 Python 之 Pygame 开发魂斗罗(十六)

学习 Python 之 Pygame 开发魂斗罗&#xff08;十六&#xff09; 完成最终的魂斗罗1. 创建Sound类2. 添加背景音乐3. 添加玩家发射子弹音效4. 增加击中boss要害音效5. 击中敌人音效6. 加入进场动画7. 解决玩家掉出地图死亡问题8. 完善玩家游戏失败函数9. 总结 完成最终的魂斗罗…

道达天际首次亮相军博会,“天网融合”引爆全场

4月23日,第十届中国指挥控制大会暨第八届中国(北京)军事智能技术装备博览会(军博会)落下帷幕。北京道达天际科技股份有限公司(简称道达天际)携DAODAJ2天网情报产品体系首次亮相,全面展示天网融合技术、产品、行业解决方案等成果,最新技术应用备受现场观众瞩目,“天网融合”创新…

性能优化对于Android程序员的重要性,看完你就明白

前言 相信我们都使用过Android手机&#xff0c;然后在使用的过程中经常会遇到手机卡顿&#xff0c;应用闪退&#xff0c;画面不流畅等问题&#xff1b;正因为如此&#xff0c;就导致用户体验非常差&#xff0c;最后选择不再使用Android手机。对此&#xff0c;很多公司对Androi…

【软考备战·希赛网每日一练】2023年4月27日

文章目录 一、今日成绩二、错题总结第一题第二题第三题第四题 三、知识查缺 题目及解析来源&#xff1a;2023年04月27日软件设计师每日一练 一、今日成绩 二、错题总结 第一题 解析&#xff1a; ADSL Modem 上网拨号方式有3种&#xff0c;即 专线方式&#xff08;静态IP&#…

半导体运动台基于dsp+fpga+ad+endac的高速数据采集FPGA设计(二)

4 系统 FPGA 程序的设计 4.1 设计方法及逻辑设计概述 4.1.1 开发环境与设计流程 Quartus II 是 Altera 公司综合开发工具&#xff0c;它集成了 FPGA/CPLD 开发过程中所设计 的所有工具和第三方软件接口&#xff0c;支持多时钟分析&#xff0c; LogicLock 基于块的…

Linux内核阅读自学精简教程目录(必读)

学习Linux内核需要一定的计算机基础知识&#xff0c;包括操作系统&#xff0c;计算机网络等。 以下是学习Linux内核的步骤&#xff1a; 了解Linux内核的基本概念和架构&#xff0c;学习Linux内核源代码的组成和结构。学习C语言和汇编语言&#xff0c;这是深入理解Linux内核的…

ssh设置别名 ,登录

1. ssh设置别名&#xff1b; 可以使用alias命令来给ssh命令起别名&#xff0c;例如&#xff1a; alias mysshssh这样就可以使用myssh命令来代替ssh命令了。如果想要永久生效&#xff0c;可以将上述命令加入到~/.bashrc文件中。 &#xff0c;如果没有 ~/.bashrc 此文件&#…

第十八章 迭代器模式

文章目录 前言一、迭代器模式基本介绍二、迭代器模式应用实例完整代码Department 系ComputerCollegeIterator 计算机学院迭代器InfoColleageIterator 信息工程学院迭代器College 学院接口ComputerCollege 计算机学院InfoCollege 信息工程学院OutPutImpl 操作迭代器Clint 测试 三…