3-5 牛顿法与阻尼牛顿法
牛顿法
简述
牛顿法（Newton's Method）是一种利用目标函数的二阶导数信息（Hessian矩阵）进行优化的方法。与最速下降法仅使用一阶梯度信息不同，牛顿法通过二次泰勒展开近似原函数，能够更准确地预测极值点的位置，因此具有更快的收敛速度（二阶收敛）。
牛顿法的基本思想是：在当前点处用二次函数近似原目标函数，然后直接求出该二次函数的极小值点作为下一次迭代点。
原理
对于多元函数
，在点
处进行二阶泰勒展开：
其中
为 Hessian矩阵（海森矩阵）：
注意：当函数的所有二阶混合偏导数连续时，有
，即 H矩阵对称。www.jlygroup.net对于多项式函数，Hessian矩阵必然对称。
对于二次函数求极值，令泰勒展开的一阶导数为零：
解得牛顿法迭代公式：
Hessian矩阵的性质
Hessian矩阵的正定性决定了极值点的性质：
正定（所有特征值 > 0）→ 局部极小点
负定（所有特征值 < 0）→ 局部极大点
不定（有正有负特征值）→ 鞍点
半定（有零特征值）→ 需要高阶检验
对于二元二次函数，Hessian矩阵为：
二阶矩阵求逆公式：对于矩阵
，当
时：
牛顿法分析
特点：
原函数为二次函数：一步求出极值点。
原函数为非二次函数：www.ycsjb.com先通过二阶泰勒展开变为二次函数，因此求出的为极值点的近似值。
优点：
收敛速度最快（二阶收敛）
对于二次函数，一次迭代即可收敛
缺点：
Hessian矩阵计算困难（需要计算所有二阶偏导数）
Hessian矩阵求逆计算量大
当Hessian矩阵不正定时，可能不收敛或收敛到鞍点
具体步骤
对于给定初始点
和收敛精度
：
① 确定初始点
，收敛精度
，令
；
② 计算梯度
和 Hessian矩阵
；
③ 判断Hessian矩阵是否正定。若不正定，可采用修正策略（见阻尼牛顿法）；
④ 构造搜索方向
；
⑤ 更新迭代点：
，
；
⑥ 重复②~⑤步，直到：
。
实例计算
这里我们举个例子，手把手带大家使用牛顿法，一步步来计算一个二元二次凸函数的最优解。
题目
解答
首先计算
的梯度函数和Hessian矩阵：
梯度：
Hessian矩阵：
注意：对于二次函数，Hessian矩阵为常数矩阵，不随
变化。
代入
，得：
计算梯度的模长：
不满足收敛条件，迭代继续。
第 1 次迭代 (
)：
计算Hessian矩阵的逆：
行列式：
计算牛顿方向：
更新迭代点：
检验收敛条件：
计算
处的梯度：
梯度的模长：
满足收敛条件，迭代结束！
最终得到最优解：
，
结论：当目标函数
为二次函数时，二阶泰勒展开是精确的，而其中的Hessian矩阵
是一个常数矩阵。因此从任意初始点进行迭代，只需一步迭代即可找到目标函数的极小值点。
代码示例
牛顿法 C++ 代码示例：
// 牛顿法 求多元函数最小值点（支持任意维度）
// func: 目标函数
// grad: 梯度函数
// hess: Hessian矩阵函数
// x0: 初始点（向量）
// eps: 收敛精度
vector<double> newton_method(
function<double(const vector<double>&)> func,
function<vector<double>(const vector<double>&)> grad,
function<vector<vector<double>>(const vector<double>&)> hess,
const vector<double>& x0,
double eps,
int max_iter = 1000
) {
vector<double> x = x0;
int n = x0.size();
cout << "\n迭代过程:" << endl;
cout << "步骤\t点 X\t\t\t函数值 f(X)\t\t梯度范数" << endl;
cout << "-----------------------------------------------------------------" << endl;
for (int iter = 0; iter < max_iter; ++iter) {
// 1. 计算当前点的梯度
vector<double> g = grad(x);
// 2. 计算梯度范数 ||∇f(X)||
double grad_norm = norm(g);
// 3. 检查收敛条件
if (grad_norm < eps) {
cout << iter << "\t";
for (size_t i = 0; i < x.size(); ++i) {
printf("%.6lf", x[i]);
cout << (i == x.size() - 1 ? "" : ", ");
}
printf("\t%.6lf\t\t%.6lf (收敛)\n", func(x), grad_norm);
break;
}
// 4. 计算Hessian矩阵
vector<vector<double>> H = hess(x);
// 5. 求解 H * d = -g （使用高斯消元法或直接求逆）
vector<double> d = solve_linear_system(H, g);
for (int i = 0; i < n; ++i) {
d[i] = -d[i]; // d = -H^(-1) * g
}
// 6. 更新迭代点 X^(k+1) = X^(k) + d^(k)
for (size_t i = 0; i < x.size(); ++i) {
x[i] += d[i];
}
// 打印迭代信息
cout << iter << "\t";
for (size_t i = 0; i < x.size(); ++i) {
printf("%.6lf", x[i]);
cout << (i == x.size() - 1 ? "" : ", ");
}
printf("\t%.6lf\t\t%.6lf\n", func(x), grad_norm);
}
return x;
}
阻尼牛顿法
问题与改进
标准牛顿法存在以下问题：
当Hessian矩阵不正定时，搜索方向可能不是下降方向
步长固定为1，可能导致函数值不降反升
初始点远离极值点时，可能不收敛
阻尼牛顿法（Damped Newton's Method）通过引入一维搜索来确定最优步长，解决了上述问题。
阻尼牛顿法迭代公式
其中：
为牛顿方向
通过一维搜索确定，使
最小，
也称为阻尼因子。
具体步骤
对于给定初始点
和收敛精度
：
① 确定初始点
，收敛精度
，令
；
② 计算梯度
和 Hessian矩阵
以及其逆矩阵
；
③ 构造搜索方向
；
④ 一维搜索：求最优步长
，使
最小；
⑤ 更新迭代点：
，
；
⑥ 重复②~⑤步，直到：
。
程序框图阻尼牛顿法的程序框图
代码示例
阻尼牛顿法 C++ 代码示例：
// 阻尼牛顿法 求多元函数最小值点
vector<double> damped_newton_method(
function<double(const vector<double>&)> func,
function<vector<double>(const vector<double>&)> grad,
function<vector<vector<double>>>(const vector<double>&)> hess,
const vector<double>& x0,
double eps,
int max_iter = 1000
) {
vector<double> x = x0;
int n = x0.size();
cout << "\n迭代过程:" << endl;
cout << "步骤\t点 X\t\t\t函数值 f(X)\t\t梯度范数\t步长" << endl;
cout << "--------------------------------------------------------------------" << endl;
for (int iter = 0; iter < max_iter; ++iter) {
// 1. 计算当前点的梯度
vector<double> g = grad(x);
// 2. 计算梯度范数
double grad_norm = norm(g);
// 3. 检查收敛条件
if (grad_norm < eps) {
cout << iter << "\t";
for (size_t i = 0; i < x.size(); ++i) {
printf("%.6lf", x[i]);
cout << (i == x.size() - 1 ? "" : ", ");
}
printf("\t%.6lf\t\t%.6lf\t\t(收敛)\n", func(x), grad_norm);
break;
}
// 4. 计算Hessian矩阵
vector<vector<double>> H = hess(x);
// 5. 求解牛顿方向 d = -H^(-1) * g
vector<double> d = solve_linear_system(H, g);
for (int i = 0; i < n; ++i) {
d[i] = -d[i];
}
// 6. 定义一维搜索函数 phi(alpha) = f(x + alpha * d)
auto phi = [&](double alpha) {
vector<double> new_x(x.size());
for (size_t i = 0; i < x.size(); ++i) {
new_x[i] = x[i] + alpha * d[i];
}
return func(new_x);
};
// 7. 使用进退法+黄金分割法确定最优步长 alpha
pair<double, double> bracket;
try {
bracket = bracket_minimum(phi, 0.0, 0.1);
} catch (...) {
printf("[Warning] bracket minimum failed! bracket = [0, 10]\n");
bracket = make_pair(0.0, 10.0);
}
double alpha_star = golden_section(phi, bracket.first, bracket.second, eps);
// 8. 更新迭代点 X^(k+1) = X^(k) + alpha* * d^(k)
for (size_t i = 0; i < x.size(); ++i) {
x[i] += alpha_star * d[i];
}
// 打印迭代信息
cout << iter << "\t";
for (size_t i = 0; i < x.size(); ++i) {
printf("%.6lf", x[i]);
cout << (i == x.size() - 1 ? "" : ", ");
}
printf("\t%.6lf\t\t%.6lf\t\t%.6lf\n", func(x), grad_norm, alpha_star);
}
return x;
}
完整示例
C++代码实现（完整）：
总结
牛顿法与阻尼牛顿法的比较：
特性 牛顿法 阻尼牛顿法
步长 固定为1 通过一维搜索确定
收敛性 局部收敛，要求初始点靠近极值点 全局收敛性更好
计算量 较小（无需一维搜索） 较大（需一维搜索）
适用场景 二次函数或初始点较好时 一般情况，更稳健
牛顿法的优缺点：
优点：
收敛速度最快（二阶收敛）
对于二次函数，一次迭代即可收敛到精确解
充分利用了函数的二阶信息
缺点：
需要计算Hessian矩阵（所有二阶偏导数）
需要求解线性方程组或矩阵求逆，计算量大
当Hessian矩阵不正定或奇异时，可能不收敛
对初始点要求较高（标准牛顿法）
阻尼牛顿法的改进：
通过一维搜索确定最优步长，保证了函数值下降
具有更好的全局收敛性
即使初始点远离极值点，也能稳定收敛
适用场景：
目标函数的二阶导数容易计算
对收敛速度要求较高
问题规模适中（Hessian矩阵不太大）
一般地，将牛顿法和阻尼牛顿法统称为牛顿型方法。牛顿型方法总体上迭代次数较少、计算速度较快。但是这类方法的主要缺点是每次迭代都要计算函数的二阶导数矩阵，并对该矩阵求逆。当目标函数的维数高时，算量和存储量大的缺点尤为明显。最速下降法的收敛速度比牛顿型方法慢，而牛顿型方法又存在上述缺点。
至此，我们已经介绍了两种基于梯度信息的优化方法：最速下降法原理简单但收敛较慢，牛顿法收敛速度快但需要计算Hessian矩阵。在接下来的文章中，我们将首先介绍坐标轮换法——一种无需计算梯度的直接搜索方法，它通过沿坐标轴方向依次搜索来寻找极值点，实现更为简单。之后，我们将介绍共轭梯度法，它巧妙地结合了最速下降法和牛顿法的优点，既避免了Hessian矩阵的计算，又克服了最速下降法的锯齿现象，具有较快的收敛速度。
