python实现FR共轭梯度法求极值
python代码实现
import numpy as np
from sympy import *
import math
# 定义符号
x1, x2, t = symbols('x1, x2, t')
def func():
# return (x1 - 1) ** 2 + (x2 - 1) ** 2
return 1.5 * x1 ** 2 + 0.5 * x2 ** 2 - x1 * x2 - 2 * x1
# return x1 ** 2 + 2 * x2 ** 2 - 2 * x1 * x2 - 4 * x1
def grad(data):
# 求梯度向量,data=[data1, data2]
f = func()
grad_vec = [diff(f, x1), diff(f, x2)] # 求偏导数,梯度向量
grad = []
for item in grad_vec:
grad.append(item.subs(x1, data[0]).subs(x2, data[1]))
return grad
def grad_len(grad):
# 梯度向量的模长
vec_len = math.sqrt(pow(grad[0], 2) + pow(grad[1], 2))
return vec_len
def zhudian(f):
# 求得min(t)的驻点
t_diff = diff(f)
t_min = solve(t_diff)
return t_min
def main(X0, theta):
f = func()
grad_vec = grad(X0)
grad_length = grad_len(grad_vec) # 梯度向量的模长
k = 0
grad_list = []
alpha_list = []
while grad_length > theta: # 迭代的终止条件
k += 1
if k == 1:
print("迭代次数:", k, end='\t')
p = -np.array(grad_vec)
grad_list.append(grad_vec)
alpha_list.append(p)
print("负梯度方向:", p, end='\t')
print("共轭方向:", p, end='\t')
# 迭代
X = np.array(X0) + t * p
t_func = f.subs(x1, X[0]).subs(x2, X[1])
t_min = zhudian(t_func)
print("精确搜索到的步长:", t_min, end='\t')
X0 = np.array(X0) + t_min * p
print("下一个迭代点:", X0, end='\t')
grad_vec = grad(X0)
grad_length = grad_len(grad_vec)
print('该点处梯度模长', grad_length)
else:
print("迭代次数:", k, end='\t')
grad_list.append(grad_vec)
alpha = grad_len(grad_list[k - 1]) ** 2 / grad_len(grad_list[k - 2]) ** 2
print("alpha", alpha)
p = -np.array(grad_list[k - 1]) + alpha * np.array(alpha_list[k - 2])
alpha_list.append(p)
print("共轭方向:", p, end='\t')
# 迭代
X = np.array(X0) + t * p
t_func = f.subs(x1, X[0]).subs(x2, X[1])
t_min = zhudian(t_func)
print("精确搜索到的步长:", t_min, end='\t')
X0 = np.array(X0) + t_min * p
print("下一个迭代点:", X0, end='\t')
grad_vec = grad(X0)
grad_length = grad_len(grad_vec)
print('该点处梯度模长', grad_length)
# 给定初始迭代点和阈值
main([-2, 4], 0.1)