Dinic’s algorithm for Maximum Flow with non-integer edge capacities

Question

There is a C++ implementation of Dinic’s algorithm for Maximum Flow Problem that I was trying to use. In that C++ code, it is assumed that all parameters are integer. I tryied to convert the code to an equivalent one in which capacity of edges can be continuous (non-integer). Here is my code:

// C++ implementation of Dinic's Algorithm 
// #include "HEADERS.h"

#include<bits/stdc++.h> 



using namespace std; 

double EPSILON = 1e-6 ;
// A structure to represent a edge between 
// two vertex 
struct Edge 
{ 
    int v ; // Vertex v (or "to" vertex) 
            // of a directed edge u-v. "From" 
            // vertex u can be obtained using 
            // index in adjacent array. 

    double flow ; // flow of data in edge 

    double C; // capacity 

    int rev ; // To store index of reverse 
            // edge in adjacency list so that 
            // we can quickly find it. 
}; 

// Residual Graph 
class Graph 
{ 
    int V; // number of vertex 
    int *level ; // stores level of a node 
    vector< Edge > *adj; 
public : 
    Graph(int V) 
    { 
        adj = new vector<Edge>[V]; 
        this->V = V; 
        level = new int[V]; 
    } 

    // add edge to the graph 
    void addEdge(int u, int v, double C) 
    { 
        // Forward edge : 0 flow and C capacity 
        Edge a{v, 0, C, static_cast<int> (adj[v].size()) }; 

        // Back edge : 0 flow and 0 capacity 
        Edge b{u, 0, 0, static_cast<int> (adj[u].size()) }; 

        adj[u].push_back(a); 
        adj[v].push_back(b); // reverse edge 
    } 

    bool BFS(int s, int t); 
    int sendFlow(int s, double flow, int t, int ptr[]); 
    double DinicMaxflow(int s, int t); 
}; 

// Finds if more flow can be sent from s to t. 
// Also assigns levels to nodes. 
bool Graph::BFS(int s, int t) 
{ 
    for (int i = 0 ; i < V ; i++) 
        level[i] = -1; 

    level[s] = 0; // Level of source vertex 

    // Create a queue, enqueue source vertex 
    // and mark source vertex as visited here 
    // level[] array works as visited array also. 
    list< int > q; 
    q.push_back(s); 

    vector<Edge>::iterator i ; 
    while (!q.empty()) 
    { 
        int u = q.front(); 
        q.pop_front(); 
        for (i = adj[u].begin(); i != adj[u].end(); i++) 
        { 
            Edge &e = *i; 
            if (level[e.v] < 0 && e.flow < e.C) 
            { 
                // Level of current vertex is, 
                // level of parent + 1 
                level[e.v] = level[u] + 1; 

                q.push_back(e.v); 
            } 
        } 
    } 

    // IF we can not reach to the sink we 
    // return false else true 
    return level[t] < 0 ? false : true ; 
} 

// A DFS based function to send flow after BFS has 
// figured out that there is a possible flow and 
// constructed levels. This function called multiple 
// times for a single call of BFS. 
// flow : Current flow send by parent function call 
// start[] : To keep track of next edge to be explored. 
//       start[i] stores count of edges explored 
//       from i. 
// u : Current vertex 
// t : Sink 
int Graph::sendFlow(int u, double flow, int t, int start[]) 
{ 
    // Sink reached 
    if (u == t) 
        return flow; 

    // Traverse all adjacent edges one -by - one. 
    for ( ; start[u] < static_cast <int> (adj[u].size()) ; start[u]++) 
    { 
        // Pick next edge from adjacency list of u 
        Edge &e = adj[u][start[u]]; 

        if (level[e.v] == level[u]+1 && e.flow < e.C) 
        { 
            // find minimum flow from u to t 
            double curr_flow = min(flow, e.C - e.flow); 

            double temp_flow = sendFlow(e.v, curr_flow, t, start); 

            // flow is greater than zero 
            if (temp_flow > 0) 
            { 
                // add flow to current edge 
                e.flow += temp_flow; 

                // subtract flow from reverse edge 
                // of current edge 
                adj[e.v][e.rev].flow -= temp_flow; 
                return temp_flow; 
            } 
        } 
    } 

    return 0; 
} 
// Returns maximum flow in graph 
double Graph::DinicMaxflow(int s, int t) 
{ 
    // Corner case 
    if (s == t) 
        return -1; 

    double total = 0; // Initialize result 

    // Augment the flow while there is path 
    // from source to sink 
    while (BFS(s, t) == true) 
    { 
        // store how many edges are visited 
        // from V { 0 to V } 
        int *start = new int[V+1]; 
        double flow = sendFlow(s, INT_MAX, t, start) ;
        // while flow is not zero in graph from S to D 
        while (flow > EPSILON )
        {
            // Add path flow to overall flow 
            total += flow;
            flow = sendFlow(s, INT_MAX, t, start) ;
        }
    } 

    // return maximum flow 
    return total; 
} 

// Driver program to test above functions 
int main() 
{ 
    Graph g(6); 
    g.addEdge(0, 1, 16 ); 
    g.addEdge(0, 2, 13 ); 
    g.addEdge(1, 2, 10 ); 
    g.addEdge(1, 3, 12 ); 
    g.addEdge(2, 1, 4 ); 
    g.addEdge(2, 4, 14); 
    g.addEdge(3, 2, 9 ); 
    g.addEdge(3, 5, 20 ); 
    g.addEdge(4, 3, 7 ); 
    g.addEdge(4, 5, 4); 

    // next exmp 
    /*g.addEdge(0, 1, 3 ); 
    g.addEdge(0, 2, 7 ) ; 
    g.addEdge(1, 3, 9); 
    g.addEdge(1, 4, 9 ); 
    g.addEdge(2, 1, 9 ); 
    g.addEdge(2, 4, 9); 
    g.addEdge(2, 5, 4); 
    g.addEdge(3, 5, 3); 
    g.addEdge(4, 5, 7 ); 
    g.addEdge(0, 4, 10); 

    // next exp 
    g.addEdge(0, 1, 10); 
    g.addEdge(0, 2, 10); 
    g.addEdge(1, 3, 4 ); 
    g.addEdge(1, 4, 8 ); 
    g.addEdge(1, 2, 2 ); 
    g.addEdge(2, 4, 9 ); 
    g.addEdge(3, 5, 10 ); 
    g.addEdge(4, 3, 6 ); 
    g.addEdge(4, 5, 10 ); */
    cout << "Maximum flow " << g.DinicMaxflow(0, 5) << endl ; 
    return 0; 
} 

In this code, I have changed type of flow and C to double. I also changed some parts of the code to adapt it to this new double type. The code works only occusionally! It either outputs Maximum flow 23 which is the right output or throws a Segmentation fault (core dumped) error. I don't really know what is wrong in this code. Any ideas?

Answer 1

I don't know if the algorithm is correct, but assuming it is, the code at the link has several issues.

1. Usage of the #include<bits/stdc++.h> header.
2. Memory leaks.

First, usage of bits/stdc++.h should be avoided, and the proper #include files should be used.

Second, the usage of std::vector would take care of the memory leaks. The code uses std::vector in some places, but totally refuses to use it in other places. For example:

int *start = new int[V+1];

should be replaced with:

std::vector<int> start(V+1);

The former was causing a memory leak due to the lack of a delete [] start; in the code. Using std::vector, the memory leak disappears.

Once std::vector is introduced, there is no need for the V member variable in the Graph class that tracks the number of vertices. The reason why is that the vector members are sized to V vertices, and a vector already knows their own size by using the vector::size() member function. So the V member variable is redundant and can be removed.

The last change that can be made is to move struct Edge within the Graph class.

---

Given all of the changes mentioned, here is a cleaned up version that returns the same result as the original code when run with the test graph set up in the main() function:

#include <vector>
#include <list>
#include <algorithm>
#include <iostream>
#include <climits>

class Graph
{
    struct Edge
    {
        int v; 
        int flow; 
        int C; 
        int rev; 
    };

    std::vector<int> level;
    std::vector<std::vector< Edge >> adj;

public:
    Graph(int V) : level(V), adj(V) {}
    void addEdge(int u, int v, int C)
    {
        adj[u].push_back({ v, 0, C, static_cast<int>(adj[v].size()) });
        adj[v].push_back({ u, 0, 0, static_cast<int>(adj[u].size()) }); 
    }

    bool BFS(int s, int t);
    int sendFlow(int s, int flow, int t, std::vector<int>& ptr);
    int DinicMaxflow(int s, int t);
};

bool Graph::BFS(int s, int t)
{
    std::fill(level.begin(), level.end(), -1);
    level[s] = 0; 
    std::list< int > q;
    q.push_back(s);
    std::vector<Edge>::iterator i;
    while (!q.empty())
    {
        int u = q.front();
        q.pop_front();
        for (i = adj[u].begin(); i != adj[u].end(); i++)
        {
            Edge &e = *i;
            if (level[e.v] < 0 && e.flow < e.C)
            {
                level[e.v] = level[u] + 1;
                q.push_back(e.v);
            }
        }
    }
    return level[t] < 0 ? false : true;
}

int Graph::sendFlow(int u, int flow, int t, std::vector<int>& start)
{
    if (u == t)
        return flow;
    for (; start[u] < static_cast<int>(adj[u].size()); start[u]++)
    {
        // Pick next edge from adjacency list of u 
        Edge &e = adj[u][start[u]];

        if (level[e.v] == level[u] + 1 && e.flow < e.C)
        {
            int curr_flow = std::min(flow, e.C - e.flow);
            int temp_flow = sendFlow(e.v, curr_flow, t, start);

            if (temp_flow > 0)
            {
                e.flow += temp_flow;
                adj[e.v][e.rev].flow -= temp_flow;
                return temp_flow;
            }
        }
    }
    return 0;
}

int Graph::DinicMaxflow(int s, int t)
{
    if (s == t)
        return -1;
    int total = 0; 
    while (BFS(s, t) == true)
    {
        std::vector<int> start(level.size() + 1);
        while (int flow = sendFlow(s, INT_MAX, t, start))
            total += flow;
    }
    return total;
}

The test function is as follows:

int main() 
{ 
    Graph g(6); 
    g.addEdge(0, 1, 16 ); 
    g.addEdge(0, 2, 13 ); 
    g.addEdge(1, 2, 10 ); 
    g.addEdge(1, 3, 12 ); 
    g.addEdge(2, 1, 4 ); 
    g.addEdge(2, 4, 14); 
    g.addEdge(3, 2, 9 ); 
    g.addEdge(3, 5, 20 ); 
    g.addEdge(4, 3, 7 ); 
    g.addEdge(4, 5, 4); 
    std::cout << "Maximum flow " << g.DinicMaxflow(0, 5); 
}

Live Example

---

Now, if you want to see what happens if you use double instead of int as the flow amount, the easiest thing to do is to create a template class based on the code above. The goal is to take where int is used and instead of replacing int with double, replace int with a template parameter (for example, T). The resulting code would be as follows:

#include <vector>
#include <list>
#include <algorithm>
#include <iostream>
#include <climits>

template <typename T>
class Graph
{
    struct Edge
    {
        int v;
        T flow;
        T C;
        int rev;
    };

    std::vector<int> level;
    std::vector<std::vector<Edge>> adj;

public:
    Graph(int V) : level(V), adj(V) {}
    void addEdge(int u, int v, T C)
    {
        adj[u].push_back({ v, T(), C, static_cast<int>(adj[v].size())});
        adj[v].push_back({ u, T(), T(), static_cast<int>(adj[u].size())}); // reverse edge 
    }
    bool BFS(int s, int t);
    T sendFlow(int s, T flow, int t, std::vector<int>& ptr);
    T DinicMaxflow(int s, int t);
};

template <typename T>
bool Graph<T>::BFS(int s, int t)
{
    std::fill(level.begin(), level.end(), -1);
    level[s] = 0; 
    std::list< int > q;
    q.push_back(s);

    typename std::vector<Edge>::iterator i;
    while (!q.empty())
    {
        int u = q.front();
        q.pop_front();
        for (i = adj[u].begin(); i != adj[u].end(); i++)
        {
            Edge &e = *i;
            if (level[e.v] < 0 && e.flow < e.C)
            {
                level[e.v] = level[u] + 1;
                q.push_back(e.v);
            }
        }
    }
    return level[t] < 0 ? false : true;
}

template <typename T>
T Graph<T>::sendFlow(int u, T flow, int t, std::vector<int>& start)
{
    if (u == t)
        return flow;

    for (; start[u] < static_cast<int>(adj[u].size()); start[u]++)
    {
        Edge &e = adj[u][start[u]];
        if (level[e.v] == level[u] + 1 && e.flow < e.C)
        {
            T curr_flow = std::min(flow, e.C - e.flow);
            T temp_flow = sendFlow(e.v, curr_flow, t, start);
            if (temp_flow > 0)
            {
                e.flow += temp_flow;
                adj[e.v][e.rev].flow -= temp_flow;
                return temp_flow;
            }
        }
    }
    return 0;
}

template <typename T>
T Graph<T>::DinicMaxflow(int s, int t)
{
    if (s == t)
        return -1;
    T total = 0; 

    while (BFS(s, t) == true)
    {
        std::vector<int> start(level.size() + 1);
        while (T flow = sendFlow(s, INT_MAX, t, start))
            total += flow;
    }
    return total;
}

Live Example

The main test example would be to simply use Graph<double> or Graph<int> instead of simply Graph -- everything else in the main function stays the same.

Now that the function is a template, any numerical type that supports the same operations as int or double can be substituted easily by just creating a Graph<numerical_type>.