方程数值解法--固定点迭代法与牛顿法原理与c++实现

方程数值求根的常用方法有二分法、固定点迭代法、牛顿法等。这里介绍后面两种方法。

首先定义问题：对于函数$y=f(x)$，在求解零点（根）时，有：
$$f(x)=0\text{\hspace{1em}(1)}$$

不动点迭代法：

我们可以将上式转化为：
$$x=g(x)\text{\hspace{1em}(2)}$$
的形式。上式即为不动点迭代的方程。不动点迭代的具体做法如下：

(1)给定初值$x_k=x_0$

(2)令$x_{k+1}=g(x_k)$

(3)如果不满足迭代条件继续执行(2)，否则终止迭代，则最终方程的根$x\ast =x_{k+1}$.

收敛性条件：

要求$g(x)$的导数在(0,1)之间。证明过程可以看下面文章：

[1](数值分析学习笔记（二） - 知乎 (zhihu.com))

2

[3](非线性方程组求解方法合集——不动点迭代法(包含算法流程及代码) - 知乎 (zhihu.com))

例子：$2x^3-x-1=0$

$g(x)$的形式是多样的，例如$x=2x^3-1$,或者$x=\sqrt[3]{\frac{x+1}{2}}$.对于第一种形式，显然不满足收敛性，所以选择第二种形式：

void fixedPointIter(const double init_value,const int max_iter, double& final_value)
{
	double x = init_value;
	int iter_count = 0;
	double x_1 = x;
	while (true)
	{
		
		//x_1 = 2.0 * std::pow(x,3)  - 1.0;
		x_1 = std::pow((x + 1) * 0.5, 1.0 / 3.0);
		
		std::cout << x_1 << "\n";

		if (iter_count>max_iter||std::abs(x_1-x)<1e-3)
		{
			break;
		}
		else
		{
			x = x_1;
			++iter_count;
		}
	}
	final_value = x_1;
}
int main()
{
	double x = 0.0;
	fixedPointIter(0.0, 50, x);
	std::cout<<x << "\n";
}

牛顿迭代法

对于式(1)，理所当然可以写成：$x=x+f(x)$,但是不满足导数小于1的收敛条件，因此，我们将其写为更好地形式，即令$g(x)=x-\frac{f(x)}{f^{\prime }(x)}$,所以有牛顿法：
$$x_{k+1}=x-\frac{f(x_k)}{f^{\prime }(x_k)}\text{\hspace{1em}(3)}$$
更多理论说明[0](02非线性方程求根：牛顿法 - 知乎 (zhihu.com))，当然也可以从几何角度解释[1]((99+ 封私信 / 80 条消息) 如何通俗易懂地讲解牛顿迭代法求开方（数值分析）？ - 知乎 (zhihu.com))

例子：$2x^3-x-1=0$，code:

void NewtonIter(const double init_value, const int max_iter, double& final_value)
{
	double x = init_value;
	int iter_count = 0;
	double x_1 = x;
	
	while (true)
	{
		double f_x = 2 * x * x * x - x - 1;
		double f_x_derivative = 6 * x * x - 1;
		x_1 = x - f_x / (f_x_derivative );
		std::cout << "newton x:" << x_1 << "\n";

		if (iter_count > max_iter || std::abs(x_1 - x) < 1e-3)
		{
			break;
		}
		else
		{
			x = x_1;
			++iter_count;
		}
	}
	final_value = x_1;
}
int main()
{
	double x = 0.0;
	NewtonIter(0.0, 50, x);
	std::cout << "newton final x:" << x << "\n";
}