最优化大作业(二): 常用无约束最优化方法

2023-12-12 14:08

本文主要是介绍最优化大作业(二): 常用无约束最优化方法,希望对大家解决编程问题提供一定的参考价值,需要的开发者们随着小编来一起学习吧!

  • 问题描述

对以下优化问题

                                                   \small minf\left ( x \right )=x_1^2+x_2^2+x_1x_2-10x_1-4x_2+60

选取初始点\small X_0=\left [ 0,0 \right ]^T,\varepsilon =10^{-2},分别用以下方法求解

(1)最速下降法;

(2)Newton法或修正Newton法;

(3)共轭梯度法。

 

  • 基本原理

(1)最速下降法

 

图1  最速下降法流程图

(2)Newton法

图2  Newton法流程图

 

(3)共轭梯度法

 

图3  共轭梯度法流程图

 

  • 实验结果

(1)最速下降法

迭代

次数

1

2

3

4

5

6

7

8

9

梯度

模值

 

5.4210

1.6680

0.9532

0.2933

0.1676

0.0516

0.0295

0.0091

搜索

方向

 

[9.00,

3.00]

[-1.71,

5.14]

[1.58,

0.53]

[-0.30,

0.90]

[0.28,

0.09]

[-0.05,

0.16]

[0.05,

0.02]

[-0.01,

0.03]

当前

迭代点

(1.00,

1.00)

(7.43,

3.14)

(6.77,

5.12)

(7.90,

5.50)

(7.78,

5.85)

(7.98,

5.91)

(7.96,

5.97)

(8.00,

5.98)

(7.99,

6.00)

当前迭代点值

47.00

14.8571

9.2057

8.2120

8.0373

8.0066

8.0012

8.0002

8.0000

表1  最速下降法迭代过程

图4  最速下降法迭代过程图

 

 

(2)Newton法

迭代次数

1

2

梯度模值

 

0.0000

搜索方向

 

[7.00,5.00]

当前迭代点

(1.00,1.00)

(8.00,6.00)

当前迭代点值

47.00

8.0000

表2  Newton法迭代过程

图5  Newton法迭代过程图

 

(3)共轭梯度法

迭代次数

1

2

3

梯度模值

 

5.4210

0.0000

搜索方向

 

[9.00,3.00]

[1.22,6.12]

当前迭代点

(1.00,1.00)

(7.43,3.14)

(8.00,6.00)

当前迭代点值

47.00

14.8571

8.0000

表3  共轭梯度法迭代过程

 

 

图6  共轭梯度法迭代过程图

 

对比结果可得,三种算法均得到同一个极值点(8, 6)。

 

  • 代码展示
import matplotlib.pyplot as plt
from sympy import *
import numpy as np
from mpl_toolkits.mplot3d import Axes3D
t = symbols('t')# 优化目标函数
def fun1():x1, x2 = symbols('x1 x2')y = np.power(x1, 2) + np.power(x2, 2) - x1*x2 -10 * x1 - 4 *x2 +60return ydef fun2(x1, x2):return np.power(x1, 2) + np.power(x2, 2) - x1 * x2 - 10 * x1 - 4 * x2 + 60# 计算当前梯度
def cal_gradient_cur(X_cur):x1, x2 = symbols('x1 x2')f = fun1()g = [diff(f, x1), diff(f, x2)]g[0] = g[0].evalf(subs={x1:X_cur[0], x2:X_cur[1]})g[1] = g[1].evalf(subs={x1:X_cur[0], x2:X_cur[1]})return np.array(g)# 计算lambda, X1: 上一轮迭代点, X2: 本次迭代点
def cal_lambda(X1, X2):g1 = np.array(cal_gradient_cur(X1))g2 = np.array(cal_gradient_cur(X2))g1_norm_2 = np.sum(g1**2, axis=0)g2_norm_2 = np.sum(g2**2, axis=0)lamda = g2_norm_2 / g1_norm_2return lamda# 更新迭代点X
def update_X(X, P):return np.array(X + t*P)# 更新迭代点X
def update_X_cur(X, t, P):return np.array(X + t*P)# 计算最优步长
def cal_best_t(X_cur):x1, x2 = symbols('x1 x2')f = fun1()f_t = f.subs({x1: X_cur[0], x2: X_cur[1]})return solve(diff(f_t, t), t)# 计算梯度模值
def cal_g_norm_cur(X):g_cur = np.array(cal_gradient_cur(X), dtype=np.float32)return np.sqrt(np.sum(g_cur**2, axis=0))def draw(X0):plt.figure()ax = plt.axes(projection='3d')xx = np.arange(-20, 20, 0.1)yy = np.arange(-20, 20, 0.1)x1, x2 = np.meshgrid(xx, yy)Z = fun2(x1, x2)ax.plot_surface(x1, x2, Z, cmap='rainbow', alpha=0.5)X = np.array([X0[0]])Y = np.array([X0[1]])X, Y = np.meshgrid(X, Y)Z = fun2(X, Y)print("初始点:(%0.2f,%0.2f,%0.2f)" % (X, Y, Z))ax.scatter(X, Y, Z, c='k', s=20)ax.set_xlabel('X')ax.set_ylabel('Y')ax.set_zlabel('Z')# ax.legend()# ax.contour(X,Y,Z, zdim='z',offset=-2,cmap='rainbow)   #等高线图,要设置offset,为Z的最小值return axclass C_gradient(object):def __init__(self, X0):self.X0 = X0# 更新搜索方向def cal_P(self, g_cur, P1, lamda):P = -1 * g_cur + lamda*P1return np.array(P)def search(self):X1 = self.X0g_norm_cur = cal_g_norm_cur(X1)  # 计算梯度模值count = 0result = []if(g_norm_cur <= 0.01):print("极值点为({:.2f},{:.2f})".format(X1[0], X1[1]))x1, x2 = symbols('x1 x2')f = fun1()min_value = f.evalf(subs={x1: X1[0], x2: X1[1]})print("极小值为{:.4f}".format(min_value))else:P = -1 * cal_gradient_cur(X1)  # 计算当前负梯度方向while True:X2 = update_X(X1, P)t_cur = cal_best_t(X2)X2 = update_X_cur(X1, t_cur, P)g_cur = cal_gradient_cur(X2)g_norm_cur = cal_g_norm_cur(X2)x1, x2 = symbols('x1 x2')f = fun1()min_value = f.evalf(subs={x1: X2[0], x2: X2[1]})result.append([float(X2[0]), float(X2[1]), float(min_value)])print("当前梯度模值为{:.4f},搜索方向为[{:.2f},{:.2f}]".format(g_norm_cur, P[0], P[1]), end=" ")print("极值点为({:.2f},{:.2f}), 极小值为{:.4f}".format(result[count][0], result[count][1], result[count][2]))if(g_norm_cur <= 0.01):return np.array(result)else:lamda = cal_lambda(X1, X2)P = self.cal_P(g_cur, P, lamda)X1 = X2count += 1def C_gradient_main():print("当前搜索方法为共轭梯度法")X0 = np.array([1, 1])ax = draw(X0)cg = C_gradient(X0)cg.ax = axresult = cg.search()ax.scatter(np.array([result[:, 0]]), np.array([result[:, 1]]), np.array([result[:, 2]]), c='k', s=20)plt.show()class steepest_gradient(object):def __init__(self, X0):self.X0 = X0def search(self):X1 = self.X0result = []count = 0while True:P = -1 * cal_gradient_cur(X1)  # 计算当前负梯度方向X2 = update_X(X1, P)t_cur = cal_best_t(X2)X2 = update_X_cur(X1, t_cur, P)g_norm_cur = cal_g_norm_cur(X2)x1, x2 = symbols('x1 x2')f = fun1()min_value = f.evalf(subs={x1: X2[0], x2: X2[1]})result.append([float(X2[0]), float(X2[1]), float(min_value)])print("当前梯度模值为{:.4f},搜索方向为[{:.2f},{:.2f}]".format(g_norm_cur, P[0], P[1]), end=" ")print("极值点为({:.2f},{:.2f}), 极小值为{:.4f}".format(result[count][0], result[count][1], result[count][2]))if(g_norm_cur <= 0.01):return np.array(result)else:X1 = X2count += 1def steepest_gradient_main():print("当前搜索方法为最速下降法")X0 = np.array([1, 1])ax = draw(X0)a = steepest_gradient(X0)a.ax = axresult = a.search()ax.scatter(np.array([result[:, 0]]), np.array([result[:, 1]]), np.array([result[:, 2]]), c='k', s=20)plt.show()class Newton(object):def __init__(self, X0):self.X0 = X0def cal_hesse(self):return np.array([[2, -1], [-1, 2]])def search(self):X1 = self.X0count = 0result = []while True:g = cal_gradient_cur(X1)g = g.reshape((1, 2)).Th = np.linalg.inv(self.cal_hesse())P = -1 * np.dot(h, g).ravel()X2 = update_X(X1, P)t_cur = cal_best_t(X2)X2 = update_X_cur(X1, t_cur, P)x1, x2 = symbols('x1 x2')f = fun1()min_value = f.evalf(subs={x1: X2[0], x2: X2[1]})g_norm_cur = cal_g_norm_cur(X2)result.append([float(X2[0]), float(X2[1]), float(min_value)])print("当前梯度模值为{:.4f},搜索方向为[{:.2f},{:.2f}]".format(g_norm_cur, P[0], P[1]), end=" ")print("极值点为({:.2f},{:.2f}), 极小值为{:.4f}".format(result[count][0], result[count][1], result[count][2]))if(g_norm_cur <= 0.01):return np.array(result)else:X1 = X2count += 1def newton_main():print("当前搜索方法为newton法")X0 = np.array([1, 1])ax = draw(X0)b = Newton(X0)result = b.search()ax.scatter(np.array([result[:, 0]]), np.array([result[:, 1]]), np.array([result[:, 2]]), c='k', s=20)plt.show()if __name__ == '__main__':steepest_gradient_main()newton_main()C_gradient_main()

 

这篇关于最优化大作业(二): 常用无约束最优化方法的文章就介绍到这儿,希望我们推荐的文章对编程师们有所帮助!



http://www.chinasem.cn/article/484903

相关文章

Python中反转字符串的常见方法小结

《Python中反转字符串的常见方法小结》在Python中,字符串对象没有内置的反转方法,然而,在实际开发中,我们经常会遇到需要反转字符串的场景,比如处理回文字符串、文本加密等,因此,掌握如何在Pyt... 目录python中反转字符串的方法技术背景实现步骤1. 使用切片2. 使用 reversed() 函

Python中将嵌套列表扁平化的多种实现方法

《Python中将嵌套列表扁平化的多种实现方法》在Python编程中,我们常常会遇到需要将嵌套列表(即列表中包含列表)转换为一个一维的扁平列表的需求,本文将给大家介绍了多种实现这一目标的方法,需要的朋... 目录python中将嵌套列表扁平化的方法技术背景实现步骤1. 使用嵌套列表推导式2. 使用itert

Python使用pip工具实现包自动更新的多种方法

《Python使用pip工具实现包自动更新的多种方法》本文深入探讨了使用Python的pip工具实现包自动更新的各种方法和技术,我们将从基础概念开始,逐步介绍手动更新方法、自动化脚本编写、结合CI/C... 目录1. 背景介绍1.1 目的和范围1.2 预期读者1.3 文档结构概述1.4 术语表1.4.1 核

在Linux中改变echo输出颜色的实现方法

《在Linux中改变echo输出颜色的实现方法》在Linux系统的命令行环境下,为了使输出信息更加清晰、突出,便于用户快速识别和区分不同类型的信息,常常需要改变echo命令的输出颜色,所以本文给大家介... 目python录在linux中改变echo输出颜色的方法技术背景实现步骤使用ANSI转义码使用tpu

Conda与Python venv虚拟环境的区别与使用方法详解

《Conda与Pythonvenv虚拟环境的区别与使用方法详解》随着Python社区的成长,虚拟环境的概念和技术也在不断发展,:本文主要介绍Conda与Pythonvenv虚拟环境的区别与使用... 目录前言一、Conda 与 python venv 的核心区别1. Conda 的特点2. Python v

Spring Boot中WebSocket常用使用方法详解

《SpringBoot中WebSocket常用使用方法详解》本文从WebSocket的基础概念出发,详细介绍了SpringBoot集成WebSocket的步骤,并重点讲解了常用的使用方法,包括简单消... 目录一、WebSocket基础概念1.1 什么是WebSocket1.2 WebSocket与HTTP

SQL Server配置管理器无法打开的四种解决方法

《SQLServer配置管理器无法打开的四种解决方法》本文总结了SQLServer配置管理器无法打开的四种解决方法,文中通过图文示例介绍的非常详细,对大家的学习或者工作具有一定的参考学习价值,需要的... 目录方法一:桌面图标进入方法二:运行窗口进入检查版本号对照表php方法三:查找文件路径方法四:检查 S

MyBatis-Plus 中 nested() 与 and() 方法详解(最佳实践场景)

《MyBatis-Plus中nested()与and()方法详解(最佳实践场景)》在MyBatis-Plus的条件构造器中,nested()和and()都是用于构建复杂查询条件的关键方法,但... 目录MyBATis-Plus 中nested()与and()方法详解一、核心区别对比二、方法详解1.and()

golang中reflect包的常用方法

《golang中reflect包的常用方法》Go反射reflect包提供类型和值方法,用于获取类型信息、访问字段、调用方法等,文中通过示例代码介绍的非常详细,对大家的学习或者工作具有一定的参考学习价值... 目录reflect包方法总结类型 (Type) 方法值 (Value) 方法reflect包方法总结

C# 比较两个list 之间元素差异的常用方法

《C#比较两个list之间元素差异的常用方法》:本文主要介绍C#比较两个list之间元素差异,本文通过实例代码给大家介绍的非常详细,对大家的学习或工作具有一定的参考借鉴价值,需要的朋友参考下吧... 目录1. 使用Except方法2. 使用Except的逆操作3. 使用LINQ的Join,GroupJoin