SHU运筹与优化上机实验

matlab实验

实验一

1.判断是否存在单位矩阵，并输出列

m=3;n=8; A = 10*rand(m,n);I = eye(m,m);
randIndex = randperm(size(A,2));
A(:,randIndex(1:m))=I;
%准备数据，想了好久这是在干嘛
count=0;
C = nchoosek(randIndex,m);%R = rref(A)
%%%
for i=1:size(C,1)if (A(:,C(i,:))==I)%disp("find eye!");count = count + 1;Column = C(i,:);end
end
%%%
if (count == 0)disp("no eye!");
else    disp("find eye!");Column
end

组合排列函数

nchoosek(n,m)

含义：从n个元素中取出m个元素的所有组合。

将矩阵化成行最简形的命令是rref或rrefmovie。
函数 rref或rrefmovie格式
R = rref(A) %用高斯—约当消元法和行主元法求A的行最简行矩阵R

2.高斯消元解方程

function x=gauss_elim(A,b)
%参量说明：A,系数矩阵；B，常数列向量；zg，增广矩阵
%将增广矩阵化为上三角，再回带求解x
%此方法较为常规，将zg(k,k)元素乘以-zg(i,k)/zg(k,k)加到第i行
%从1:n-1列，主对角元素的以下行，通过两层循环来遍历
zg=[A,b];
zj=rref(zg);
n=length(b);
ra=rank(A);
rz=rank(zg);
temp=rz-ra;
if temp>0disp('无解');return;
end
if ra==rzif ra==nx=zeros(n,1);for p=1:n-1for k=p+1:nm=zg(k,p)/zg(p,p);zg(k,p:n+1)=zg(k,p:n+1)-m*zg(p,p:n+1);endend%用第p层（从1到n-1）消元b=zg(1:n,n+1);A=zg(1:n,1:n);x(n)=b(n)/A(n,n);for q=n-1:-1:1x(q)=(b(q)-sum(A(q,q+1:n)*x(q+1:n)))/A(q,q);endend
end

find（）

查找非零元素的索引和值

k = find(X)
k = find(X,n)
k = find(X,n,direction)

返回一个包含数组 X 中每个非零元素的线性索引的向量。

实验二

实现函数function [xs,Bs,x_num]=BFS(A,b)，求基本可行解和对应的基（矩阵）。

测试程序

%BFS_test.m
A=[2,1,1,0,0;1,1,0,1,0;0,1,0,0,1;]
b=[10,8,7]'[xs,Bs,x_num]=BFS(A,b)

function [xs,Bs,x_num]=BFS(A,b)
%参量说明：A,系数矩阵；b，常数列向量；xs,基本可行解，Bs对应基矩阵，x_num，个数
clc %清屏
%clear all % 清除变量
close all % 关掉图像窗口
format RAT % 显示分数形式%m = 3; %矩阵A的行数
%n =5; %矩阵A的列数% 随机生成一个m*n行满秩矩阵A测试,这里按照ppt给的例子设置
[m,n]=size(A);x_num = 0;idxs = nchoosek(1:n,m); %所有可能的基阵的下标组合
count = size(idxs,1);for k = 1:1:countidx = idxs(k,:);Bs_ = A(:,idx);if det(Bs_) == 0%fprintf(['\n1.B不能构成基阵.\n']); % B是奇异矩阵，不能构成基阵else%fprintf(['\n1. B能构成基阵.\n']); % B不是奇异矩阵，能构成基阵% 令X_N=0，解方程B_X*B+N*X_N=b,计算基本解为X_B=B^(-1)b,XN=0;xs_ = zeros(1,n);X_B = inv(Bs_)*b;X_N = 0;xs_(idx)=X_B;%fprintf('2. 该基阵对应的基本解：xs= [%s].\n',rats(xs));% 判断是否为基本可行解[r,c]=size(xs_);flag=0;%判断基本解是否为基本可行解for i = 1:cif (xs_(1,i)<0)flag=flag+1;endendif (flag == 0)%fprintf(['\n3. xs是基本可行解.\n']);x_num =x_num + 1;xs(x_num,:)=xs_;Bs(:,:,x_num)=Bs_;else%fprintf(['\n3. xs不是基本可行解.\n']);endend
end

实验三

测试程序

%test
A=[ 2 -3 2 1 0;1/3 1 5 0 1];
b=[15 20]';
c=[1 2 1 0 0];format rat %元素使用分数表示
% [x_opt,fx_opt,iter] = Simplex_eye(A,b,c);

函数实现

function [x_opt,fx_opt,iter] = Simplex_eye(A,b,c)
%参量说明：x_opt为最求解，fx_opt为最优函数值，iter为迭代次数
% max z=(c)'x
% Ax=b,x>=0
% r(A)=m,A大小(m*n),(n>=m),b>=0(m*1),c(1*n)这里有问题输入可以是，
%可能是方便输入,直接输入的c其实是是转置后（c'）的行向量
%A是化作标准型后的系数矩阵，
format rat %元素使用分数表示
[m,n]=size(A)%m约束条件个数, n决策变量数
I = eye(m,m);
randIndex = randperm(size(A,2));%A的列数整数的随机排列
C = nchoosek(randIndex,m);%列数的排列组合的集合
ind_B=[];
for i=1:size(C,1)%遍历每一种组合if (A(:,C(i,:))==I)%找到单位矩阵disp("find eye!");ind_B = C(i,:);%找到基变量的列数序号end
end
ind_N = setdiff(1:n, ind_B);  %非基变量的索引，从小到大排序
iter=0;%记录迭代次数%开始循环迭代
while(true)%取基变量，非基变量取值为0%选择初始可行基x0=zeros(n,1);x0(ind_B) = b;%初始基可行解cB = c(ind_B);%计算检验数，选一个检验数最大的非基变量作为换入变量cN=c(ind_N)%注意基变量检验数为0B=A(:,ind_B)%基本矩阵Pjlist=A(:,ind_N)wlist=inv(B)*Pjlist;%wlist=y_kzlist=(cB)*inv(B)*Pjlist;clist=cN-zlist;clist;[maxVal maxInd] = max(clist);%找到换入变量xk，k=maxIndk=maxInd;%判断是否已经是最优解if all(clist(:)<=0)%clist中检验值是否都<=0x_opt = x0;fx_opt = c*x_opt;returnend%判断是否已经是无界if all(wlist(:)<=0)%clist中元素是否都<=0x_opt = [];fx_opt = [];disp('问题无界')returnend%选一个基变量作为换出变量Pk=A(:,k)yk=inv(B)*Pk%     yk(yk<=0)=0.0000001;thetalist=b./yk;%a_ik>0，b一定大于0，如果a_ik<=0就让b_i/a_ik变得无穷大thetalist(thetalist<=0) = 10000;[~,minId] = min(thetalist);%确定出基变量xr索引r,这里写错好多次改好久r=ind_B(minId)% 换基ind_B(ind_B == r) = k;        %新的基变量索引ind_N = setdiff(1:n, ind_B);  %非基变量索引% 更新A和b,化最简型A(:,ind_N) = A(:,ind_B) \ A(:,ind_N)%基矩阵的逆乘以非基矩阵b = A(:,ind_B) \ b %基矩阵的逆乘以bA(:,ind_B) = eye(m,m) %基矩阵更新为单位矩阵iter=iter+1;
end
end

花了好久把单纯形法和单纯形表代码写出来，结果正确！

1.setdiff函数
set difference.

C=setdiff(A,B) for vector A and B, return the values in A that are not in B with no repetitions. C will be sorted.

对于向量A,向量B，C=setdiff(A,B)函数返回在向量A中却不在向量B中的元素，并且C中不包含重复元素，并且从小到大排序。
————————————————
版权声明：本文为CSDN博主「风景不在对岸wj」的原创文章，遵循CC 4.0 BY-SA版权协议，转载请附上原文出处链接及本声明。
原文链接：https://blog.csdn.net/u011089523/article/details/79986511

参考资料

运筹学——matlab实现单纯形法

CSDN一搜全是上大智科的学长，学长太强了。我的代码水平还是不够。

以收集上大智科学长的运筹学博客为乐，这些代码可以一直传承下去。

单纯形法以及对偶单纯形法的Matlab实现

“最近在上《运筹与优化》这门课，讲到了单纯形法这部分，老师让我们上机用Matlab实现单纯形法，数学公式实现过程看不懂的我就暴力模拟了单纯形法的整个过程，包含无界解等情况。”——努力变成大白的小白学长

学长写的代码好复杂，不如上面18级学长写的清楚。

这个17级的学长写了120行，其实只要写50行不到（20行注释），可以年年优化，减少代码行数。

pth_1-utm_source=distribute.pc_relevant.none-task-blog-2%7Edefault%7ECTRLIST%7Edefault-2.highlightwordscore)

“最近在上《运筹与优化》这门课，讲到了单纯形法这部分，老师让我们上机用Matlab实现单纯形法，数学公式实现过程看不懂的我就暴力模拟了单纯形法的整个过程，包含无界解等情况。”——努力变成大白的小白学长

学长写的代码好复杂，不如上面18级学长写的清楚。

这个17级的学长写了120行，其实只要写50行不到（20行注释），可以年年优化，减少代码行数。

纸质作业周四交

SHU运筹与优化上机实验相关推荐

运筹学matlab实验报告,运筹学上机实验报告利用Matlab求解整数线性规划
四川师范大学数学与软件科学学院运筹学上机实验报告. 学期:__2011_至__2012__ 第___一__ 学期 2011年11月9日课程名称:__ 运筹学 ________ 专业:_信息与计算 ...
java编写程序上机实验,《Java程序设计》上机实验
<<Java程序设计>上机实验>由会员分享,可在线阅读,更多相关<<Java程序设计>上机实验(19页珍藏版)>请在技术文库上搜索. 1.tor的安装及 ...
派件系统c语言实验报告,物流规划与优化选址实验报告.doc
物流规划与优化选址实验报告 <物流规划与优化> 课程实验报告一.实验任务与要求重心法是根据待选物流配送中心的数量,将各起迄点预先分配给各个物流配送中心,从而形成个数等于物流配送中心数量 ...
C语言上机报告例文,c语言上机实验报告_大一c语言上机实验报告_c语言实验报告怎么写...
计算机的同学会进行上机实验,包括ERP,JA,C语言等等.下面是出国留学网为大家整理的上机实验心得体会,供大家参考. 上机实验心得体会(一) 通过该实验,对所学的知识有了进一步的了解.在实验的过程中, ...
《汇编语言》上机实验内容//理解
[实验目标要求] <汇编语言>是计算机科学与技术专业必修的专业基础课程.汇编语言程序设计实验的目标是学习汇编语言程序设计的基本方法和技能,熟练掌握用汇编语言设计.编写.调试和运行程序的方法 ...
数值分析上机题matlab线性方程组,数值分析上机实验报告 - 线性方程组部分实验题1...
s=A(i,(i+1):n)*x((i+1):n,1); else s=0; end x(i,1)=(b(i)-s)/A(i,i);end %Cholosky分解方法***************** ...
通信系统计算机仿真上机实验报告,昆明理工大学计算机仿真实验.docx
文档介绍: <计算机仿真>上机实验报告姓名: 学号:-专业:-测控技术与仪器班级:_121-班实验一常微分方程的求解及系统数学模型的转换实验目的通过实验熟悉计算机仿真中常用到的Matl ...
素数c语言程序解题思路,C语言上机实验题目解题思路.doc
上机实验题目解题思路目录第十三次实验:指针之一2 2453:步骤:2 2454:步骤:2 3575:步骤:方法同24543 3576:步骤:3 3580:步骤:3 3582:步骤:3 第十二次实验 ...
西工大matlab计算机实验题,西工大信号系统上机实验一实验二
上机实验1 连续时间信号的时域分析一.实验目的 (1)掌握连续时间信号的时域运算的基本方法: (2)掌握相关函数的调用格式及作用: (3)掌握连续信号的基本运算: (4)掌握利用计算机进行卷积运算的 ...

SHU运筹与优化上机实验

matlab实验

实验一

实验二

实验三

SHU运筹与优化上机实验相关推荐

最新文章

热门文章