Python求解无约束最优化问题:牛顿法和最小二乘法

上学期上了一门课叫做数值计算,需要用python写一些迭代解方程的代码。百度后的结果都不尽如人意,因此自己写了一些,这里将比较重要的两个放上来(考试的时候考到了)。

无约束优化问题的牛顿法
import numpy as np
import sympy as sym

#设置函数并求梯度矩阵和Hesse矩阵
vars = sym.symbols('x1 x2')
f = sym.sympify(['100 * (x2 - x1 ** 2) ** 2 + (1 - x1) ** 2'])
GradientMatrix = sym.zeros(len(vars), 1)
HessiansMatrix = sym.zeros(len(vars), len(vars))

for i,fi in enumerate(f):
    for j,r in enumerate(vars):
        GradientMatrix[j, 0] = sym.diff(fi,r)
        for k,s in enumerate(vars):
            HessiansMatrix[j,k]  =  sym.diff(sym.diff(fi,r),s)

HessiansMatrix = sympy.lambdify([x1, x2], HessiansMatrix)
GradientMatrix = sympy.lambdify([x1, x2], GradientMatrix)

x0 = np.array([[-1.2], [1]])
e = 1e-3

#无约束优化问题的牛顿法
def Newton(x, e):
    dx = np.linalg.solve(HessiansMatrix(x[0, 0], x[1, 0]), GradientMatrix(x[0, 0],x[1, 0]))
    print(dx)
    print(x)
    while(np.sum(abs(dx)) > e):
        x = x - dx
        dx = np.linalg.solve(HessiansMatrix(x[0, 0], x[1, 0]), GradientMatrix(x[0, 0],x[1, 0]))
        print(dx)
        print(x)
        
Newton(x0, e)
无约束优化问题的最小二乘法
import numpy as np
import sympy as sym

#设置函数并求雅可比矩阵
vars = sym.symbols('f1 f2 x1 x2')
f1 = 10 * (x2 - x1 ** 2)
f2 = 1 - x1
funcs = sympy.Matrix([f1, f2])
args = sympy.Matrix([x1, x2])
jac = funcs.jacobian(args)
funcs = sympy.lambdify([x1, x2], funcs)
jac = sympy.lambdify([x1, x2], jac)

x0 = np.array([[-1.2], [1]])
e = 1e-3

#无约束优化问题的最小二乘法
def GaussNewton(x, e):
    while(np.sum(abs(jac(x[0, 0], x[1, 0]).T @ funcs(x[0, 0],x[1, 0]))) > e):
        dx = np.linalg.solve(jac(x[0, 0], x[1, 0]).T @ jac(x[0, 0], x[1, 0]), -jac(x[0, 0], x[1, 0]).T @ funcs(x[0, 0],x[1, 0]))
        x = x + dx
        print(dx)
        print(x)
        
GaussNewton(x0, e)

留下评论