给定一个关于浮点数 x 的函数 f(x) 以及三个初始的不同猜测值来寻找函数的根,我们需要求出该函数的根。这里,f(x) 可以是代数函数或超越函数。
示例:
输入:一个函数 f(x) = x^3 + 2x^2 + 10x - 20 以及三个初始猜测值 - 0, 1 和 2
输出:根的值为 1.3688 或在允许偏差范围内的任何其他值。
输入:一个函数 f(x) = x^5 - 5x + 2 以及三个初始猜测值 - 0, 1 和 2
输出:根的值为 0.4021 或在允许偏差范围内的任何其他值。
Muller 方法是一种寻找形如 f(x)=0 的方程根的算法。它由 David E. Muller 于 1956 年发现。
该方法首先假设根的三个初始值,然后通过这三个点构造一条抛物线,并将该抛物线与 x 轴的交点作为下一次的近似值。这个过程不断重复,直到找到具有所需精度级别的根为止。
为什么要学习 Muller 方法?
Muller 方法是求根方法之一,其他方法还包括二分法、试位法、割线法和牛顿-拉夫森法。但是,它相较于这些方法具有特定的优势——
- 收敛速度,即我们在每一步中向根靠近的程度,在 Muller 方法中大约是 1.84,而对于割线法是 1.62,对于试位法和二分法则是线性的,即 1。因此,Muller 方法比二分法、试位法和割线法更快。
- 虽然它比收敛速度为 2 的牛顿-拉夫森法慢,但它克服了牛顿-拉夫森法最大的缺点之一,即每一步都需要计算导数。
因此,这表明 Muller 方法是一种计算函数根的有效方法。
算法及其工作原理
- 假设函数的任意三个不同的初始根,设为 x0、x1 和 x2。
- 现在,针对这些点 —— x0、x1 和 x2 的函数 f(x) 值,画一个二次多项式,即一条抛物线。
通过这些点的抛物线 p(x) 的方程如下——
p(x) = c + b(x – x2) + a(x – x2)^2 ,其中 a、b 和 c 是常数。
- 画出抛物线后,找出该抛物线与 x 轴的交点,我们设为 x3。
- 寻找抛物线与 x 轴的交点,即 x3:
– 为了求出 x3(即 p(x) 的根),其中 p(x) = c + b(x – x2) + a(x – x2)^2,使得 p(x3) = c + b(x3 – x2) + a(x3 – x2)^2 = 0,我们需要对 p(x) 应用二次公式。由于会有两个根,但我们需要取离 x_2 更近的那一个。为了避免由于相近数相减而产生的舍入误差,请使用以下方程:
x3 – x2 = \frac{-2c}{b\pm \sqrt{b^{2}-4ac}}
现在,因为 p(x) 的根必须更靠近 x_2,所以我们必须从上述方程可能得出的两个值中,取分母较大的那个值。
– 为了求上述方程中的 a、b 和 c,将 x 在 p(x) 中分别设为 x0、x1 和 x2,并设这些值为 p(x0)、p(x1) 和 p(x2),如下所示——
p(x0) = c + b(x0 – x2) + a(x0 – x2)^2 = f(x0)。
p(x1) = c + b(x1 – x2) + a(x1 – x2)^2 = f(x1)。
p(x2) = c + b(x2 – x2) + a(x2 – x2)^2 = c = f(x2)。
– 因此,我们得到了三个方程和三个变量 —— a、b、c。在解出这些变量的值后,我们得到 a、b 和 c 的如下值——
c = p(x_2) = f(x_2)
b = (d_2*(h_1)^2 - d_1*(h_2)^2 ) / ( h_1h_2 * (h_1 - h_2))
a = (d_1*(h_2) - d_2*(h_1)) / ( h_1h_2 * (h_1 - h_2))
– 其中,
d1 = p(x0) – p(x2) = f(x0) – f(x_2)
d2 = p(x1) – p(x2) = f(x1) – f(x_2)
h1 = x0 – x_2
h2 = x1 – x_2
– 现在,将这些值代入 x3 – x2 的表达式中,求出 x_3。
这就是 p(x) = x_3 的根的求法。
- 如果 x3 非常接近 x2(在允许误差范围内),那么 x3 就成为 f(x) 的根;否则,继续重复寻找下一个 x3 的过程,并将之前的 x1、x2 和 x3 作为新的 x0、x1 和 x2。
代码实现
以下是上述算法的 C++ 实现。
C++
// C++ Program to find root of a function, f(x)
#include
using namespace std;
const int MAX_ITERATIONS = 10000;
// Function to calculate f(x)
float f(float x)
{
// Taking f(x) = x ^ 3 + 2x ^ 2 + 10x - 20
return 1*pow(x, 3) + 2*x*x + 10*x - 20;
}
void Muller(float a, float b, float c)
{
int i;
float res;
for (i = 0;;++i)
{
// Calculating various constants required
// to calculate x3
float f1 = f(a);
float f2 = f(b);
float f3 = f(c);
float d1 = f1 - f3;