Pages

Showing posts with label Network Flow. Show all posts
Showing posts with label Network Flow. Show all posts

Tuesday, September 15, 2009

Maximum/Minimum Weighted Bipartite Matching using Cycle Cancelling

Problem:

This is an extension to our maximum cardinality bipartite matching problem we introduced earlier. Imagine the same situation, we are given a bipartite graph G = (V,E) in which the vertices can be separated into two disjoint sets such that there are no edges between vertices belonging in the same set. Instead of weightless edges like our previous problem - we are now faced with weighted edges. The problem now is to match the vertices such that the total weight of the matched edges is as small (or as large) as possible whilst still retaining the maximum number of matchings. The figure below depicts this situation:



Applications:

Like the weightless version, we can compute the most optimal pairings between two disjoint sets. For example, let X represent a set of workers and Y to be a set of jobs. There is an edge (u,v) in E of worker u in X to job v in Y of weight W, which gives the productivity of the worker u given job v. If we were tasked to assign the workers to the set of jobs such that each worker is given exactly 1 job and we were also told to maximise the total productivity of the assignments then the answer is simply the maximum weighted matching in the bipartite graph. We can also define the edge weights to mean the converse, for example, the weight for edge (u,v) could mean how difficult worker u would find job v – in which our task could be to minimise the total difficulty whilst ensuring all workers are assigning to exactly 1 job. This does not complicate things as we simply need to reverse our operations to derive the opposite result.

Algorithm:

We need a new approach to tackle this problem. Augmenting the paths isn't enough as it doesn't take into account the weightings on the edge. One idea is to arbitrarily assign the maximum matching without consideration of the weights and then try to successively improve (i.e. maximise) the cost by finding a negative cycle an augmenting the matchings along this cycle. To represent this we calculate a modified graph G’, with edge weights representing the cost difference we will gain by matching a vertex x with another vertex y. In other words, for any u and v, we put a weight W on the edge from u to v, where W is defined as the cost we will gain if we gave the current partner of v to u, i.e. W = cost(u, matching[v]) – cost(u, matching[u]). Now we simply look for a negative cycle in G’ and augment the partners along the cycle to produce a better matching. How does this work? If we can find a negative cycle from a vertex u to itself, it means we can alternate partners to reduce the total cost whilst also keeping the same number of matched partners.

Consider the following bipartite graph:



And the matchings below:

 

The matching on the left has a weight of 15; the optimal minimum matching has a weight of 8 (shown on the right). Algorithmically, let’s assign the vertices on the top-left and bottom-left (A and B respectively) and the vertices on the top-right and bottom-right (C and D respectively). We compute the modified G’ edge weights of (A,B) and (B,A) – the cost of switching the partners around (take note their matched counterparts – i.e. the right side, is calculated implicitly). Proceeding to do this we yield:

W(A,B) = cost(A, matching[B]) – cost(A, matching[A]) = 6 – 5 = 1
W(B,A) = cost(B, matching[A]) – cost(B, matching[B]) = 2 – 10 = -8

Therefore a possible cycle could be from A -> matching[B] -> B -> matching[A] -> A. This yields a cost of W(A,B) + W(B,A) = 1 + (-8) = -7. Since this is a negative cycle we can save a total cost of up to 7 (you can verify this by hand) if we augment along this cycle. From this, A becomes matched to matching[B] (i.e. D) and B becomes matched to matching[A] (i.e. C). The matching on the right hand side diagram highlights this minimal matching. To actually find a negative cycle, an easy way is to use Floyd-Warshall’s algorithm which is illustrated in our implementation below. Additional book-keeping is required to keep track of the parents so we can back-track and augment the negative cycle path.

To summarise, our algorithm is as follows:

- Start with an initial perfect matching M
- While there is a negative cycle C on G’, augment along the negative cycle C to produce a better matching M’
- Return the matching M

Other Notes:

This problem can also be solved using the Hungarian algorithm (a famous combinatorial optimisation algorithm) or using minimum-cost flows. In fact, our cycle cancelling implementation is a specific subset of the minimum cost flow algorithm. The difference between our cycle cancelling algorithm and the one used for the general minimum cost flow problem is that we need to run a maximum flow algorithm to derive an initial network (here we simply assign the vertices to each other as we are guaranteed a perfect matching it’s a complete bipartite graph).

There are many variations on maximum (or minimum) weighted bipartite matching. The first variation is called the assignment problem where we are given an equal number of vertices on each side and the graph itself is complete (i.e. all the vertices link to each other between the two disjoint sets). The second variation we are given an uneven number of vertices on each side but the graph itself is still complete, here we proceed to reduce the problem into the first variation by adding dummy vertices to make the sides even. Any edges which go to these vertices have a weight of 0 and hence do not have an effect on the final answer. Another variation is where there are an uneven number of vertices and the graph isn’t complete – hence a perfect matching may not be possible, this is where the general minimum cost flow algorithm comes in which we will go over later.

Applying the Algorithm:

We briefly demonstrate how the algorithm is used to solve two algorithmic problems.


ACM ICPC Pacific Northwest Regionals 2004 (GoingHome)

URL: http://acm.tju.edu.cn/toj/showp1636.html

The problem can be easily reduced to the assignment problem by constructing a bipartite graph of houses and people. The edge weights are simply the Manhattan distance between house i and person j. We then just run our algorithm over the graph and we obtain our answer!

#include <iostream>
#include <string>
#include <vector>
#include <algorithm>

using namespace std;

int cost[111][111];      // cost of assigning house i to person j
int curMatchings[111];   // house i matched to person curMatchings[i]
int delta[111][111];   // the delta/cost graph
int parent[111][111];   // used for backtracking the cycle
int cycle[111];         // path of the cycle
int totalLen;         // length of the cycle

bool augmentingCycle(int N) {
   memset(delta,0,sizeof(delta));
   memset(parent,0,sizeof(parent));
   memset(cycle,-1,sizeof(cycle));
   // derive the delta graph
   for (int i = 0; i < N; i++) {
      for (int j = 0; j < N; j++) {
         delta[i][j] = cost[i][curMatchings[j]] - cost[i][curMatchings[i]];
         parent[i][j] = j;
      }
   }
   for (int k = 0; k < N; k++) {
      for (int i = 0; i < N; i++) {
         for (int j = 0; j < N; j++) {
            if (delta[i][k] + delta[k][j] < delta[i][j]) {
               // minimum cost matching
               delta[i][j] = delta[i][k] + delta[k][j];
               // book-keep the optimal path so far for i->j
               parent[i][j] = parent[i][k];
               // detect for cycle (if vertex i and j are the same vertex)
               if (i == j) {
                  totalLen = 0;
                  // backtrack and construct the cycle path
                  do {
                     cycle[totalLen++] = i;
                     i = parent[i][j];
                  } while (i != j);
                  return true;
               }
            }
         }
      }
   }
   // did not find an augmenting path
   return false;
}

int cycleCancel(int N) {
   int res = 0;
   // pseudo max-flow for assignment problem
   for (int i = 0; i < N; i++) curMatchings[i] = i;
   while (augmentingCycle(N)) {
      // augment the negative cycle
      int r = curMatchings[cycle[0]];
      for (int i = 0; i < totalLen - 1; i++) {
         curMatchings[cycle[i]] = curMatchings[cycle[i+1]];
      }
      curMatchings[cycle[totalLen-1]] = r;
   }
   // compute the final cost
   for (int i = 0; i < N; i++) {
      res += cost[i][curMatchings[i]];
   }
   return res;
}

int main() {
   int N, M;
   while (cin >> N >> M) {
      if (N == 0 && M == 0) break;
      char ch;
      vector<pair<int,int> > houses;
      vector<pair<int,int> > people;
      for (int i = 0; i < N; i++) {
         for (int j = 0; j < M; j++) {
            cin >> ch;
            if (ch == 'H') houses.push_back(make_pair(i,j));
            else if (ch == 'm') people.push_back(make_pair(i,j));
         }
      }
      // compute the weights associated with matching house i to person j
      for (int i = 0; i < houses.size(); i++) {
         for (int j = 0; j < people.size(); j++) {
            cost[i][j] = abs(houses[i].first - people[j].first) + 
               abs(houses[i].second - people[j].second);
         }
      }
      cout << cycleCancel(houses.size()) << "\n";
   }
   return 0;
}

TC SRM387 D2 1000 (MarblesRegroupingHard)
URL: http://www.topcoder.com/stat?c=problem_statement&pm=8538

This is a re-visit of a previous problem we solved using DP. This is a slightly less obvious reduction than the previous problem; we need to compute the total number of colours which exists in the whole dataset. Then we can draw an edge from box i to colour j with a total weight of all the other j-coloured marbles from other boxes minus the ones in box i with colour j (as they don’t need to be moved). Take note that there may be a different number of colours to the number of boxes but we are always guaranteed that the number of colours is less than or equal to the number of boxes, hence a matching is always possible based on the pigeonhole principle. We can add dummy colours to the graph to reduce it to the assignment problem – which is what is done below:

class MarblesRegroupingHard {
public:
   int minMoves(vector <string>);
};

int totColours[15];
int B[51][15];

int cost[55][55];      // cost of assigning house i to person j
int curMatchings[55];   // house i matched to person curMatchings[i]
int delta[55][55];   // the delta/cost graph
int parent[55][55];   // used for backtracking the cycle
int cycle[55];         // path of the cycle
int totalLen;         // length of the cycle

bool augmentingCycle(int N) {
   memset(delta,0,sizeof(delta));
   memset(parent,0,sizeof(parent));
   memset(cycle,-1,sizeof(cycle));
   // derive the delta graph
   for (int i = 0; i < N; i++) {
      for (int j = 0; j < N; j++) {
         delta[i][j] = cost[i][curMatchings[j]] - cost[i][curMatchings[i]];
         parent[i][j] = j;
      }
   }
   for (int k = 0; k < N; k++) {
      for (int i = 0; i < N; i++) {
         for (int j = 0; j < N; j++) {
            if (delta[i][k] + delta[k][j] < delta[i][j]) {
               // minimum cost matching
               delta[i][j] = delta[i][k] + delta[k][j];
               // book-keep the optimal path so far for i->j
               parent[i][j] = parent[i][k];
               // detect for cycle (if vertex i and j are the same vertex)
               if (i == j) {
                  totalLen = 0;
                  // backtrack and construct the cycle path
                  do {
                     cycle[totalLen++] = i;
                     i = parent[i][j];
                  } while (i != j);
                  return true;
               }
            }
         }
      }
   }
   // did not find an augmenting path
   return false;
}

int cycleCancel(int N) {
   int res = 0;
   // pseudo max-flow for assignment problem
   for (int i = 0; i < N; i++) curMatchings[i] = i;
   while (augmentingCycle(N)) {
      // augment the negative cycle
      int r = curMatchings[cycle[0]];
      for (int i = 0; i < totalLen - 1; i++) {
         curMatchings[cycle[i]] = curMatchings[cycle[i+1]];
      }
      curMatchings[cycle[totalLen-1]] = r;
   }
   // compute the final cost
   for (int i = 0; i < N; i++) {
      res += cost[i][curMatchings[i]];
   }
   return res;
}


int MarblesRegroupingHard::minMoves(vector <string> boxes) {
   int res = 0;
   memset(totColours,0,sizeof(totColours));
   memset(B,0,sizeof(B));
   int numColours = 0;
   for (int i = 0; i < boxes.size(); i++) {
      istringstream iss(boxes[i]);
      int colours;
      int k = 0;
      while (iss >> colours) {
         totColours[k] += colours;
         B[i][k++] = colours;
      }
      numColours = k;
   }
   // make the weighted bipartite graph
   // use the cost matrix to fill in
   for (int i = 0; i < boxes.size(); i++) {
      for (int j = 0; j < boxes.size(); j++) {
         if (j >= numColours) cost[i][j] = 0;
         else cost[i][j] = totColours[j] - B[i][j];
      }   
   }
   return res = cycleCancel(boxes.size());
}

Friday, August 28, 2009

TC SRM447 D1 500 (PeopleYouMayKnow)

The problem statement can be accessed via:
http://www.topcoder.com/stat?c=problem_statement&pm=10580

The problem itself is pretty simple: it asks us to given a friends graph determine the number of friends to remove such that there are no paths from friend A to friend B which has a length less than 3.

The key observation is to remove all friends which have a shortest path between A and B which exceeds 3 as these will never conflict with our objective. To determine intermediate distances we can just use Floyd-Warshall's algorithm to calculate all-pairs shortest paths between friends. Then we can discard any nodes which violate dist[person1][i] + dist[i][person2] <= 3. Note that the problem fixes the constraint to not allow any paths with a length of less than 3 - this allows us to use a simple bipartite matching to determine the minimum cut of the network. For the general problem of having a length less than n, we can simply use a network flow algorithm (which we use for the implementation here).

We can divide the friends into two (one used as an in-node, the other used as an out-node) and run the maximum flow algorithm over the network. This concept is very similar to the one used in TC SRM360 D1 500 (Prince of Persia).

Implementation below:


class PeopleYouMayKnow {
public:
int maximalScore(vector <string>, int, int);
};

int augment(int cur, int flow, int ret);

int adjmat[201][201];
int used[201][201];
bool visited[201];
int t;

int maxflow() {
int flow = 0, thisFlow = 0;
memset(visited,false,sizeof(visited));
while ((thisFlow = augment(0, INT_MAX-1, 0)) != 0) {
memset(visited,false,sizeof(visited));
flow += thisFlow;
}
return flow;
}

int augment(int cur, int flow, int ret) {
if (cur == t-1) return flow;
visited[cur] = true;
for (int i = 0; i < t && ret == 0; i++) {
if (adjmat[cur][i] > used[cur][i] && !visited[i] &&
(ret = augment(i, min(flow,adjmat[cur][i]-used[cur][i]), 0)) > 0)
used[cur][i] += ret;
else if (used[i][cur] > 0 && !visited[i] &&
(ret = augment(i,min(flow, used[i][cur]), 0)) > 0)
used[i][cur] -= ret;
}
return ret;
}

int distTable[51][51];
int validNodes[51];
#define INF (1<<20)

int PeopleYouMayKnow::maximalScore(vector <string> friends, int A, int B) {
int res = 0;
int sz = friends.size();
for (int i = 0; i < sz; i++)
for (int j = 0; j < sz; j++)
if (i == j) distTable[i][j] = 0;
else distTable[i][j] = friends[i][j] == 'Y' ? 1 : INF;
for (int k = 0; k < sz; k++)
for (int i = 0; i < sz; i++)
for (int j = 0; j < sz; j++)
distTable[i][j] <?= distTable[i][k] + distTable[k][j];
// build the network - allowing infinite capacities between
// out -> in nodes and restricting capacities between nodes with 1 capacity
for (int i = 0; i < sz; i++)
for (int j = 0; j < sz; j++)
if (friends[i][j] == 'Y')
adjmat[i+1+sz][j+1] = INF;
for (int i = 0; i < sz; i++)
if (distTable[person1][i] + distTable[i][person2] > 3)
for (int j = 0; j < sz; j++) adjmat[j+1][i+1] = adjmat[i+1][j+1] = 0;
else
validNodes[i] = 1;
for (int i = 0; i < sz; i++)
if (validNodes[i]) adjmat[i+1][i+1+sz] = 1; else adjmat[i+1][i+1+sz] = 0;
adjmat[0][person1+1] = INF;
adjmat[person1+1][person1+1+sz] = INF;
adjmat[person2+1][person2+1+sz] = INF;
adjmat[person2+1+sz][2*sz+1] = INF;
t = 2*sz+2;
return res = maxflow();
}

Friday, August 7, 2009

TC SRM358 D1 1000 (SharksDinner)

The problem statement can be accessed via:
http://www.topcoder.com/stat?c=problem_statement&pm=7834&rd=10768

The problem states that there are sharks that can eat another provided that all their statistics (skill, speed and intelligence) are at least as good as another. We are asked to find the minimum number of sharks that will remain at the end, in other words, this can be viewed as finding the maximum number of sharks being eaten. This problem easily lends itself to a matching problem, if it weren't for the fact that a shark can eat at most 2 other sharks (instead of just one other shark) then this problem is simply a maximum cardinality bipartite matching problem. However, since the constraint is rather low (allowing a shark to eat at most 2 sharks) we can divide a shark into multiple nodes and solve it as a bipartite matching problem.

However we'll solve it using the maximum flow algorithm (which can also be used to solve the bipartite matching) by making a slight modification. For each shark we give it an "in" node and an "out" node and the capacity between these two nodes is exactly 2 (corresponding to the fact that a shark can only eat at most 2). Using this method allows us to scale the maximum number a shark can eat without increasing the number of nodes (beyond the in and out node expansion). One small note we need to be careful of is sharks which can eat each other. To handle this, a simple modification towards matching edges (i.e. shark A can eat shark B) is made - if shark A and shark B has the same statistics then shark A can eat shark B iff it has a lower index. So you resolve the problem of sharks eating each other cyclically.

Now it's just a matter of constructing the network and running the maximum flow algorithm over it. The final answer is simply the number of sharks in the data set minus the ones that got eaten from our matching.


class SharksDinner {
public:
int minSurvivors(vector <int>, vector <int>, vector <int>);
};

class Shark {
public:
int size;
int speed;
int intel;
Shark(int a, int b, int c) : size(a), speed(b), intel(c) { }
};

bool operator<(const Shark& lhs, const Shark& rhs) {
if (lhs.size != rhs.size) return lhs.size < rhs.size;
if (lhs.speed != rhs.speed) return lhs.speed < rhs.speed;
return lhs.intel < rhs.intel;
}

int augment(int cur, int flow, int ret);

class node {
public:
int num;
int cap;
node() { }
node(int n, int c) : num(n), cap(c) { }
};

bool operator<(const node& lhs, const node& rhs) {
if (lhs.num != rhs.num) return lhs.num < rhs.num;
if (lhs.cap != rhs.cap) return lhs.cap < rhs.cap;
return false;
}

#define INF 9999999
int sz; // the size of the graph (vertices)
int cap[202][202]; // capacity matrix
int used[202][202]; // network flow array
bool visited[202];
vector<vector<node> > adjlist;

int maxflow() {
int flow = 0;
int thisFlow = 0;
memset(visited,false,sizeof(visited));
while ((thisFlow = augment(0, INT_MAX-1, 0)) != 0) {
memset(visited,false,sizeof(visited));
flow += thisFlow;
}
return flow;
}

int augment(int cur, int flow, int ret) {
if (cur == sz-1) return flow;
visited[cur] = true;
for (int k = 0; k < adjlist[cur].size() && ret == 0; k++) {
int i = adjlist[cur][k].num;
int c = adjlist[cur][k].cap;
if (c > used[cur][i] && !visited[i] && (ret = augment(i, min(flow,
c-used[cur][i]), 0)) > 0) used[cur][i] += ret;
else if (used[i][cur] > 0 && !visited[i] && (ret = augment(i,
min(flow, used[i][cur]), 0)) > 0)
used[i][cur] -= ret;
}
return ret;
}

int SharksDinner::minSurvivors(vector <int> size, vector <int> speed,
vector <int> intelligence) {
int res = 0;
vector<Shark> sharkData;
for (int i = 0; i < size.size(); i++) {
sharkData.push_back(Shark(size[i],speed[i],intelligence[i]));
}
adjlist = vector<vector<node> >(sharkData.size() * 3 + 2, vector<node>());
// source
for (int i = 0; i < sharkData.size(); i++) {
adjlist[0].push_back(node(i+1,INF));
}
// sub-nodes
for (int i = 0; i < sharkData.size(); i++) {
adjlist[i+1].push_back(node(i+1+sharkData.size(),2));
}
// end-nodes
for (int i = 0; i < sharkData.size(); i++) {
for (int j = 0; j < sharkData.size(); j++) {
if (i == j) continue;
if (sharkData[i].size < sharkData[j].size ||
sharkData[i].speed < sharkData[j].speed ||
sharkData[i].intel < sharkData[j].intel) continue;
if (sharkData[i].size == sharkData[j].size &&
sharkData[i].speed == sharkData[j].speed &&
sharkData[i].intel == sharkData[j].intel && i > j) continue;
adjlist[i+1+sharkData.size()].push_back(node(j+1+(sharkData.size()*2),INF));
adjlist[j+1+(sharkData.size()*2)].push_back(node(i+1+sharkData.size(),0));
}
}
// sink
for (int i = 0; i < sharkData.size(); i++) {
adjlist[i+1+(sharkData.size()*2)].push_back(node(3*sharkData.size()+1,1));
}
sz = 3 * sharkData.size() + 2;
// run flow
int flow = maxflow();
return res = sharkData.size() - flow;
}

Friday, July 31, 2009

TC SRM360 D1 500 (PrinceOfPersia)

The problem statement can be accessed via:
http://www.topcoder.com/stat?c=problem_statement&pm=7876&rd=10772

This problem lends itself to a graph problem as can be observed by the fact that, the input format is in a grid format and there are obstacles (forbidden cells) and its related to path accessibility between the prince and princess. One way to solve this problem is to use the Maximum Flow algorithm. Take note, although they are required to "meet" in any cell - it is sufficient enough for the prince (or princess) to go to princess (or prince) as it is not limited by time.

How does it relate to the Maximum Flow problem? We are basically asked to find the minimum cut and using the max-flow min-cut theorem we can calculate this by running a maximum flow algorithm through the graph. To see that what we really want is the minimum cut, let's revise the definition of a cut (in terms of graph theory of course): A cut is a partition of vertices of a graph into two disjoint subsets. So if we imagine the princess in one subset and the prince in another subset, the minimum cut between these two disjoint subsets is the number of accessible ways in which one can go to another.

So now we move on to applying the actual algorithm to the problem. If we let each empty cell be a node (or vertex) on the network and assign 1-capacity edges to adjacent empty nodes and apply flow from one "P" (source) to another "P" (sink) do we obtain the correct answer? Well, unfortunately no. The problem lies in which one "bottleneck" node (on the minimum cut set) can have an outgoing flow of more than 1, so it may over-calculate the answer. To see this, consider the example:

# . P
. . .
P . #

The central node in the snippet can have an incoming flow of 2 if both its adjacent nodes gives it a flow of 1. The problem with this is that the bottom-left P (assume it to be the sink) will have a maximum flow of 2 (and hence a minimum cut of 2). So our problem now is to ensure that each empty node will only give out a maximum of 1 out-going flow.

One wrong approach is to ensure that there is only an incoming flow of maximum 1 (i.e. 1 edge going into the node). The problem with this is that it totally breaks down the network layout. The solution is to further break each empty cell into 2 nodes: let's call them A and B. These nodes are used as "filters" or more precisely input and output filters. The edge from A->B will have a capacity of 1 which limits the maximum flow across the node to be 1 whereas the incoming and outgoing arcs of A and B respectively are infinite. This allows the network to retain the correct layout (and hence the correct answer) as well as fixing the original problem. We just need to ensure that all input arcs to the node use the sub-node A, and all output arcs from the node use the sub-node B.

Implementation-wise it's pretty standard if you follow the ideas discussed earlier. A sample implementation from the practice room is given below:


// the maximum number of vertices
#define NN 256

// adjacency matrix (fill this up)
int cap[NN][NN];

// flow network
int fnet[NN][NN];

// BFS
int q[NN], qf, qb, prev[NN];

int fordFulkerson( int n, int s, int t )
{
memset( fnet, 0, sizeof( fnet ) );
int flow = 0;
// find an augmenting path of at least 1
while( true )
{
memset( prev, -1, sizeof( prev ) );
qf = qb = 0;
prev[q[qb++] = s] = -2;
while( qb > qf && prev[t] == -1 )
for( int u = q[qf++], v = 0; v < n; v++ )
if( prev[v] == -1 && fnet[u][v] - fnet[v][u] < cap[u][v] )
prev[q[qb++] = v] = u;
if( prev[t] == -1 ) break;
// get the bottleneck capacity
int bot = 0x7FFFFFFF;
for( int v = t, u = prev[v]; u >= 0; v = u, u = prev[v] )
bot <?= cap[u][v] - fnet[u][v] + fnet[v][u];
// update the flow network
for( int v = t, u = prev[v]; u >= 0; v = u, u = prev[v] )
fnet[u][v] += bot;
flow += bot;
}
return flow;
}

int dx[] = {-1, 1, 0, 0};
int dy[] = {0, 0, -1, 1};

#define INF 999999

int PrinceOfPersia::minObstacles(vector <string> maze) {
int res = 0;
bool gotP = false;
int source = 0, sink = 0;
for (int i = 0; i < maze.size(); i++) {
for (int j = 0; j < maze[i].size(); j++) {
// assign the source and sink
if (maze[i][j] == 'P' && !gotP) {
source = (i * maze[0].size() + j)*2; gotP = true;
}
else if (maze[i][j] == 'P' && gotP) {
sink = (i * maze[0].size() + j)*2;
}
// build connections
if (maze[i][j] == '.') {
cap[(i * maze[0].size() + j)*2][(i * maze[0].size() + j)*2+1] = 1;
}
else if (maze[i][j] == 'P') {
cap[(i * maze[0].size() + j)*2][(i * maze[0].size() + j)*2+1] = INF;
}
// iterate through adjacent cells
for (int k = 0; k < 4; k++) {
int mi = i + dx[k];
int mj = j + dy[k];
if (mi < 0 || mj < 0 || mi >= maze.size() || mj >= maze[0].size())
continue;
// impossible case
if (maze[mi][mj] == 'P' && maze[i][j] == 'P') return -1;
if (maze[mi][mj] != '#') {
cap[(i * maze[0].size() + j)*2+1][(mi * maze[0].size() + mj)*2] = INF;
}
}
}
}
int flow = fordFulkerson(maze[0].size()*maze.size()*2, source, sink);
return res = flow;
}