火星探险问题

题意：

​ 火星的部分地形可以用一个 P*Q 的网格表示。登陆舱位于方格 (1,1)，目标点位于方格 (P,Q)。现有 k 辆探测车从登陆舱出发前往目标点。探测车的坐标不能减小 (即 x,y 均不降)，且不能到达有障碍的地面。探测车可以收集经过的岩石样本，且每个方格的样本只能收集一次。求每个探测车的路线 (可以有重合部分) 使到达目标的车最多的同时这些车收集到样本最多。
​ 输入为 k,P,Q 以及一副地图，0 表示空地，2 表示岩石样本，1 表示障碍。

​ k,P,Q<=35

分析：

​ 对于首要最优化问题：到达目标的车最多，我们可以直接采用最大流。
​ 为了满足次要最优化问题：到达的车收集的样本最多，我们将最大流改进为最大费用最大流即可。

​ 但是一个严峻的问题是：此题需要输出每个探测车的路线。对于这个问题，我们可以从以下几个方面思考，找到合理的方法。
​ (1) 在增广的同时记录每一次的路线。
​ (2) 尝试从最后的残量网络上找路线。
​ (3) 放弃这道题，去做 “深海机器人问题”。

​ 第一个方法很容易被否决了。因为增广时我们会遇到反向边，这样记录的每条路线都不是真正的探测车道路。如果可以利用什么方法将新的増广路参考旧的修改出真正对应的线路，这个方法依旧可行。但是太麻烦。
​ 第三个方法也被否决了。轻易放弃并不是我们应有的精神。所以我们还是在方法 (2) 的基础上加以分析。
​ 对于第二个方法，我们先明确这样一个事实：最后的残量网络中正向边的剩余流量 $\not =$容量，则这条边被某辆车经过了。根据这个事实，我们可以想出以下结论：(1) 容量和剩余流量的差值为经过这条边的车数 (2) 一定可以找到[ 流量] 条路径覆盖每条被车经过了的边。由于我们并不在乎每辆车的身份，所以两辆车经过同一条边后分道扬镳时，我们可以随便指定其中一辆车的航向。那么，我们只需要随便找到一条路径，组成它的每一条边都被车经过，将这条路径输出，同时路径上每一条边经过的车数减一，继续重复，就可以得到[ 流量] 条路径了。

代码：

#include <iostream>
#include <cstring>
#include <cstdio>

#define MX 40004
#define oo 12312312
#define S 0
#define T 10001

using namespace std;

typedef struct edge_t
{
    int u,v,c,w;
}edge;
edge e[MX];
int fst[MX],nxt[MX],lnum;
int mp[101][101];
int n,m,k;

inline void addeg(int nu,int nv,int nc,int nw)
{
    nxt[++lnum]=fst[nu];
    fst[nu]=lnum;
    e[lnum]=(edge){nu,nv,nc,nw};
}

void init()
{
    memset(fst,0xff,sizeof(fst));
    lnum=-1;
}

void input()
{
    scanf("%d%d%d",&k,&m,&n);
    for(int i=1;i<=n;i++)
        for(int j=1;j<=m;j++)
            scanf("%d",&mp[i][j]);
}

int getp(int l,int x,int y)
{
    return l*m*n+(x-1)*m+y;
}

void build()
{
    int nx,ny;
    for(int i=1;i<=n;i++)
    {
        for(int j=1;j<=m;j++)
        {
            if(mp[i][j]==0)addeg(getp(0,i,j),getp(1,i,j),k,0),addeg(getp(1,i,j),getp(0,i,j),0,0);
            else if(mp[i][j]==2)
            {
                addeg(getp(0,i,j),getp(1,i,j),1,1),addeg(getp(1,i,j),getp(0,i,j),0,-1);
                addeg(getp(0,i,j),getp(1,i,j),k,0),addeg(getp(1,i,j),getp(0,i,j),0,0);
            }
            if(i+1<=n)addeg(getp(1,i,j),getp(0,i+1,j),k,0),addeg(getp(0,i+1,j),getp(1,i,j),0,0);
            if(j+1<=m)addeg(getp(1,i,j),getp(0,i,j+1),k,0),addeg(getp(0,i,j+1),getp(1,i,j),0,0);
        }
    }           //if(y==x-n*m+1||y==x-n*m+m)outp(x),cout<<"->",outp(y),cout<<endl;
    addeg(S,getp(0,1,1),k,0),addeg(getp(0,1,1),S,0,0);
    addeg(getp(1,n,m),T,k,0),addeg(T,getp(1,n,m),0,0);
}

int dis[MX],inq[MX],pre[MX],q[MX];
int spfa(int *flow,int *cost)
{
    int x,y,h=0,t=1,mfl=+oo;
    memset(dis,0xaf,sizeof(dis));
    memset(pre,0,sizeof(pre));
    dis[S]=0;
    inq[S]=1;
    pre[S]--;
    q[++h]=S;
    while(h>=t)
    {
        x=q[(t++)%MX];
        inq[x]=0;
        for(int i=fst[x];i!=-1;i=nxt[i])
        {
            y=e[i].v;
            if(e[i].c&&dis[y]<dis[x]+e[i].w)
            {
                dis[y]=dis[x]+e[i].w;
                pre[y]=i;
                if(!inq[y])
                {
                    q[(++h)%MX]=y;
                    inq[y]=1;
                }
            }
        }
    }
    if(dis[T]<-oo)return 0;
    x=pre[T];
    while(x!=-1)mfl=min(mfl,e[x].c),x=pre[e[x].u];
    x=pre[T];
    while(x!=-1)e[x].c-=mfl,e[x^1].c+=mfl,x=pre[e[x].u];
    *flow+=mfl;
    *cost+=mfl*dis[T];
    return 1;
}

int mcf()
{
    int flow=0,cost=0;
    while(spfa(&flow,&cost));
    return cost;
}

void outp(int x)
{
    cout<<(((x-1)%(n*m))/m+1)<<" "<<((x-1)%m+1);
}

void out(int x,int car)
{
    if(x==getp(0,n,m))return;
    for(int i=fst[x];i!=-1;i=nxt[i])
    {
        int y=e[i].v;
        if((i&1)==0&&e[i].c!=k)
        {
            if(y==x-n*m+1)cout<<car<<" "<<1<<endl;
            else if(y==x-n*m+m)cout<<car<<" "<<0<<endl;
            e[i].c++;
            out(y,car);
            return ;
        }
    }
}

int main()
{
    init();
    input();
    build();
    mcf();
    for(int i=1;i<=k;i++)out(getp(0,1,1),i);
    return 0;
}