线性规划-单纯形法学习笔记

线性规划其实和方程组很像，只不过其中的 $=$ 换成了 $\leq$。

即线性规划就是。

$$\forall i\in[1,m]\ ,\ \sum_{j=1}^{n}a_{i,j}\times x_j\leq c_i\\ \forall i\in [1,n]\ ,\ x_i\geq 0$$

$$\text{最大化}\sum_{j=1}^{n}w_{j}\times x_j$$

不等号不怎么好处理，考虑转成等号。具体来说，可以直接在左式加上一些非负实数，即

$$\forall i\in[1,m]\ ,\ \sum_{j=1}^{n}a_{i,j}\times x_j+x_{n+j}= c_i\\ \forall i\in [1,2n]\ ,\ x_i\geq 0$$

$$\text{最大化}\sum_{j=1}^{2n}w_{j}\times x_j\ (\forall j>n , w_j=0)$$

考虑我们如何处理这个式子：显然，如果我们通过一些等效操作把 $w_j$ 变成非正的，那么最后的答案显然就是常数项。

首先我们考虑如何进行转换操作。为了方便，我们将式子转换一下：

$$\forall i\in[1,m]\ ,\ x_{n+i}=c_i-\sum_{j=1}^{n}a_{i,j}\times x_j$$

这里我们把在要求的式子中左边的变量（一开始是 $x_{n+1\dots 2n}$）称为基，其余的称为非基。

然后我们令 $\forall j\leq n,x_j=0$。那么显然有 $x_{n+i}=c_i$。这就是初值。当然初值可能不符合约束，这个后面会提到。

首先我们先找到一个待换的变量 $x_j$，要求 $w_j<0$，然后我们找到一个可以使答案变大的约束，具体来说，我们应该找到一个 $i$，使 $\frac{c_i}{a_{i,j}}$ 最大。

找到这个数之后，我们考虑换基，即把一个基替换式子中的非基。注意这里需要保证一个式子的左边是基，右边是非基。

这部分转轴可以通过类似于高斯消元的步骤完成。重复上述步骤直到所有 $w_j<0$ 即可。

接下来考虑一件事：如果我们找不到一开始的初值，即一开始 $\exists\ i\leq n, c_i<0$ 怎么办？

最正经的做法是再新建一个变量，然后进行换基。但是这样很麻烦。

考虑一种可能有问题的做法，即每次找到一个 $c_i<0$，然后取一个 $x_j$ 满足 $a_{i,j}>0$，进行转轴操作。

但是这样可能会导致时间炸掉，所以我们每次找的过程中，如果有多个 $x_j>0$，随机取一个即可。

这样UOJ只有97分，不过应该能应付大部分线性规划。

线性规划模板

#include<iostream>
#include<cstdio>
#include<cstring>
#include<cmath>
#include<cstdlib>
#define N 510
#define db double
#define eps 1e-8
using namespace std;
int n,m,id[N],ans[N];
db num[N][N];
void pivot(int x,int y)
{
    swap(id[x+n],id[y]);
    db res=-num[x][y];
    num[x][y]=-1;
    for(int j=0;j<=n;j++) num[x][j]/=res;
    for(int i=0;i<=m;i++)
    if(i!=x && fabs(num[i][y])>eps)
    {
        res=num[i][y],num[i][y]=0;
        for(int j=0;j<=n;j++) num[i][j]+=res*num[x][j];
    }
}
int solve(int t)
{
    for(int i=1;i<=n;i++) id[i]=i;
    while(true)
    {
        int x=0,y=0;
        for(int i=1;i<=m;i++) if(num[i][0]<-eps && (!x || rand()&1)) x=i;
        if(!x) break;
        for(int i=1;i<=n;i++) if(num[x][i]>eps && (!y || rand()&1)) y=i;
        if(!y) return -2;
        pivot(x,y);
    }
    while(true)
    {
        int x=0,y=0;
        db res=0;
        for(int i=1;i<=n;i++) if(num[0][i]>res+eps) res=num[0][i],y=i;
        if(!y) break;
        res=1e10;
        for(int i=1;i<=m;i++) if(num[i][y]<-eps && -num[i][0]/num[i][y]<res-eps) res=-num[i][0]/num[i][y],x=i;
        if(!x) return -1;
        pivot(x,y);
    }
    printf("%.9f\n",num[0][0]);
    for(int i=n+1;i<=n+m;i++) ans[id[i]]=i-n;
    if(t)
    for(int i=1;i<=n;i++)
        printf("%.9f ",ans[i]?num[ans[i]][0]:0);
    return 1;
}
int main()
{
    int t;
    scanf("%d%d%d",&n,&m,&t);
    for(int i=1;i<=n;i++) scanf("%lf",&num[0][i]);
    for(int i=1;i<=m;i++)
    {
        for(int j=1;j<=n;j++) scanf("%lf",&num[i][j]),num[i][j]*=-1;
        scanf("%lf",&num[i][0]);
    }
    int p=solve(t);
    if(p==-1) puts("Unbounded");
    else if(p==-2) puts("Infeasible");
    return 0;
}

CF375E Red and Black Tree

题目大意

给出一棵 $n$ 个节点的树，树上的节点有红黑两种颜色。每次操作可以交换两个节点的颜色，问最少需要多少次操作可以使得树上任意一个点均存在与它距离 $\leq x$ 的黑点，在这里认为树上两个节点的距离为它们之间的最短路径长度。

题解

这题的一个做法是 $O(n^3)$ 的树形dp，但是跑得贼慢。

考虑如何用线性规划解决这个问题，显然有：

$$\sum_{i=1}^{n}x_i=\text{黑点总个数}\\ \forall i\in[1,n]\ ,\ \sum_{dis(i,j)<k}x_j\geq 1 \\ \forall i\in[1,n]\ ,\ x_i\geq 0$$

最后求

$$\sum_{[i\text{是白色}]}{x_i}$$

显然符合线性规划的要求。

注意这里是 $\geq$，所以要把模板中的东西反一下即可。

时间复杂度 $O(\text{能过})$。

#include<iostream>
#include<cstdio>
#include<cstring>
#include<cmath>
#include<cstdlib>
#define N 1010
#define db double
#define eps (1e-9)
using namespace std;
namespace LW{
    int n,m;
    db num[N][N];
    void pivot(int x,int y)
    {
        db res=num[x][y];
        num[x][y]=1;
        for(int j=0;j<=n;j++) num[x][j]/=res;
        for(int i=0;i<=m;i++)
        if(i!=x && fabs(num[i][y])>eps)
        {
            res=num[i][y],num[i][y]=0;
            for(int j=0;j<=n;j++) num[i][j]-=res*num[x][j];
        }
    }
    int solve()
    {
        while(true)
        {
            int x=0,y=0;
            for(int i=1;i<=m;i++) if(num[i][0]<-eps && (!x || rand()&1)) x=i;
            if(!x) break;
            for(int i=1;i<=n;i++) if(num[x][i]<-eps && (!y || rand()&1)) y=i;
            if(!y) return -1;
            pivot(x,y);
        }
        while(true)
        {
            int x=0,y=0;
            db res=0;
            for(int i=1;i<=n;i++) if(num[0][i]>eps){res=num[0][i],y=i;break;}
            if(!y) break;
            res=1e10;
            for(int i=1;i<=m;i++) if(num[i][y]>eps && num[i][0]/num[i][y]<res-eps) res=num[i][0]/num[i][y],x=i;
            if(!x) return -1;
            pivot(x,y);
        }
        printf("%d\n",(int)(num[0][0]+0.5));
        return 1;
    }
}

int nxt[N<<1],to[N<<1],head[N],cnt;
int p[N<<1],val[N];
int n,k;
void add(int u,int v,int w)
{
    nxt[++cnt]=head[u];
    to[cnt]=v;
    p[cnt]=w;
    head[u]=cnt;
}
void dfs(int u,int fa,int dis,int rt)
{
	if(dis>k)return;
    LW::num[rt][u]=-1;
    for(int i=head[u];i;i=nxt[i])
    {
        int v=to[i];
        if(v!=fa) dfs(v,u,dis+p[i],rt);
    }
}
int main()
{
    scanf("%d%d",&n,&k);
    LW::n=n,LW::m=n+2;
    int s=0;
    for(int i=1;i<=n;i++) scanf("%d",&val[i]),s+=val[i];
    for(int i=1;i<n;i++)
    {
        int u,v,w;
        scanf("%d%d%d",&u,&v,&w);
        add(u,v,w);
        add(v,u,w);
    }
    for(int i=1;i<=n;i++) LW::num[0][i]=val[i]-1;
    LW::num[n+1][0]=s;LW::num[n+2][0]=-s;
    for(int i=1;i<=n;i++)
    {
        LW::num[n+1][i]=1;
        LW::num[n+2][i]=-1;
        dfs(i,0,0,i);
        LW::num[i][0]=-1;
    }
	if(LW::solve()==-1) puts("-1");
    return 0;
}