BZOJ3876: [Ahoi2014]支线剧情 线性规划

来源:互联网 发布:微博上最恶心的公知 编辑:程序博客网 时间:2024/06/11 01:08

题意:给一个DAG,每条边有费用,每选一次就要付出这个费用,求从1号点出发的若干条路径将所有边都覆盖的最小费用。
线性规划方程:设每条边的选取次数为变量。
对于每个变量:大于等于1
对于每个点:所有连向这个点的边的经过次数之和大于等于这个点连出的所有边经过次数之和。
最小化:sigma(每条边经过次数*费用)
可以用对偶原理获取初始解,但这样矩阵是m*m的,听说过不了。
对于每个变量x,令x’=x-1,则x’>=0,用x=x’+1替换原式中所有的x,可以将矩阵缩小到n*m,但是没有初始解了。
线性规划解初始解的过程:
首先对于每个小于等于的约束条件,如果我们人为地在左边减去一个数,当这个数很大时所有条件肯定是都能满足的。那么设辅助变量x0(x0>=0),在原线性规划每个式子左边减去x0,得到的线性规划必定有解。
而且若x0可以取0,说明并不需要强行减去一个变量,也就是原线性规划有解。同理,x0大于0说明原线性规划无解。那么我们求-x0最小值。首先找到最小的bl,执行pivot(l,0),就可以将所有的b全变成非负数,也就得到了初始解。接下来解max(-x0)。若结束后x0变为了非基变量,由于其为0,可以全部消去,若是基变量还需找个该行任意的非0变量执行pivot后再消去x0。然后将现在每个变量带入原表达式,再求解最大值即可。
代码:

#include<cstdio>#include<cmath>#include<cstring>#include<algorithm>#include<cassert>using namespace std;const double eps=1e-8;struct simplex{    int n,m;    double a[301][5001],b[301],c[5001],v;    int pos[5301],nex[5001];    int B[301],N[5001];    void pivot(int l,int e)    {        int fr=m+1;        swap(pos[B[l]],pos[N[e]]);        swap(B[l],N[e]);        b[l]/=a[l][e];        for(int i=0;i<=m;++i)        if(i!=e&&fabs(a[l][i])>eps)        {            nex[i]=fr,fr=i;            a[l][i]/=a[l][e];        }        a[l][e]=1/a[l][e];        for(int i=1;i<=n;++i)        if(i!=l&&fabs(a[i][e])>eps)        {            b[i]-=b[l]*a[i][e];            for(int j=fr;j<=m;j=nex[j])            {                a[i][j]-=a[l][j]*a[i][e];            }            a[i][e]=-a[i][e]*a[l][e];        }        v+=b[l]*c[e];        for(int i=fr;i<=m;i=nex[i])        c[i]-=a[l][i]*c[e];        c[e]=-c[e]*a[l][e];    }    void init()    {        for(int i=0;i<=m;++i)        N[i]=pos[i]=i;        for(int i=1;i<=n;++i)        B[i]=pos[m+i]=m+i;        double minn=1e100;        int l=0,e;        for(int i=1;i<=n;++i)        if(b[i]<minn)        minn=b[i],l=i;        if(minn>=-eps) return;        c[0]=-1;        for(int i=1;i<=n;++i)        a[i][0]=-1;        pivot(l,0);        double res=max();        assert(fabs(res)<eps);        if(pos[0]>m)        {            l=pos[0]-m;            e=-1;            for(int i=0;i<=m;++i)            if(fabs(a[l][i])>eps)            {                e=i;break;            }            if(e==-1) return;            pivot(l,e);        }        e=pos[0];        for(int i=1;i<=n;++i)        a[i][e]=0;        memset(c,0,m+1<<3);        v=0;    }    double max()    {        int l,e;        while(1)        {            e=-1;            for(int i=0;i<=m;++i)            if(c[i]>eps)            {                e=i;break;            }            if(e==-1) return v;            double minn=1e100;            for(int i=1;i<=n;++i)            if(a[i][e]>eps&&b[i]/a[i][e]<minn)            {                l=i,minn=b[i]/a[i][e];            }            pivot(l,e);        }    }    void receive(double r[])    {        for(int i=1;i<=m;++i)        {            if(pos[i]<=m) c[pos[i]]+=r[i];            else            {                int l=pos[i]-m;                v+=r[i]*b[l];                for(int j=0;j<=m;++j)                c[j]-=r[i]*a[l][j];            }        }    }}s;int n;int cost[5001];int tot;struct e{    int no;    e *n;    e(int no,e *n):no(no),n(n){}}*f[301],*g[301];double a[5001];bool final[301];int main(){    scanf("%d",&n);    for(int i=1;i<=n;++i)    {        int k,b,t;        scanf("%d",&k);        if(!k) final[i]=1;        while(k--)        {            scanf("%d%d",&b,&t);            cost[++tot]=t;            f[i]=new e(tot,f[i]);            g[b]=new e(tot,g[b]);        }    }    int &cnt=s.n=0;    s.m=tot;    for(int i=2;i<=n;++i)    {        if(final[i]) continue;        ++cnt;        for(e *j=f[i];j;j=j->n)        {            ++s.a[cnt][j->no];            --s.b[cnt];        }        for(e *j=g[i];j;j=j->n)        {            --s.a[cnt][j->no];            ++s.b[cnt];        }    }    s.init();    for(int i=1;i<=tot;++i)    {        a[i]=-cost[i];        s.v+=a[i];    }    s.receive(a);    printf("%.0lf\n",-s.max());    return 0;}
0 0
原创粉丝点击