牛顿迭代法及其实际应用(附C++代码)

牛顿迭代法介绍

     牛顿迭代法主要用于求解对于不便于计算的非线性方程实根的问题,利用牛顿迭代法求出近似解。

原理介绍

对于非线性方程f(x)x_{0}附近进行Tayor展开: 

f(x) = f({x_0}) + (x - {x_0}){f^{(1)}}({x_0}) + {(x - {x_0})^2}\frac{​{​{f^{(2)}}({x_0})}}{​{2!}} + ....

取其线性部分作为非线性方程f(x) = 0的近似方程 

f(x) = f({x_0}) + (x - {x_0}){f^{(1)}}({x_0})

{f^{(1)}}({x_0}) \ne 0,则解为:

{x_1} = {x_0} - \frac{​{f({x_0})}}{​{​{f^{(1)}}({x_0})}}

再把f(x)x_{1}附近展开Tayor级数,取其线性部分作为f(x) = 0的近似方程,若{f^{(1)}}({x_1}) \ne 0,则有:

{x_2} = {x_1} - \frac{​{f({x_1})}}{​{​{f^{(1)}}({x_1})}}

如此可以推导出Newton法:

 {x_{n + 1}} = {x_n} - \frac{​{f({x_n})}}{​{​{f^{(1)}}({x_n})}}

实际应用举例

为了更好理解牛顿迭代法,利用求解y = {x^3} + {e^x}方程根进行举例

求导方程为:y = 3{x^2} + {e^x}

C++代码实现

#include<iostream>
#include<cmath>
using namespace std;
//实例:y=x^3+e^x
double x;        //定义变量x
double X;        //临时变量,记录第k次的x值,用于判断精度,第k次与第k+1次做差比较精度             
double y(double x)          //定义关于x的表达式
{
	return -x*x*x + exp(x);
}
double dy(double x)         //定义导数公式
{
	return -3 * x*x + exp(x);
}

bool accuracy(double x0)
{
	if (fabs(X - x0) > 1e-5)          //定义精度10的负5次方
	{
		return 1;
	}
	else
		return 0;
}
void ND(double x0)    //牛顿迭代法
{  
	int n = 0;            //记录迭代次数
	do{
		double _y = y(x0);
		double _dy = dy(x0);
		X = x0;
		x0 = x0 - _y / _dy;
		n++;
		cout << "第" << n << "次迭代结果为:" ;
		printf("%.5f\n", x0);
	} while (accuracy(x0));
	cout << "共迭代" << n << "次,最终近似解为:" ;
	printf("%.5f\n", x0);
}

void main()
{
	x = 4.5;    //定义x0初始值
	ND(x);
}

故由此可知,函数 y = {x^3} + {e^x}函数近似解为4.53640

补充优化

 在实际求解非线性放出实根时,部分方程可能会出现求导困难的情况,故可以利用差商代替求导,则公式为:

{x_{n + 1}} = {x_n} - \frac{​{f({x_n})}}{​{f({x_n}) - f({x_{n - 1}})}}({x_n} - {x_{n - 1}})

C++代码为:

#include<iostream>
#include<cmath>
using namespace std;

//实例:y=x^3+e^x
double x;        //定义变量x
double X;        //临时变量,记录第k次的x值,用于判断精度,第k次与第k+1次做差比较精度             
double y(double x)          //定义关于x的表达式
{
	return -x*x*x + exp(x);
}

bool accuracy(double x0)
{
	if (fabs(X - x0) > 1e-5)          //定义精度10的负5次方
	{
		return 1;
	}
	else
		return 0;
}
void ND_optimize(double x1,double x0)    //牛顿迭代法差商改进
{
	double _y1 = y(x0);
	double _y2 = y(x1);
	int n = 0;            //记录迭代次数
	do{
		double _y1 = y(x0);
		double _y2 = y(x1);

		X = x0;
		x0 = x0 - (_y1 / (_y1-_y2)) * (x0-x1);
		n++;
		cout << "第" << n << "次迭代结果为:";
		printf("%.5f\n", x0);
	} while (accuracy(x0));
	cout << "共迭代" << n << "次,最终近似解为:";
	printf("%.5f\n", x0);
}

void main()
{
	double x1 = 4.5;    //定义x1初始值
	double x2 = 4;      //定义x2初始值
	ND_optimize(x1,x2);
}

### 高斯牛顿迭代 C++ 实现 以下是基于高斯牛顿迭代C++ 示例代码实现,用于解决非线性最小二乘问题。此方通过不断更新参数估计值来逼近最优解。 #### 示例代码 ```cpp #include <iostream> #include <cmath> #include <vector> using namespace std; // 定义残差函数及其雅可比矩阵计算 double residual(double x, double a, double b) { return exp(a * x + b) - (a * x + b); // 假设的目标函数形式 } void jacobian(vector<double>& J, double x, double a, double b) { J[0] = x * exp(a * x + b) - x; J[1] = exp(a * x + b) - 1; } int main() { const int max_iter = 100; // 最大迭代次数 const double tol = 1e-6; // 收敛精度 vector<double> data_x = {0.0, 0.5, 1.0}; // 输入数据点 vector<double> params(2, 0.0); // 初始参数猜测 [a, b] for (int iter = 0; iter < max_iter; ++iter) { vector<vector<double>> J(data_x.size(), vector<double>(params.size())); vector<double> r(data_x.size()), delta_params(params.size()); // 构建雅可比矩阵和残差向量 for (size_t i = 0; i < data_x.size(); ++i) { r[i] = residual(data_x[i], params[0], params[1]); jacobian(J[i], data_x[i], params[0], params[1]); } // 解正规方程 JTJΔp = -JT*r vector<vector<double>> JTJ(params.size(), vector<double>(params.size(), 0)); vector<double> JTr(params.size(), 0); for (size_t p = 0; p < params.size(); ++p) { for (size_t q = 0; q <= p; ++q) { for (size_t i = 0; i < data_x.size(); ++i) { JTJ[p][q] += J[i][p] * J[i][q]; } JTJ[q][p] = JTJ[p][q]; // 对称性质 } for (size_t i = 0; i < data_x.size(); ++i) { JTr[p] -= J[i][p] * r[i]; } } // 使用高斯消元或其他数值方求解 Δp for (size_t k = 0; k < params.size(); ++k) { for (size_t i = k + 1; i < params.size(); ++i) { double factor = JTJ[i][k] / JTJ[k][k]; for (size_t j = k; j < params.size(); ++j) { JTJ[i][j] -= factor * JTJ[k][j]; } JTr[i] -= factor * JTr[k]; } } for (int i = static_cast<int>(params.size()) - 1; i >= 0; --i) { delta_params[i] = JTr[i]; for (size_t j = i + 1; j < params.size(); ++j) { delta_params[i] -= JTJ[i][j] * delta_params[j]; } delta_params[i] /= JTJ[i][i]; } bool converged = true; for (size_t i = 0; i < params.size(); ++i) { params[i] += delta_params[i]; if (fabs(delta_params[i]) > tol) { converged = false; } } cout << "Iteration " << iter + 1 << ": "; for (auto param : params) { cout << param << " "; } cout << endl; if (converged) break; } cout << "Final parameters: "; for (auto param : params) { cout << param << " "; } cout << endl; return 0; } ``` 上述代码实现了高斯牛顿算法的核心逻辑[^3],包括构建雅可比矩阵、计算残差以及通过正规方程求解参数增量的过程。 --- ### 相关理论说明 高斯牛顿是一种优化技术,适用于非线性最小二乘问题。其基本思想是利用泰勒展开近似目标函数,并通过迭代逐步改进参数估计值。每次迭代过程中,需要计算当前点处的雅可比矩阵并求解正规方程 \( \mathbf{J}^\top\mathbf{J}\Delta\mathbf{p} = -\mathbf{J}^\top\mathbf{r} \)[^1]。 ---
评论 2
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值