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

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

相关文章

Linux换行符的使用方法详解

《Linux换行符的使用方法详解》本文介绍了Linux中常用的换行符LF及其在文件中的表示,展示了如何使用sed命令替换换行符,并列举了与换行符处理相关的Linux命令,通过代码讲解的非常详细,需要的... 目录简介检测文件中的换行符使用 cat -A 查看换行符使用 od -c 检查字符换行符格式转换将

SpringBoot实现数据库读写分离的3种方法小结

《SpringBoot实现数据库读写分离的3种方法小结》为了提高系统的读写性能和可用性,读写分离是一种经典的数据库架构模式,在SpringBoot应用中,有多种方式可以实现数据库读写分离,本文将介绍三... 目录一、数据库读写分离概述二、方案一:基于AbstractRoutingDataSource实现动态

Java中的String.valueOf()和toString()方法区别小结

《Java中的String.valueOf()和toString()方法区别小结》字符串操作是开发者日常编程任务中不可或缺的一部分,转换为字符串是一种常见需求,其中最常见的就是String.value... 目录String.valueOf()方法方法定义方法实现使用示例使用场景toString()方法方法

Java中List的contains()方法的使用小结

《Java中List的contains()方法的使用小结》List的contains()方法用于检查列表中是否包含指定的元素,借助equals()方法进行判断,下面就来介绍Java中List的c... 目录详细展开1. 方法签名2. 工作原理3. 使用示例4. 注意事项总结结论:List 的 contain

macOS无效Launchpad图标轻松删除的4 种实用方法

《macOS无效Launchpad图标轻松删除的4种实用方法》mac中不在appstore上下载的应用经常在删除后它的图标还残留在launchpad中,并且长按图标也不会出现删除符号,下面解决这个问... 在 MACOS 上,Launchpad(也就是「启动台」)是一个便捷的 App 启动工具。但有时候,应

SpringBoot日志配置SLF4J和Logback的方法实现

《SpringBoot日志配置SLF4J和Logback的方法实现》日志记录是不可或缺的一部分,本文主要介绍了SpringBoot日志配置SLF4J和Logback的方法实现,文中通过示例代码介绍的非... 目录一、前言二、案例一:初识日志三、案例二:使用Lombok输出日志四、案例三:配置Logback一

Python实现无痛修改第三方库源码的方法详解

《Python实现无痛修改第三方库源码的方法详解》很多时候,我们下载的第三方库是不会有需求不满足的情况,但也有极少的情况,第三方库没有兼顾到需求,本文将介绍几个修改源码的操作,大家可以根据需求进行选择... 目录需求不符合模拟示例 1. 修改源文件2. 继承修改3. 猴子补丁4. 追踪局部变量需求不符合很

mysql出现ERROR 2003 (HY000): Can‘t connect to MySQL server on ‘localhost‘ (10061)的解决方法

《mysql出现ERROR2003(HY000):Can‘tconnecttoMySQLserveron‘localhost‘(10061)的解决方法》本文主要介绍了mysql出现... 目录前言:第一步:第二步:第三步:总结:前言:当你想通过命令窗口想打开mysql时候发现提http://www.cpp

Mysql删除几亿条数据表中的部分数据的方法实现

《Mysql删除几亿条数据表中的部分数据的方法实现》在MySQL中删除一个大表中的数据时,需要特别注意操作的性能和对系统的影响,本文主要介绍了Mysql删除几亿条数据表中的部分数据的方法实现,具有一定... 目录1、需求2、方案1. 使用 DELETE 语句分批删除2. 使用 INPLACE ALTER T

MySQL INSERT语句实现当记录不存在时插入的几种方法

《MySQLINSERT语句实现当记录不存在时插入的几种方法》MySQL的INSERT语句是用于向数据库表中插入新记录的关键命令,下面:本文主要介绍MySQLINSERT语句实现当记录不存在时... 目录使用 INSERT IGNORE使用 ON DUPLICATE KEY UPDATE使用 REPLACE