最优化方法Python计算：连续函数的单峰区间计算

我们知道，闭区间上的一元连续函数必在区间上取得最大值和最小值。实践中我们需要能数值地确定含有$f(x)$的唯一最优解$x_0$的区间$[a,b]$。这里介绍寻求连续函数$f(x)$在一点$x^{\ast }$附近单峰区间的包围算法及其Python实现。
算法的思想是从$x^{\ast }$开始沿着$f(x)$下降的方向逐步探索，直至第一次遇到上升：
(1)设定初始步长$s$和缩放系数$\lambda$；
(2)设定$a,c$为$x^{\ast }$及$x^{\ast }+s$，并确保$f(a)>f(c)$。必要时需调整搜索方向；
(3)取$b$为$c+s$。比较$f(c)$与$f(b)$，若$f(c)<f(b)$，意味着$f(a)>f(c)<f(b)$，即可断定$x_0\in (a,b)$，$(a,b)$即为所求；
(4)否则，令$a$为$c$，$c$为$b$，并调整步长$s$为$\lambda \times s$，转(3)。

通常，将步长$s$初始化为$|s|=0.01$，$\lambda$初始化为2。当$s$初始化为正数时，首次探索自左向右。反之，$s$取负数初始值则首次探索自右向左。上图给出了一个按此思想探寻包围函数$f(x)$极小值点的区间$(a,b)$的实例。初始时，步长$s<0$，以起始点$x^{\ast }$为$a$，$c=a+s<a$位于$a$的左边。如图中(a)所示。由于$x^{\ast }$位于$f(x)$的下降区间，故$f(c)>f(a)$。交换$a$和$c$，如图中(b)所示。令$s=-s$为下降方向，则此后$s>0$。取$b=c+s$，如图中(c )所示。由于$f(b)<f(c)$，故将$a$置为$c$，$c$置为$b$，如图中(d)所示。将步长扩大$s=\lambda s$，并置$b=c+s$，由于$f(b)$仍然小于$f(c)$（见图中(e)），故再次将$a$和$c$移至$c$，$b$，扩大步长$s$为$\lambda s$并置$b=c+s$。此时，$f(a)>f(c)<f(b)$，见图中(f)。可见$f(x)$的极小值点被区间$(a,b)$所包围，即$f(x)$在$(a,b)$为一单峰函数。故$(a,b)$即为所求。
下列代码实现算法。

def myBracket(f,xstar,s=1e-2,lamd=2):
   a=xstar					#a初始化为xstar
   ya=f(a)
   c=a+s					#c初始化为a+s
   yc=f(c)
   if yc>ya:				#s为上升方向
      a,c=c,a				#交换a,c
      ya,yc=yc,ya
      s=-s					#调整s为下降方向
   b=c+s					#b置为c+s
   yb=f(b)
   while yb<=yc:			#b同a、c处于同一下降区间
      a,ya=c,yc				#a置为c
      c,yc=b,yb				#c置为b
      s*=lamd				#扩大步长s
      b=c+s					#重置b为c+s
      yb=f(b)
   if a>b:					#若a大b小
      a,b=b,a				#交换a,b
   return a,b

函数myBracket有4个参数：f表示函数$f(x)$。xstar表示初始点$x^{\ast }$。s表示步长$s$，缺省值为0.01。lamd表示放大系数$\lambda$，缺省值为2。函数体中第2～5行分别初始化点$a$和$c$以及对应的函数值$y_a$和$y_c$。第6～9行的if语句矫正a，c保证ya>yc，及步长方向与函数值的下降方向一致。第10～11行初始化点$b$及对应的函数值$y_c$。第10～15行的while循环执行区间$(a,b)$的迭代，直至$y_c<y_b$。第16～17行的if语句确保$a<b$。
例1 设函数$f(x)=\sin x$，用myBracket分别对$x^{\ast }=0$、$x^{\ast }=\frac{4\pi }{3}$，计算单峰区间。
解：下列代码完成本例计算。

import numpy as np					#导入numpy
xstar=0								#设置xstar为0
print(myBracket(np.sin, xstar))		#计算附近的单峰区间
xstar=4*np.pi/3						#设置xstar为4pi/3
print(myBracket(np.sin, xstar))		#计算附近的单峰区间

利用行内注释信息，不难理解本程序。运行程序，输出

(-2.55, -0.6300000000000001)
(4.50879020478639, 5.46879020478639)

第1行输出包含距$x^{\ast }=0$最近的局部最小值点$x_0=-\frac{\pi }{2}$的区间$(-2.55,-0.63)$，如下图(a)所示。第2行输出包含距$x^{\ast }=\frac{4\pi }{3}$最近的局部最小值点$x_0=\frac{3\pi }{2}$的区间$(4.51,5.47)$，如下图(b)所示。

写博不易，敬请支持：
如果阅读本文于您有所获，敬请点赞、评论、收藏，谢谢大家的支持！