Data Recovery (UVa Live Archive North America - Southeast - 2010/2011)

http://acmicpc-live-archive.uva.es/nuevoportal/data/problem.php?p=4865

問題

N*Mの行列がある。各成分は[0,100]の値になっているか空欄になっている。各列の和と各行の和を与えるので元の行列を復元せよ。ただし、値が一意に定まらない箇所は-1と答えよ。
1<=N,M<=50

解法

最大流で解く。
二部グラフの片方の集合をy座標ノード、もう一方の集合をx座標ノードとする。その間を結ぶエッジは、行列のその箇所が空欄であれば容量が100のエッジを張る。srcとy座標ノード、x座標ノードとdestの容量は列(行)の和から確定した箇所の数値を引いた物とする。
このグラフ上で最大流を流すと、y座標ノードからx座標ノードへの流量が解となる。あとはこの解が一意に定まるかを調べれば良い。考え方としては、(x,y)の箇所の値を±1した解が存在しなければそこの場所の値は一意に定まる。これは最大流を複数回やれば求めれるが、それだと遅すぎるので、最大流を流した残余グラフ上で(x,y)のエッジを通る閉路が存在するかチェックをすれば良い。
計算量はO(N^2M^2)で計算可能。

ソースコードは続きから

#include <stdio.h>
#include <string.h>
#include <algorithm>
#include <iostream>
#include <math.h>
#include <assert.h>
#include <vector>
#include <queue>

using namespace std;
typedef long long ll;
typedef unsigned int uint;
typedef unsigned long long ull;
static const double EPS = 1e-9;
static const double PI = acos(-1.0);

#define REP(i, n) for (int i = 0; i < (int)(n); i++)
#define FOR(i, s, n) for (int i = (s); i < (int)(n); i++)
#define FOREQ(i, s, n) for (int i = (s); i <= (int)(n); i++)
#define FORIT(it, c) for (__typeof((c).begin())it = (c).begin(); it != (c).end(); it++)
#define MEMSET(v, h) memset((v), h, sizeof(v))

typedef int Weight;
struct Edge {
  int src;
  int dest;
  Weight weight;
  Edge(int src, int dest, Weight weight) : src(src), dest(dest), weight(weight) {;}
  bool operator<(const Edge &rhs) const {
    if (weight != rhs.weight) { return weight > rhs.weight; }
    if (src != rhs.src) { return src < rhs.src; }
    return dest < rhs.dest;
  }
};
typedef vector<Edge> Edges;
typedef vector<Edges> Graph;
typedef vector<Weight> Array;
typedef vector<Array> Matrix;

void PrintMatrix(const Matrix &matrix) {
  for (int y = 0; y < (int)matrix.size(); y++) {
    for (int x = 0; x < (int)matrix[y].size(); x++) {
      printf("%d ", matrix[y][x]);
    }
    puts("");
  }
}

Weight augment(const Graph &g, Matrix &capacity, Matrix &flow, const vector<int> &level, vector<bool> &finished, int from, int t, Weight cur) {
  if (from == t || cur == 0) { return cur; }
  if (finished[from]) { return 0; }
  for (Edges::const_iterator it = g[from].begin(); it != g[from].end(); it++) {
    int to = it->dest;
    if (level[to] <= level[from]) { continue; }
    Weight f = augment(g, capacity, flow, level, finished, to, t, min(cur, capacity[from][to]));
    if (f > 0) {
      capacity[from][to] -= f;
      capacity[to][from] += f;
      flow[from][to] += f;
      flow[to][from] -= f;
      return f;
    }
  }
  finished[from] = true;
  return 0;
}

Weight MaxFlow(const Graph &g, Matrix &capacity, Matrix &flow, int s, int t) {
  int n = g.size();
  capacity = Matrix(n, Array(n));
  flow = Matrix(n, Array(n, 0));
  for (int from = 0; from < n; from++) {
    for (Edges::const_iterator it = g[from].begin(); it != g[from].end(); it++) {
      int to = it->dest;
      capacity[from][to] += it->weight;
    }
  }
  int ans = 0;
  while (true) {
    vector<int> level(n, -1);
    level[s] = 0;
    queue<int> que;
    que.push(s);
    for (int d = n; !que.empty() && level[que.front()] < d; ) {
      int from = que.front();
      que.pop();
      if (from == t) { d = level[from]; }
      for (Edges::const_iterator it = g[from].begin(); it != g[from].end(); it++) {
        int to = it->dest;
        if (capacity[from][to] > 0 && level[to] == -1) {
          que.push(to);
          level[to] = level[from] + 1;
        }
      }
    }
    vector<bool> finished(n);
    bool end = true;
    while (true) {
      Weight f = augment(g, capacity, flow, level, finished, s, t, 2000000000LL);
      if (f == 0) { break; }
      ans += f;
      end = false;
    }
    if (end) { break; }
  }
  return ans;
}

void addEdge(Graph &g, int from, int to, int c) {
  g[from].push_back(Edge(from, to, c));
  g[to].push_back(Edge(to, from, 0));
}

bool visit[110];
bool dfs(const Graph &g, const Matrix &capacity, int from, int target, bool first) {
  if (from == target) { return true; }
  if (first) { visit[from] = true; }
  FORIT(it, g[from]) {
    int to = it->dest;
    if (first && to == target) { continue; }
    if (visit[to] || capacity[from][to] == 0) { continue; }
    visit[to] = true;
    if (dfs(g, capacity, to, target, false)) { return true; }
  }
  return false;
}

int width, height;
bool fix[60][60];
int ans[60][60];
int matrix[60][60];
int row[60];
int colum[60];

inline int Y(int y) { return y; }
inline int X(int x) { return height + x; }
inline int SOURCE() { return width + height; }
inline int DEST() { return width + height + 1; }
inline int SIZE() { return width + height + 2; }

int main() {
  while (scanf("%d %d", &height, &width), width|height) {
    MEMSET(fix, false);
    MEMSET(ans, -1);
    MEMSET(matrix, -1);
    // get input
    REP(y, height) {
      REP(x, width) {
        scanf("%d", &matrix[y][x]);
      }
    }
    REP(y, height) {
      scanf("%d", &row[y]);
    }
    REP(x, width) {
      scanf("%d", &colum[x]);
    }

    // make graph
    Graph g(SIZE());
    REP(y, height) {
      REP(x, width) {
        if (matrix[y][x] != -1) {
          fix[y][x] = true;
          ans[y][x] = matrix[y][x];
          row[y] -= matrix[y][x];
          colum[x] -= matrix[y][x];
        } else {
          addEdge(g, Y(y), X(x), 100);
        }
      }
    }
    REP(y, height) {
      addEdge(g, SOURCE(), Y(y), row[y]);
    }
    REP(x, width) {
      addEdge(g, X(x), DEST(), colum[x]);
    }

    // calc maxflow & answer
    Matrix capacity;
    Matrix flow;
    MaxFlow(g, capacity, flow, SOURCE(), DEST());
    REP(y, height) {
      REP(x, width) {
        if (!fix[y][x]) {
          ans[y][x] = flow[Y(y)][X(x)];
        }
      }
    }

    // check fix
    REP(y, height) {
      REP(x, width) {
        if (fix[y][x]) { continue; }
        MEMSET(visit, false);
        if (capacity[Y(y)][X(x)] >= 1 && dfs(g, capacity, X(x), Y(y), true)) { continue; }
        MEMSET(visit, false);
        if (flow[Y(y)][X(x)] >= 1 && dfs(g, capacity, Y(y), X(x), true)) { continue; }
        fix[y][x] = true;
      }
    }

    // print ans
    REP(y, height) {
      REP(x, width) {
        if (x != 0) { putchar(' '); }
        printf("%d", fix[y][x] ? ans[y][x] : -1);
        //printf("%d", ans[y][x]);
      }
      puts("");
    }
  }
}