无约束最优化牛顿法二次凸函数

### 牛顿法应用于二次凸函数的原理

对于无约束最优化问题，当目标函数为二次凸函数时，牛顿法能够有效地找到全局最小值。通过二阶泰勒级数展开来逼近原函数，在每一步迭代中构建一个局部二次模型并求其极小值作为下一个迭代点。

具体来说，给定当前点 \( x_k \)，则有：

\[ q_k(x) = f_k + g_k^T (x - x_k) + \frac{1}{2} (x - x_k)^T G_k (x - x_k) \]

其中，

- \( f_k = f(x_k) \)
- \( g_k = \nabla f(x_k) \) 是梯度向量
- \( G_k = \nabla^2 f(x_k) \) 是 Hessian 矩阵[^3]。

为了使上述表达式取得最小值，令一阶导等于零得到更新公式:

\[ x_{k+1} = x_k - [G_k^{-1}]g_k \]

这便是牛顿法的核心思想所在[^1]。

### Python 实现示例

下面是针对特定形式的二维二次凸函数 `z = 0.2 * x**2 + 9 * x + 0.5 * y ** 2 + 3 * y` 使用 SymPy 进行自动微分以及应用牛顿法寻找最优解的例子[^4]:

import numpy as np
from sympy import symbols, diff, Matrix


def newton_method(f_str, var_list, init_point=np.array([0., 0.]), tol=1e-8, max_iter=100):
    """
    :param f_str: 函数字符串表示
    :param var_list: 变量列表
    :param init_point: 初始猜测点，默认为[0, 0]
    :param tol: 收敛精度阈值
    :param max_iter: 最大迭代次数
    """

    # 定义符号变量
    vars_sym = symbols(' '.join(var_list))
    
    # 构建函数及其各阶导数
    func = eval(f_str.replace('^', '**'))
    grad = Matrix([diff(func, v) for v in vars_sym])
    hess = Matrix([[diff(grad[i], j) for i in range(len(vars_sym))] for j in vars_sym])

    point = init_point.copy()
    iter_num = 0
    
    while True:
        current_grad = np.array([float(g.subs(dict(zip(vars_sym, point)))) for g in grad]).reshape(-1, 1)
        current_hess_inv = np.linalg.inv(np.array(hess).astype(float))

        delta_x = -(current_hess_inv @ current_grad)

        if all(abs(dx)[0] < tol for dx in delta_x): break
        
        point += delta_x.flatten()

        iter_num += 1
        if iter_num >= max_iter: raise Exception("达到最大迭代次数")

    return {'optimal_solution': tuple(point), 'function_value_at_optimum': float(func.subs(dict(zip(vars_sym, point))))}


if __name__ == "__main__":
    result = newton_method(
        "0.2*x**2 + 9*x + 0.5*y**2 + 3*y", ["x", "y"]
    )
    print(result)

此代码片段展示了如何定义一个通用版本的新顿算法，并将其应用于具体的二次多项式上以获取最佳解决方案。