单纯形 Simplex 学习笔记

前言:这些是一些大神的博客
https://www.cnblogs.com/ECJTUACM-873284962/p/7097864.html
https://ldxcaicai.github.io/LP/


线性规划

标准型:
max(∑i=1ncixi)max(\sum_{i=1}^n c_ix_i)max(i=1ncixi)
limit:∑j=1nai,j∗xj≤bi, x≥0limit:\sum_{j=1}^na_{i,j}*x_j\le b_i,\ x\ge 0limit:j=1nai,jxjbi, x0
松弛型:
max(∑i=1ncixi)max(\sum_{i=1}^n c_ix_i)max(i=1ncixi)
limit:∑j=1nai,j∗xj=bi, x≥0limit:\sum_{j=1}^na_{i,j}*x_j= b_i,\ x\ge 0limit:j=1nai,jxj=bi, x0


标准型与松弛型的转换:
xi+m=bi−∑j=1nai,j∗xjx_{i+m}=b_i-\sum_{j=1}^na_{i,j}*x_jxi+m=bij=1nai,jxj,那么有 ∑j=1nai,j∗xj+xi+m=bi\sum_{j=1}^na_{i,j}*x_j+x_{i+m}= b_ij=1nai,jxj+xi+m=bi

单纯形求解线性规划:

一直重复以下操作:

  1. 在目标函数中找到一个 ce≥0c_e\ge 0ce0eee
  2. 在约束条件中找到一个al,e>0a_{l,e}>0al,e>0且使blal,e\frac{b_l}{a_{l,e}}al,ebl最大的 lll
  3. 将第 lll 个等式进行移项,xe=b[l]−∑j!=ea[l][j]−xl+ma[l][e]x_e=\frac{b[l]-\sum_{j!=e}a[l][j]-x_{l+m}}{a[l][e]}xe=a[l][e]b[l]j!=ea[l][j]xl+m
  4. 带到其他行,消去 xex_exe,新代入 xl+mx_{l+m}xl+m
  5. 当无法找到 eee 时,说明已经最优
  6. 当无法找到合法的 lll 时,说明该线性规划无界

令:当存在 b[i]<0b[i]<0b[i]<0 的时候可能会出锅,可以选择对 b[i]<0b[i]<0b[i]<0 的列随机扰动
具体实现参考代码


满足线性规划问题约束条件的所有点组成的集合就是线性规划的可行域。若可行域有界(以下主要考虑有界可行域),线性规划问题的目标函数最优解必然在可行域的顶点上达到最优。

单纯形法就是通过设置不同的基向量,经过矩阵的线性变换,求得基可行解(可行域顶点),并判断该解是否最优,否则继续设置另一组基向量,重复执行以上步骤,直到找到最优解。所以,单纯形法的求解过程是一个循环迭代的过程


对偶原理:

最大化和最小化可以通过正负来转换,而大于小于号则不行

考虑如下的一个线性规划:

min(∑i=1ncixi)min(\sum_{i=1}^n c_ix_i)min(i=1ncixi)
limit:∑j=1nai,j∗xj≥bi, x≥0limit:\sum_{j=1}^na_{i,j}*x_j\ge b_i,\ x\ge 0limit:j=1nai,jxjbi, x0

利用待定系数法的思想可以知道,该线性规划的任意可行解 {x1,x2,...,xn}\{x_1,x_2,...,x_n\}{x1,x2,...,xn} 一定能表示成如下形式:

∑i=1ncixi≥∑i=1myi(∑j=1nai,j∗xj)≥∑i=1myi∗bi\sum_{i=1}^nc_ix_i\ge \sum_{i=1}^my_i(\sum_{j=1}^na_{i,j}*x_j)\ge\sum_{i=1}^my_i*b_ii=1ncixii=1myi(j=1nai,jxj)i=1myibi
∑i=1myi∗aj,i≤cj,yi≥0\sum_{i=1}^my_i*a_{j,i}\le c_j,y_i\ge 0i=1myiaj,icj,yi0

于是可以转换为一个新的线性规划
找到 ∑i=1myi∗bi\sum_{i=1}^my_i*b_ii=1myibi 的上界就相当于找到 ∑i=1nci∗xi\sum_{i=1}^nc_i*x_ii=1ncixi
但是我们还要证明当 ∑i=1myi∗bi\sum_{i=1}^my_i*b_ii=1myibi 取到最优解即最大值时,∑i=1nci∗xi\sum_{i=1}^nc_i*x_ii=1ncixi 也取到最优解即最小值

证明:
可以参考
https://blog.youkuaiyun.com/itnerd/article/details/83352441

将线性规划表示成矩阵的形式就是要证这么一个东西:
min{cTx∣Ax≥b,x≥0}⇔max⁡{bTy∣ATy≤c,y≥0}min\{c^Tx|Ax\ge b,x\ge 0\}\Leftrightarrow \max\{b^Ty|A^Ty\le c,y\ge 0\}min{cTxAxb,x0}max{bTyATyc,y0}

对于第一个线性规划,构造一个函数 f(y)=cTx+yT(b−Ax),y≥0f(y)=c^Tx+y^T(b-Ax),y\ge 0f(y)=cTx+yT(bAx),y0
g(y)=min{f(y)∣Ax≥b,x≥0}g(y)=min\{f(y)|Ax\ge b,x\ge0\}g(y)=min{f(y)Axb,x0}
∀g(y)≤cTx\forall g(y)\le c^Txg(y)cTxmax(g(y))=cTxmax(g(y))=c^T xmax(g(y))=cTx
于是有
min(cTx)=max(g(y))min(c^Tx)=max(g(y))min(cTx)=max(g(y))
=max(yTb+min((cT−yTA)x))=max(y^Tb+min((c^T-y^TA)x))=max(yTb+min((cTyTA)x))
=max(bTy+min((cT−yTA)x))=max(b^Ty+min((c^T-y^TA)x))=max(bTy+min((cTyTA)x))

cT−yTA≥0c^T-y^TA\ge 0cTyTA0,那么 min((cT−yTA)x)=0min((c^T-y^TA)x)=0min((cTyTA)x)=0
反之 min((cT−yTA)x)=−∞min((c^T-y^TA)x)=-\inftymin((cTyTA)x)=
所以
=max(bTy+min((cT−yTA)x))=max(bTy)=max(b^Ty+min((c^T-y^TA)x))=max(b^Ty)=max(bTy+min((cTyTA)x))=max(bTy)
令:翻到一个感性的证明,有点像经济学上买家卖家的博弈
在这里插入图片描述

线性规划互相松弛定理

在这里插入图片描述


这里给出 NOI 2008 志愿者招募的题解及代码

ai,ja_{i,j}ai,j 表示第 iii 类志愿者在第 jjj 天是否工作,xix_ixi 表示第 iii 类志愿者的个数

那么有线性规划

min(∑i=1mxi∗ci)min(\sum_{i=1}^mx_i*c_i)min(i=1mxici)
limit:∑j=1naj,i∗xj≥bi, xi≥0limit:\sum_{j=1}^na_{j,i}*x_j\ge b_i,\ x_i\ge 0limit:j=1naj,ixjbi, xi0

对偶过后发现
max(∑i=1nyi∗bi)max(\sum_{i=1}^ny_i*b_i)max(i=1nyibi)
limit:∑j=1mai,j∗yj≤cilimit:\sum_{j=1}^ma_{i,j}*y_j\le c_ilimit:j=1mai,jyjci
已经是标准形式了

#include<bits/stdc++.h>
#define cs const
using namespace std;
int read(){
	int cnt = 0, f = 1; char ch = 0;
	while(!isdigit(ch)){ ch = getchar(); if(ch == '-') f = -1; }
	while(isdigit(ch)) cnt = cnt*10 + (ch-'0'), ch = getchar();
	return cnt * f;
}
cs int N = 1e3 + 5, M = 1e4 + 5;
cs double eps = 1e-8, INF = 1e15;
int n, m; double a[M][N];
void pivot(int l, int e){
	double t = a[l][e]; a[l][e] = 1;
	for(int i = 0; i <= n; i++) a[l][i] /= t;
	for(int i = 0; i <= m; i++) if(l!=i && fabs(a[i][e]) > eps){
		t = a[i][e]; a[i][e] = 0;
		for(int j = 0; j <= n; j++) a[i][j] -= a[l][j] * t;
	}
}
bool simplex(){
	while(true){
		int l = 0, e = 0; double mi = INF;
		for(int i = 1; i <= n; i++) if(a[0][i] > eps) { e = i; break; } if(!e) break;
		for(int i = 1; i <= m; i++) if(a[i][e] > eps && a[i][0] / a[i][e] < mi) mi = a[i][0] / a[i][e], l = i;
		if(!l) return false; pivot(l, e);
	} return true;
}
int main(){
	n = read(), m = read();
	for(int i = 1; i <= n; i++) a[0][i] = read();
	for(int i = 1; i <= m; i++){
		int l = read(), r = read(), c = read();
		for(int j = l; j <= r; j++) a[i][j] = 1;
		a[i][0] = c;
	} if(simplex()) cout << (long long)(-a[0][0] + 0.5);
	return 0;
}

附:随机扰动的代码

    bool init(){
        while(true){
            int l = 0, e = 0;
            for(int i = 1; i <= m; i++)
                if(a[i][0] < -eps && (!l || (rand()&1))) l = i;
            if(!l) break;
            for(int j = 1; j <= n; j++) 
                if(a[l][j] < -eps && (!e || (rand()&1))) e = j;
            if(!e) return false; pivot(l, e);
        } return true;
    }
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包

打赏作者

FSYo

你的鼓励将是我创作的最大动力

¥1 ¥2 ¥4 ¥6 ¥10 ¥20
扫码支付:¥1
获取中
扫码支付

您的余额不足,请更换扫码支付或充值

打赏作者

实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值