单纯形法的matlab实现(极小化问题)

实验报告

实验题目: 单纯形法的matlab实现学生:

学号:

实验时间: 2013-4-15

一．实验名称: 单纯形法的MATLAB 实现

二．实验目的及要求:

1. 了解单纯形算法的原理及其matlab 实现.

2. 运用MATLAB 编辑单纯形法程序解决线性规划的极小化问题, 求出最优解及目标函数值.

三．实验容:

1. 单纯形方法原理:

单纯形方法的基本思想, 是从一个基本可行解出发, 求一个使目标函数值有所改善的基本可行解; 通过不断改进基本可行解, 力图达到最优基本可行解. 对问题

.

0 ,A s.t. def min ≥=x b x cx f 其中A 是一个m ×n 矩阵, 且秩为m, c 为n 维行向量, x 为n 维列向量, b 为m 维非负列向量. 符号“”表示右端的表达式是左端的定义式, 即目标函数f 的具体形式就是cx . 记

),...,,(n 21p p p A =

令A =(B,N), B 为基矩阵, N 为非基矩阵, 设

⎥⎦

⎤

⎢⎣⎡=0B 1-)

0(b x

是基本可行解, 在)

0(x

处的目标函数值

b c b c c cx

f 1

-B 1-N B )

0(0B 0B ),(=⎥⎦

⎤⎢⎣⎡==,

其中B c 是c 中与基变量对应的分量组成的m 维行向量; N c 是c 中与非基变量对应的分量组成的n-m 维行向量. 现由基本可行解)

0(x 出发求解一个改进的基本可行解.

设⎥⎦

⎤

⎢

⎣⎡=N B x x x 是任一可行解, 则由b Ax =得到 N 1-1-B N B B x b x -=,

在点x 处的目标函数值

∑∈--=⎥⎦

⎤

⎢⎣⎡==R j j j j f x x c c cx f x )c z (),(0N B N B ,

其中R 是非基变量下标集,

j j p c 1-B B z =.

2. 单纯形方法计算步骤:

首先给定一个初始基本可行解, 设初始基为B, 然后执行下列主要步骤: （1）解b x B =B , 求得_

1

b b B x B ==-, 令0=N x , 计算目标函数值B B x

c f =.

（2）求单纯形乘子w , 解B c wB =, 得到1

-=B c w B . 对于所有非基变量, 计算判别数

j j j j j c c -z -=p w . 令

}

c -{z max c -z j j R

j k k ∈=.

若0c -z k k ≤, 则对于所有非基变量0c -z j j ≤, 对应基变量的判别数总是为零, 因此停止计算, 现行基本可行解是最优解. 否则, 进行下一步.

（3）解k k p By =, 得到k 1

k p B y -=, 若0k ≤y , 即k y 的每个分量均非正数, 则停止计算,

问题不存在有限最优解. 否则进行步骤（4）. （4）确定下标r, 使

x k =⎪⎭

⎪⎬⎫⎪⎩⎪⎨⎧>=0min ik ik i k y y b y b r r , r B x 为离基变量, k x 为进基变量. 用k p 替换r B p , 得到新的基矩阵B, 返回步骤（1）.

3. 单纯形方法表格形式:

表 3.1.1

表 3.1.2（3.1.1略去左端列后的详表）

假设0B 1

-≥=b b , 由上表得0,N B ==x b x . 若0c -N B c N -1

B ≤, 则现行基本可行解是最优解.

若0c -N B c N -1

B >, 则用主元消去法求改进的基本可行解. 先根据

}c -{z max c -z j j R j k k ∈=选择主列, 再根据⎪⎭

⎪⎬⎫⎪⎩⎪⎨⎧>=0min ik

ik i k y y b y b r r 找主行, 主元为rk y , 然后进行主元消去, 得到新单纯形表. 表的最后一行是判别数和函数目标值.

四．实验流程图及其MATLAB实现: 1. 流程图:

end

t;

f=y*z;

%%========选取主元==========%%

%---------选取主列---------%

[alpha,q]=max(t);

q;

W(k)=q;%x下标矩阵

%-------------------------%

%--------选取主元----------%

for p=1:m

if B(p,q)<=0

r(p)=N;

else r(p)=z(p)/B(p,q);

end

end

[beta,p]=min(r);

p;

y(p)=C(q);

%-------------------------%

%%==========================%%

B(p,:)=B(p,:)/B(p,q);

for i=1:m

if i~=p

B(i,:)=B(i,:)-B(p,:)*B(i,q);

end

end

if max(t)<=0

break;

end

B;

end

%++++++++++++++++++++++++++++++++++++++% Z=B(:,end);

if length(x(W))~=length(Z)

x=char(' NONE');

f=char(' NONE');

disp(' 不存在有限最优解');

else x(W)=Z';

end