前言

本系列涉及线性与非线性规划中的几种规划算法

1.本节介绍复合形方法（单纯形搜索方法）
2.参考：陈宝林-最优化理论与算法
3.采用python编程实现，已测试，代码可行

一、复合形方法

这里的复合形方法，我觉得和参考书中的单纯形搜索方法是一个意思。基本思想是先生成k个可行点，得到一个初始的复合形，然后改变复合形的形状（反射、扩张、收缩、压缩）。
这考虑反射这一种手段的复合形方法的具体细节参见参考书和下述。

二、算法实现

1.算法流程

只是反射方式下的复合形方法的流程如下图：

2.算法细节

3.代码

#复合型方法，随机初始区间X0_R=[[a,b],[],[]...]未给出，Xi自变量，X0给定一个初始可行点初值
#定点数n+1<=k<=2n,收敛精度delta,a0每次缩减的倍数，lamda,a0的最小限制a0_limit=e-10,反射a0=1.3，b0，扩张a1=1,b1,缩小a2=0.7,b2,压缩a3=0.5
import math
from sympy import *
import random
import numpy as np
class Complex_method:def __init__(self, F, G, Xi, X0_R, k, X0=None):self.f = Fself.G = Gself.X0_R = X0_Rself.i = len(X0_R)# 自变量个数self.Xi = Xi  # 定义自变量self.k = k  # 产生的顶点数if not X0:self.Xk = np.array(self.random_X(self.k,[]))#self.Xk = np.array([[7,8],[6,8],[8,7]],dtype=float)else:Xk = []Xk.append(X0)self.Xk = np.array(self.random_X(k-1,Xk))# 生成自变量字典，如何实现字典的numpy向量计算def x_dict(self,X):x = {}for index,item in enumerate(self.Xi):x[item] = X[index]self.X_dict = x#随机顶点Xk是存放顶点的集合def random_X(self,k,Xk):n = 0while n < k:X = []for index,item in enumerate(self.Xi):X.append(random.uniform(self.X0_R[index][0],self.X0_R[index][1]))if self.condition_G(X):Xk.append(X)n += 1return Xk#约束条件def condition_G(self,X):self.x_dict(X)n = Truefor i in self.G:if i.subs(self.X_dict) > 0:n = falsereturn n#计算目标函数值def f_value(self,X):self.x_dict(X)self.fX = self.f.subs(self.X_dict)return self.fX#求各顶点的函数值并排序def fXk_value(self):fXk = []for i in self.Xk:fXk.append(self.f_value(i))self.fXk = np.array(fXk)#顶点的值self.fmeans = self.fXk.mean()#平均值max_fxk = max(fXk)min_fxk = min(fXk)acc = math.sqrt(sum((np.array(fXk)-min_fxk)**2)/(self.k-1))#self.fXl = min_fxkself.fXH = max_fxkmax_index = fXk.index(max_fxk)min_index = fXk.index(min_fxk)self.XL = self.Xk[min_index]#最小(好)self.XL_index = min_indexself.XH = self.Xk[max_index]#最大(坏)self.XH_index = max_indexdel fXk[max_index]max_1 = max(fXk)self.fXG = max_1max1_index = self.fXk.tolist().index(max_1)self.XG = self.Xk[max1_index]#次大(坏)self.XG_index = max1_indexreturn acc#计算中心点和def center(self):Xcenter = (np.sum(self.Xk,axis=0) - self.XH) /(self.k-1)return Xcenter#计算反射点def reflect(self,XC,a0):XR = XC + a0 * (XC-self.XH)#self.Xk[self.XH_index] = XRreturn XR#扩张a1def again_reflect(self,XR,XC,a1):n = falseXE = XR + a1 * (XR-XC)if self.condition_G(XE):  # 扩张后的点可行if self.f_value(XE) < self.f_value(XR):self.Xk[self.XH_index] = XEn = Truereturn  n#缩小def reduct(self,XC,a2):n = falseXK = self.XH + a2 * (XC - self.XH)if self.condition_G(XK):  # 扩张后的点可行if self.f_value(XK) < self.f_value(self.XH):self.Xk[self.XH_index] = XKn = Truereturn n#压缩def compress(self,XL,Xk,a3):Xk = XL - a3 * (XL - Xk)Xk[self.XL_index] = XLself.Xk = Xkreturn Xk#开始计算#收敛精度delta, 每次缩减的倍数lamda,反射a0 = 1.3, a0的最小限制a0_limit = e - 10，b0，扩张a1 = 1, b1, 缩小a2 = 0.7, b2, 压缩a3 = 0.5def find_(self,delta,lamda,a0,a0_limit):a1 = a0aa = 0bb = 0ci = 0while True:acc = self.fXk_value()#是否为最优解print(f'顶点为Xk={self.Xk}')print(f'XL={self.XL},fXL={self.fXl},XG={self.XG},fXG={self.fXG},XH={self.XH},fXH={self.fXH}')print(f'精度acc={acc}')if  acc<= delta:print(f'最优解为x*={self.XL},极小值为f*={self.fXl}')breakelse:Xcenter = self.center()  # 中心可行否？a0 = a1#每次重新反射都要还原a0while True:if aa == 1:#为了aa = 0  # 退出到accbreak#即用来控制elif self.condition_G(Xcenter)的退出elif self.condition_G(Xcenter) :print(f'可行中心点={Xcenter}')while True:XR = self.reflect(Xcenter,a0)#反射点行否？if self.condition_G(XR) :while True:bb = 0if self.f_value(XR) < self.f_value(self.XH):#是否下降ci += 1self.Xk[self.XH_index] = XRprint(f'反射点XR={XR}')aa = 1break#控制从if self.f_value(XR) < self.f_value(self.XH):的退出else:#不下降if a0 <= a0_limit:self.XH = self.XG###关于XG是否也要换为XH？？交换的效果更好print(f'a0达到极限')break#控制从if a0 <= a0_limit:退出else:a0 *= lamda#print(f'不下降的反射点XR={XR},a0减小，a0={a0}')bb = 1breakif bb == 1:continue#控制从a0 *= lamda#的退出break#从反射点循环退出，控制从if self.condition_G(XR)的退出else:a0 *= lamda#反射点不可行else:#中心点不可行Xa = self.XLXb = Xcenterfor i in range(self.i):self.X0_R[i][0] = Xa[i]self.X0_R[i][1] = Xb[i]self.Xk = self.random_X(self.k,[])print(f'中心点不可行，更新随机区间={self.X0_R}')break#控制else:#中心点不可行的退出print(f'迭代次数k={ci}')
xi = symbols('x1 x2')
F = (xi[0]-5)**2+4*(xi[1]-6)**2
g1 = 64-xi[0]**2-xi[1]**2
g2 = xi[1]-xi[0]-10
g3 = xi[0]-10
G = [g1,g2,g3]
X0_R = [[4,10],[5,15]]
k = 3
X0 = [7,9]
q1 = Complex_method(F,G,xi,X0_R,k,X0)
#收敛精度delta, 每次缩减的倍数lamda,反射a0 = 1.3, a0的最小限制a0_limit = e - 10，b0，扩张a1 = 1, b1, 缩小a2 = 0.5, b2, 压缩a3 = 0.5
q1.find_(1E-9,0.5,1,1E-10)

本算法只设计了反射的情况，对于其他的几种方式有待完善

4.示例

step2:求得最优解为x*=[5.21789682 6.06412585],极小值为f*=0.0639275221179644

总结

更加复杂的复合形方法有待完善

线性与非线性规划：复合形方法相关推荐

线性与非线性规划：黄金分割方法
文章目录前言一.黄金分割法二.进退法 1.进退法介绍 2.代码三.黄金分割代码 1.代码 2.仿真示例总结前言本系列涉及线性与非线性规划中的几种规划算法 1.本节介绍黄金分割(0.618 ...
线性与非线性规划：混合罚函数方法
文章目录前言一.混合罚函数二.算法实现 1.算法步骤 2.算法细节 3.代码 4.示例总结前言本系列涉及线性与非线性规划中的几种规划算法 1.本节介绍混合罚函数方法 2.参考:高立-数值最 ...
线性与非线性规划：随机方向法
文章目录前言一.随机方向法二.算法实现 1.算法步骤 2.算法细节 3.代码示例 4.示例总结前言本系列涉及线性与非线性规划中的几种规划算法 1.本节介绍随机方向法 2.参考:陈宝林-最优 ...
线性与非线性规划：三次拟合
文章目录前言一.三次拟合方法二.代码实现 1.代码 2.示例总结前言本系列涉及线性与非线性规划中的几种规划算法: 1.本节介绍三次拟合方法 2.参考:陈宝林-最优化理论与算法 3.采用py ...
线性与非线性规划——单纯形法python实现
线性与非线性规划--单纯形法python实现单纯形法原理及python实现,时间:2022/1/11 本文不对单纯形方法进行讲解,若有需要请参照<最优化理论与算法第二版>(陈宝林著) ...
【控制】《鲁棒控制-线性矩阵不等式处理方法》-俞立老师-第4章-控制系统综合
第3章回到目录第5章 <鲁棒控制-线性矩阵不等式处理方法>-俞立老师-目录 4.1 H∞H_\inftyH∞ 控制 4.1.1 状态反馈 H∞H_\inftyH∞ 控制 4.1.2 ...
【控制】《鲁棒控制-线性矩阵不等式处理方法》-俞立老师-第3章-系统性能分析
第1章回到目录第3章 <鲁棒控制-线性矩阵不等式处理方法>-俞立老师-目录第1章引言第1章引言 clear clcsyms xi omega k xi = 0.1; omega ...
【控制】《鲁棒控制-线性矩阵不等式处理方法》-俞立老师-第2章-线性矩阵不等式
第1章回到目录第3章 <鲁棒控制-线性矩阵不等式处理方法>-俞立老师-目录第2章引言第2章引言
【控制】《鲁棒控制-线性矩阵不等式处理方法》-俞立老师-第1章-引言
无回到目录第2章跳转链接章节跳转链接第1章引言 <鲁棒控制-线性矩阵不等式处理方法>-俞立老师-目录跳转链接第1章引言第1章引言

线性与非线性规划：复合形方法

文章目录

前言