Corei13's blog

By Corei13, 9 years ago, In English

The main challenges with compile time computations are

  • we can only work with constants (so no looping is possible)
  • maximum recursion depth is very limited compared to runtime computations (for example, maximum recursion depth is only 900 in g++)

Here is an example of how to use template meta programming for compile time computations. The example shows how to implement primality checking in compile time, but the same principle can be applied to other problems. For instance, in the main function, replacing Generate <level, Prime>::value by Generate <level, Sqrt>::value will generate integer square roots for all numbers between [0, 2^level).

#include <bits/stdc++.h>

using namespace std;

namespace Helpers {
    template <int lhs, int rhs>
    struct Plus {
        static const int value = lhs + rhs;
    };
    struct True {
        static const bool value = true;
    };
    struct False {
        static const bool value = false;
    };
    // Array <T, T args...> contains the array {args...} of size ::size at ::value
    // For example if A = Array <int, 1, 2, 3>, then A::value[] = {1, 2, 3} and A::size = 3
    template <typename T, T... args>
    struct Array {
        static const T value[sizeof...(args)];
        static const int size = sizeof...(args);
    };
    template <typename T, T... args>
    const T Array <T, args...>::value[sizeof...(args)] = {args...};
}

// Sqrt <N>::value is the integer square root of N
template <int N>
struct Sqrt {
    template <int lo, int hi>
    struct Helper {
        static const int mid = (lo + hi + 1) / 2;
        static const int value = conditional < (N / mid < mid), Helper < lo, mid - 1 >, Helper <mid, hi> >::type::value;
    };
    template <int n>
    struct Helper <n, n> {
        static const int value = n;
    };
    static const int value = Helper <0, N>::value;
};


// Prime <N>::value is true if and only if N is a prime
template <int N>
struct Prime {
    // Helper::value is true iff any prime in [lo, hi] divides N
    template <int lo, int hi>
    struct Helper {
        static const int mid = (lo + hi) / 2;
        static const bool value = conditional < Helper <lo, mid>::value, Helpers::True, Helper < mid + 1, hi >>::type::value;
    };
    template <int d>
    struct Helper <d, d> {
        static const bool value = Prime <d>::value and N % d == 0;
    };
    static const bool value = not conditional < (N <= 1), Helpers::True, Helper <1, Sqrt <N>::value> >::type::value;
};

/*
    Generate <N, F> contains the Helpers::Array with ::value = {F<0>::value, ..., F<2 ^ N - 1>::value} at ::value
    Example:
        using A = Generate<5, Sqrt>::value;
        for (int i = 0; i < A::size; ++i) {
            cout << A::value[i] << endl;
        }
*/

template <int N, template <int> class F>
struct Generate {
    template <int n, int... args>
    struct Helper {
        using value = typename Helper < n - 1, args..., Helpers::Plus <args, sizeof...(args)>::value... >::value;
    };

    template <int... args>
    struct Helper <0, args...> {
        using value = Helpers::Array <decltype(F <0>::value), F <args>::value...>;
    };
    using value = typename Helper <N, 0>::value;
};


int main(int argc, char const *argv[]) {
    const int level = 10; // for primes upto 2^level
    using A = Generate <level, Prime>::value;
    cout << A::size << endl;
    for (int i = 0; i < A::size; ++i) {
        cout << "> " << boolalpha << i << ' ' << A::value[i] << endl;
    }
    return 0;
}

Full text and comments »

  • Vote: I like it
  • +60
  • Vote: I do not like it

By Corei13, 9 years ago, In English

For those who are interested, a C++11 implementation of highest-label push relabel maximum flow algorithm. Style and format is taken from here.

Uses just gap relabeling heuristic. Global relabeling heuristic was implemented but removed because of poor performance.

/*
    Implementation of highest-label push-relabel maximum flow
    with gap relabeling heuristic.

    Running time:
        O(|V|^2|E|^{1/2})

    Usage:
        - add edges by AddEdge()
        - GetMaxFlow(s, t) returns the maximum flow from s to t
    
    Input:
        - graph, constructed using AddEdge()
        - (s, t), (source, sink)

    Output:
        - maximum flow value

    Todo:
        - implement Phase II (flow network from preflow network)
        - implement GetMinCut()
*/

template <class T> struct Edge {
    int from, to, index;
    T cap, flow;

    Edge(int from, int to, T cap, T flow, int index): from(from), to(to), cap(cap), flow(flow), index(index) {}
};

template <class T> struct PushRelabel {
    int n;
    vector <vector <Edge <T>>> adj;
    vector <T> excess;
    vector <int> dist, count;
    vector <bool> active;
    vector <vector <int>> B;
    int b;
    queue <int> Q;

    PushRelabel (int n): n(n), adj(n) {}

    void AddEdge (int from, int to, int cap) {
        adj[from].push_back(Edge <T>(from, to, cap, 0, adj[to].size()));
        if (from == to) {
            adj[from].back().index++;
        }
        adj[to].push_back(Edge <T>(to, from, 0, 0, adj[from].size() - 1));

    }

    void Enqueue (int v) {
        if (!active[v] && excess[v] > 0 && dist[v] < n) {
            active[v] = true;
            B[dist[v]].push_back(v);
            b = max(b, dist[v]);
        }
    }

    void Push (Edge <T> &e) {
        T amt = min(excess[e.from], e.cap - e.flow);
        if (dist[e.from] == dist[e.to] + 1 && amt > T(0)) {
            e.flow += amt;
            adj[e.to][e.index].flow -= amt;
            excess[e.to] += amt;    
            excess[e.from] -= amt;
            Enqueue(e.to);
        }
    }

    void Gap (int k) {
        for (int v = 0; v < n; v++) if (dist[v] >= k) {
            count[dist[v]]--;
            dist[v] = max(dist[v], n);
            count[dist[v]]++;
            Enqueue(v);
        }
    }

    void Relabel (int v) {
        count[dist[v]]--;
        dist[v] = n;
        for (auto e: adj[v]) if (e.cap - e.flow > 0) {
            dist[v] = min(dist[v], dist[e.to] + 1);
        }
        count[dist[v]]++;
        Enqueue(v);
    }

    void Discharge(int v) {
        for (auto &e: adj[v]) {
            if (excess[v] > 0) {
                Push(e);
            } else {
                break;
            }
        }

        if (excess[v] > 0) {
            if (count[dist[v]] == 1) {
                Gap(dist[v]); 
            } else {
                Relabel(v);
            }
        }
    }

    T GetMaxFlow (int s, int t) {
        dist = vector <int>(n, 0), excess = vector<T>(n, 0), count = vector <int>(n + 1, 0), active = vector <bool>(n, false), B = vector <vector <int>>(n), b = 0;
        
        for (auto &e: adj[s]) {
            excess[s] += e.cap;
        }

        count[0] = n;
        Enqueue(s);
        active[t] = true;
        
        while (b >= 0) {
            if (!B[b].empty()) {
                int v = B[b].back();
                B[b].pop_back();
                active[v] = false;
                Discharge(v);
            } else {
                b--;
            }
        }
        return excess[t];
    }

    T GetMinCut (int s, int t, vector <int> &cut);
};

Full text and comments »

  • Vote: I like it
  • +27
  • Vote: I do not like it