计算数学

数值优化

一维搜索之连续抛物线插值

在黄金分割法中，除了比较函数值大小之外，函数值本身并未在最优求解中用到。这说明，我们并未充分利用现有的信息来加速迭代收敛。这一节介绍连续抛物线插值法求解一维搜索问题，这是一种较少浪费函数值的方法。
下面阐述连续抛物线插值的基本思想，先在最小值附件选3个点$r,s,t$。计算这3个点的目标函数值$f$，且可以给出过这3个点的抛物线 \[P\left( x \right) = f\left( r \right) + {d_1}\left( {x - r} \right) + {d_3}\left( {x - r} \right)\left( {x - s} \right)\] 其中 \[{d_1} = \frac{{f\left( s \right) - f\left( r \right)}}{{s - r}}\] \[{d_2} = \frac{{f\left( t \right) - f\left( s \right)}}{{t - s}}\] \[{d_3} = \frac{{{d_2} - {d_1}}}{{t - r}}\] 令$P\left( x \right)$可以给出抛物线的极小值公式 \[x = \frac{{r + s}}{2} - \frac{{\left( {f\left( s \right) - f\left( r \right)} \right)\left( {t - r} \right)\left( {t - s} \right)}}{{2\left[ {\left( {s - t} \right)\left( {f\left( t \right) - f\left( s \right)} \right) - \left( {f\left( s \right) - f\left( r \right)} \right)\left( {t - s} \right)} \right]}}\] 此点即为极小值新的近似。有些教材给的算法到此为止，这样的话只能给出较粗的极小值。在实际应用时，可以重复以上过程，用新的x代替r，s，t中最远的或次最远的点，反复迭代，即可以给出需求精度下的极小值。
下面给出连续抛物线求一维搜索算法：

[Input] 最大允许迭代次数、误差容限，初始猜想极小点r，s，t 。
[Output] 极值点的自变量与函数值。
[STEP1] 对$i = 1,2,3, \cdots $，作如下循环 \[x = \frac{{r + s}}{2} - \frac{{\left( {f\left( s \right) - f\left( r \right)} \right)\left( {t - r} \right)\left( {t - s} \right)}}{{2\left[ {\left( {s - t} \right)\left( {f\left( t \right) - f\left( s \right)} \right) - \left( {f\left( s \right) - f\left( r \right)} \right)\left( {t - s} \right)} \right]}}\] \[t = s\] \[s = r\] \[r = x\] 如果$i$超过与设定的最大允许迭代次数，停止循环，并给出告警信息。如果已经满足预设的精度需求，则停止循环。
[STEP2] 输出计算结果，并停止该算法。

一般而言，连续抛物线插值法要快于黄金分割方法，因为它更有效的利用了函数值的信息。用连续抛物线插值法求解函数 \[f\left( x \right) = {x^6} - 11{x^3} + 17{x^2} - 7x + 1\] 在区间$\left[ {0,1} \right]$上的极小值。 go语言版本代码如下。


package main
import (
	"fmt"
	."math"
)
func main() {
/*------------------------------------------------------
  Author  : Song Yezhi                       
  verison : 2021-10-7 18:25
  go build -gcflags "-N -l"  
     
------------------------------------------------------*/
	r :=0.0
	s :=0.7
	t :=1.0
	_ , _ = paraInterp(r,s,t,myfx)

}

func myfx(x float64) float64 {
/*------------------------------------------------------
  Author  : Song Yezhi                       
  verison : 2021-10-7 18:26
  go build -gcflags "-N -l"  
------------------------------------------------------*/
	fx := Pow(x,6)-11.0*Pow(x,3)+17.0*x*x - 7.0*x +1.0
	return fx
}

func paraInterp(r,s,t float64,funcX func(float64) float64) (x, fx float64) {
/*------------------------------------------------------
  Author  : Song Yezhi                       
  verison : 2021-10-7 18:28
  go build -gcflags "-N -l"  
     连续抛物型法一维搜索
  -----------------------------------------------------  
  Input  Parameters :
     r,s,t  ------初值
     funcX -----  目标函数
  Output Parameters :
     x---极小点
     fx--极小点函数值
------------------------------------------------------*/
	var imax int = 20
	var tol float64 = 1e-8
	
	var fs,ft,fr float64
	var tmp1,tmp2 float64
	
	fs = funcX(s)
	fr = funcX(r)
	ft = funcX(t)

	for i := 0; i < imax; i++ {		
		tmp1 = (fs -fr) * (t-r) * (t-s)
		tmp2 = 2.0*((s-r)*(ft-fs) -(fs-fr)*(t-s))
		
		x = (r+s)/2.0 - tmp1/tmp2
		
		if Abs(fs-fr) < tol {
		    break
		}
		
	    t = s
		s = r
		r = x
		
		ft = fs
		fs = fr
		fr = funcX(r)
		// *+++single function evaluation +++*   
		fmt.Printf("i= %4d  x= %12.7f  f(x)=%12.7f \n", i, x, ft)
	}
	return x,fx
}

Astrodynamics

计算数学

目录

数值优化

一维搜索之连续抛物线插值