Python布伦特求根

基本思想

如果给定求解区间，而且区间内存在零点，则可根据二分法逐渐逼近，最终至少能找到一个零点。

如果给定一个点，同时可以利用该点的斜率，那么可以来通过切线和横轴的交点来不断确定坐标位置，此即牛顿法。

如果有三个点，就可以做出一个二次函数，通过这个二次函数来确定迭代步长的方法，叫做逆二次插值法。

这三种方法在应付不同的求根问题时各有优劣，布伦特方法的思想，便是找到一个判定条件，来决定迭代求根时到底应该采用哪种方案。

记$b_k$为第$k$次的根估计值，那么$b_k$存在一个对位点，满足$|f(a_k)<|f(b_k)|$，且$f(a_k)f(b_k)<0$。则根据$b_k$的取值不同，有下面四个不等式

1. $|\delta |<|b_k-b_{k-1}|$
2. $|\delta |<|b_k-b_{k-1}-b_{k-2}|$
3. $|s-b_k|<\frac{1}{2}|b_k-b_{k-1}|$
4. $|s-b_k|<\frac{1}{2}|b_{k-1}-b_{k-2}|$

其中，$\delta$为步长容忍度，$s$为当前假设的值。

则当下列条件满足任意一个的时候，下一次迭代采用二分法

• 上次迭代为二分法，且不等式1不成立
• 上次迭代为二分法，且不等式3不成立
• 上次迭代为插值法，且不等式2不成立
• 上次迭代为插值法，且不等式4不成立
• 通过插值法得到的临时值，不在区间$[\frac{3a_l+b_k}{4},b_k]$中间

如果采用插值法，若$b_k,b_{k-1},b_{k-1}$各不相同，则采用二次插值，否则用线性插值。

实现

scipy提供了两种布伦特方法，分别是

brenth(f, a, b, args=(), xtol=2e-12, rtol=8.881784197001252e-16, maxiter=100, full_output=False, disp=True)
brentq(f, a, b, args=(), xtol=2e-12, rtol=8.881784197001252e-16, maxiter=100, full_output=False, disp=True)

二者参数基本一致，其中$a,b$为求根区域，$f$为待求根方程，args为除了未知数外，f的其他参数。xtol, rtol为误差相关量，maxiter为最大迭代次数。

这两个函数的区别是，brenth采用了双曲外推。

下面对布伦特求根和二分求根做一下比较，看看速度上孰优孰劣。

from scipy.optimize import bisect, brenth, brentq
import timeit

def func(x):
    return x**3


timeit.timeit(lambda:bisect(func, -5, 5), number=10000)
timeit.timeit(lambda:brenth(func, -5, 5), number=10000)
timeit.timeit(lambda:brentq(func, -5, 5), number=10000)

结果表明，对于$x^3$这个函数的零点，这三个算法并没有体现出时间上的优势。