每日 29

Newton法程序设计
（1）求解无约束优化问题：；
（2）终止准则取；
（3）完成Newton法（牛顿法）的MATLAB编程、调试；
（4）选取几个与实验二中相同的初始点，并给出相关实验结果的对比及分析（从最优解、最优值、收敛速度（迭代次数）等方面进行比较）；
python：
import numpy as np

目标函数：f(x) = 100*(x1² - x2)² + (x1 - 1)²

def objective(x):
x1, x2 = x
return 100 * (x12 - x2)2 + (x1 - 1)**2

梯度函数：返回梯度向量 ∇f(x)

def gradient(x):
x1, x2 = x
dx1 = 200 * (x12 - x2) * 2 * x1 + 2 * (x1 - 1) # ∂f/∂x1
dx2 = -200 * (x12 - x2) # ∂f/∂x2
return np.array([dx1, dx2])

海森矩阵函数：返回Hesse矩阵 H(x)

def hessian(x):
x1, x2 = x
h11 = 200 * (6 * x1**2 - 2 * x2) + 2 # ∂²f/∂x1²
h12 = -400 * x1 # ∂²f/∂x1∂x2
h22 = 200 # ∂²f/∂x2²
return np.array([[h11, h12], [h12, h22]])

牛顿法主函数

def newton_method(initial_x, tol=1e-5, max_iter=1000):
x = initial_x.copy()
iter_count = 0
history = [x.copy()] # 记录迭代点

while iter_count < max_iter:
    grad = gradient(x)
    grad_norm = np.linalg.norm(grad)

    # 终止条件
    if grad_norm <= tol:
        break

    try:
        H = hessian(x)
        # 求解海森矩阵线性方程组 H*p = -grad 得到搜索方向
        p = np.linalg.solve(H, -grad)
        x = x + p  # 牛顿步长通常为1（阻尼牛顿法可加线搜索，此处简化）
    except np.linalg.LinAlgError:
        print("警告：海森矩阵奇异，迭代终止")
        break

    history.append(x.copy())
    iter_count += 1

return x, objective(x), iter_count, history

测试不同初始点

initial_points = [
np.array([0, 0]), # 初始点1：原点
np.array([2, 3]), # 初始点2：较远点
np.array([-1, 1]) # 初始点3：负半轴附近
]

对比实验二的最速下降法结果

print("--------------------- 牛顿法实验结果 ---------------------")
for idx, x0 in enumerate(initial_points):
x_opt, f_opt, n_iter, _ = newton_method(x0)
print(f"初始点 {idx+1}😊
print(f"最优解 x = {x_opt.round(6)}")
print(f"最优函数值 f(x) = {f_opt.round(6)}")
print(f"迭代次数: {n_iter}")
print("-" * 50)

理论最优解验证（已知最优解为(1,1)，函数值为0）

print("\n--------------------- 理论最优解 ---------------------")
print("理论最优解: x = [1, 1], f(x) = 0")