码迷,mamicode.com
首页 > 其他好文 > 详细

单纯形法与线性规划入门

时间:2020-07-18 22:38:37      阅读:62      评论:0      收藏:0      [点我收藏+]

标签:最大化   ret   register   href   simple   线性规划   and   code   c++   

前言

说实话其实这个算法真的不难,只要你有耐心把它啃完。。。

线性规划

线性规划的标准型长这个样子:

\[\begin{align}max\ \ &\sum_{j=1}^nc_jx_j&\\s.t.\ \ &\sum_{j=1}^na_{i,j}x_j\le b_i&i=1,2,...,m\\&x_j\ge0&j=1,2,...,n\end{align} \]

它的松弛型则是长这个样子:

\[\begin{align}max\ \ &\sum_{j=1}^nc_jx_j&\\s.t.\ \ &x_{n+i}=b_i-\sum_{j=1}^na_{i,j}x_j&i=1,2,...,m\\&x_j\ge0&j=1,2,...,n+m\end{align} \]

考虑和原本的标准型相比,松弛型有什么优势。

容易发现,原本限制条件中的不等号,被我们转化成了等号。

这一步转化是之后所有操作的基础。

\(Pivot\)(转轴)

考虑到我们需要最大化的式子\(\sum_{j=1}^nc_jx_j\)

如果我们能够通过某些手段把\(x\)的系数全部都变成负数,那么因为\(x_j\ge0\),只要令\(x_{1\sim n}=0\)即可得到最大值。

但我们该怎样把\(x\)的系数都变成负数呢?

其实,只要通过等量代换就可以了。

考虑\(x_j\ge 0\)并不是仅针对\(x_{1\sim n}\),而同时要求\(x_{n+1\sim n+m}\)也满足这个条件。

而它们之间存在一定的等量关系,就可以试着通过等量代换来调整所需最大化的式子。

我们规定把\(x_e\)\(x_{n+l}\)代换的过程叫做\(Pivot(l,e)\)

已知:

\[x_{n+l}=b_l-\sum_{j=1}^na_{l,j}x_j \]

我们移出右式中的\(x_e\),将它和\(x_{n+l}\)调换位置得到:

\[x_e=\frac{b_l}{a_{l,e}}-\frac1{a_{l,e}}x_{n+l}-\sum_{j=1}^n[j\not=e]\frac{a_{l,j}}{a_{l,e}}x_{j} \]

然后我们再把其余式子中的\(x_e\)全部用右式代换,就彻底抹消了\(x_e\)的存在,完成了一次转轴。

\(Init\)(初始化)

考虑如果存在某个\(b_i<0\),我们就无法令\(x_{1\sim n}=0\),因为这样就不满足限制了。

仔细观察转轴操作,可发现这一操作过程中并不会改变\(b_i\)的正负性,因此只要我们在一开始把所有\(b_i\)都强制修改为正数即可。

为了达成这一目标,每次我们找出一个\(b_l<0\),并找到一个\(a_{l,e}<0\),然后只要执行\(Pivot(l,e)\)即可。

如果找不到,显然无解。

\(Simplex\)(单纯形法)

单纯形法主过程的核心问题,就是每次应该对哪一对\(l,e\)执行转轴操作:

  • 首先,需要满足\(c_e>0\),且转轴后可以让新的\(c_e<0\)(这是我们执行转轴操作的初衷)。

  • 其次,在所有可选的\(l\)中,我们要使\(\frac{b_l}{a_{l,e}}\)最小(要选择最紧的限制)。

只要不断进行\(Pivot(l,e)\),直到无法继续就可以了。

代码

【UOJ179】线性规划(据说因为出题人卡常卡精度,只要拿到\(97\)分就可以了?)

#include<bits/stdc++.h>
#define Tp template<typename Ty>
#define Ts template<typename Ty,typename... Ar>
#define Reg register
#define RI Reg int
#define Con const
#define CI Con int&
#define I inline
#define W while
#define N 20
#define DB double
#define eps 1e-8
using namespace std;
int n,m,op;DB a[N+5][N+5];
namespace SimplexMethod//单纯形法
{
	int id[2*N+5];DB s[N+5];
	I void Pivot(CI l,CI e)//转轴
	{
		RI i,j;DB t=a[l][e];swap(id[n+l],id[e]);//交换编号
		for(a[l][e]=1,i=0;i<=n;++i) a[l][i]/=t;//对这个式子变形
		for(i=0;i<=m;++i) if(i^l&&fabs(a[i][e])>eps)//枚举其他式子
			for(t=a[i][e],a[i][e]=j=0;j<=n;++j) a[i][j]-=t*a[l][j];//等量代换
	}
	I bool Init()//初始化,使所有b[i]为正
	{
		RI i,l,e;W(1)//不断操作
		{
			for(l=e=0,i=1;i<=m;++i) a[i][0]<-eps&&(!l||rand()&1)&&(l=i);if(!l) break;//找到b[l]<0
			for(i=1;i<=n;++i) a[l][i]<-eps&&(!e||rand()&1)&&(e=i);//找到a[l][e]<0
			if(!e) return puts("Infeasible"),0;Pivot(l,e);//找不到则无解,否则转轴
		}return 1;
	}
	I bool Simplex()//单纯形法主过程
	{
		RI i,l,e;DB Mn;W(1)//不断操作
		{
			for(l=e=0,Mn=1e9,i=1;i<=n;++i) if(a[0][i]>eps) {e=i;break;}if(!e) break;//找到c[e]>0
			for(i=1;i<=m;++i) a[i][e]>eps&&a[i][0]/a[i][e]<Mn&&(Mn=a[i][0]/a[i][e],l=i);//找到限制最紧的l
			if(!l) return puts("Unbounded"),0;Pivot(l,e);//如果所有a[l][e]<0无穷大,否则转轴
		}return 1;
	}
	I void Solve()//单纯形法
	{
		RI i;for(srand(20050521),i=1;i<=n;++i) id[i]=i;if(Init()&&Simplex())
		{
			if(printf("%.10lf\n",-a[0][0]),!op) return;//输出答案
			for(i=1;i<=m;++i) s[id[n+i]]=a[i][0];for(i=1;i<=n;++i) printf("%.10lf ",s[i]);//输出方案
		}
	}
}
int main()
{
	RI i,j;for(scanf("%d%d%d",&n,&m,&op),i=1;i<=n;++i) scanf("%lf",a[0]+i);//读入
	for(i=1;i<=m;++i) {for(j=1;j<=n;++j) scanf("%lf",a[i]+j);scanf("%lf",a[i]);}//读入
	return SimplexMethod::Solve(),0;//单纯形法
}

单纯形法与线性规划入门

标签:最大化   ret   register   href   simple   线性规划   and   code   c++   

原文地址:https://www.cnblogs.com/chenxiaoran666/p/SimplexMethod.html

(0)
(0)
   
举报
评论 一句话评论(0
登录后才能评论!
© 2014 mamicode.com 版权所有  联系我们:gaon5@hotmail.com
迷上了代码!