单纯型法MATLAB编程求解
单纯型法编程求解,本程序主要参考了大佬博主“天涯铭”在其博文《MATLAB多维极值之单纯形法》中的编程,个人只是根据自己需要进行了些许修改,比如加入了一些矩阵用于记录迭代运行数据,然后删除了一小段程序段,但整体并没有过多的创新。
大家可以去“天涯铭”大佬的博客去看看,干货很多的!!
附大佬博文连接如下:https://blog.csdn.net/STM89C56/article/details/105725205
写这篇博文是因为一些相关的程序编程用到了个人改动后的函数作为子函数运行,在个人其他博文介绍时会涉及到,所以附程序如下,分享出来方便大家交流。
clear all;clc
%% 初始参数设定,求f最小值
syms x1 x2 x
f=x1^2+2*x2^2-4*x1-2*x1*x2;% f 函数符号表达式
x0=[0 0];% 给定初始点
h=4; % 初始单纯形构造因子
alpha=1; % alpha 反射因子
gama=2; % gama 扩张因子
beta=0.5;% beta 收缩因子
eps=0.5;% eps 精度值
%% 构造初始单纯形
x0=x0';
X(:,1)=x0;
Num_var_f=length(symvar(f));%二维函数
for i = 1:Num_var_f %二维函数,需要三个顶点坐标
x = zeros(Num_var_f,1);%预分配内存
x(i) = x0(i) + h;
X(:,i+1)=x;% 将单纯形存储
end
n = 0; %迭代循环次数
m = 1; %精度差值
k = 0; %数据记录
g = 0; %操作记录
%% 计算函数值
for i = 1:Num_var_f+1
F_val(i)=double(subs(f,symvar(f),X(:,i)'));
end
[F_sort,F_index]=sort(F_val);%从小到大排列函数
% 最好点
f_best_0=F_sort(1);
X_best_0=X(:,F_index(1));
% 最差点
f_bad_n=F_sort(end);
X_bad_n=X(:,F_index(end));
% 次差点
f_nextbad_1n=F_sort(end-1);
X_nextbad_1n=X(:,F_index(end-1));
% 计算形心
Xc_n1=1/Num_var_f*sum(X(:,F_index(1:end-1)),2);
f_Xc_n1=double(subs(f,symvar(f),Xc_n1'));
%% 循环迭代语句
while m > eps
%% 反射
flag = 0;
X_reflect_n2=Xc_n1+alpha*(Xc_n1-X_bad_n);
f_reflect_n2=double(subs(f,symvar(f),X_reflect_n2'));
%% 反射成功
if f_reflect_n2<f_best_0 % 反射点R<最好点
%% 进行扩张
X_reflect_expand_n3=Xc_n1+gama*(Xc_n1-X_bad_n);%扩张
f_reflect_expand_n3=double(subs(f,symvar(f),X_reflect_expand_n3'));
if f_reflect_expand_n3 < f_best_0%扩张点继续缩小
X(:,F_index(end)) = X_reflect_expand_n3; %将其作为单纯型的新顶点
g=1;%扩张点
else
X(:,F_index(end)) = X_reflect_n2; %采用扩张前的点作为单纯形的新顶点
g=0;%反射点
end
flag = 1; %反射成功标志位
end
%% 反射失败
if flag == 0 %说明没有反射成功
if f_reflect_n2 < f_nextbad_1n %反射点<次差点
X(:,F_index(end)) = X_reflect_n2;%采用该反射点作为新的顶点
g=0;%反射点
%% 进行收缩
else % f_reflect >= f_nextbad
X_reflect_shrink_n4=Xc_n1+beta*(X_reflect_n2-Xc_n1);% 收缩点
f_reflect_shrink_n4=double(subs(f,symvar(f),X_reflect_shrink_n4'));
X(:,F_index(end)) = X_reflect_shrink_n4;
g=2;%收缩点
end
end
%% 计算迭代完后函数序列,并计算精度
for i = 1:Num_var_f+1
F_val(i)=double(subs(f,symvar(f),X(:,i)'));
end
[F_sort,F_index]=sort(F_val);%从小到大排列函数
% 最好点
f_best_0=F_sort(1);
X_best_0=X(:,F_index(1));
% 最差点
f_bad_n=F_sort(end);
X_bad_n=X(:,F_index(end));
% 次差点
f_nextbad_1n=F_sort(end-1);
X_nextbad_1n=X(:,F_index(end-1));
% 计算形心
Xc_n1=1/Num_var_f*sum(X(:,F_index(1:end-1)),2);
f_Xc_n1=double(subs(f,symvar(f),Xc_n1'));
% 计算精度与迭代次数
m=abs(double(subs(f,symvar(f),X(:,F_index(end))'))-f_Xc_n1);
n=n+1;
% 运算数据记录
k=k+1;
A(k)=g;
B(:,k+1)=X(:,F_index(end))';
M(k)=m;
end
%% 结果输出
minx=X_best_0;
A=A';B=B';M=M';
A
B
M
n
minx
m