// Flujo máximo: la clase que se usa en todo el capítulo. Las aristas se guardan en pares: la arista
// e y su inversa e ^ 1, que nace con capacidad 0. Guardar la capacidad que queda es guardar la red
// residual.
#pragma once

#include <algorithm>
#include <cassert>
#include <cstddef>
#include <limits>
#include <queue>
#include <vector>

class MaxFlow {
public:
    explicit MaxFlow(int nodes) : adjacency_(static_cast<std::size_t>(nodes)) {}

    // Devuelve el número de la arista, para preguntar después cuánto flujo lleva.
    int add_edge(int from, int to, long long capacity) {
        assert(capacity >= 0);
        int id = static_cast<int>(arcs_.size());
        arcs_.push_back({to, capacity});
        arcs_.push_back({from, 0});
        adjacency_[index(from)].push_back(id);
        adjacency_[index(to)].push_back(id + 1);
        return id;
    }

    // Camino de aumento con menos aristas, hallado con BFS. Se detiene apenas llega a t.
    long long edmonds_karp(int s, int t) {
        long long total = 0;
        while (true) {
            std::vector<int> arrived_by(adjacency_.size(), -1);
            std::vector<bool> seen(adjacency_.size(), false);
            std::queue<int> pending;
            seen[index(s)] = true;
            pending.push(s);
            while (!pending.empty() && !seen[index(t)]) {
                int v = pending.front();
                pending.pop();
                for (int id : adjacency_[index(v)]) {
                    const Arc& arc = arcs_[index(id)];
                    if (arc.residual > 0 && !seen[index(arc.to)]) {
                        seen[index(arc.to)] = true;
                        arrived_by[index(arc.to)] = id;
                        pending.push(arc.to);
                    }
                }
            }
            if (!seen[index(t)]) return total;
            total += push_along(arrived_by, s, t);
        }
    }

    // Camino de aumento cualquiera, el primero que encuentra un DFS.
    long long ford_fulkerson(int s, int t) {
        long long total = 0;
        while (true) {
            std::vector<int> arrived_by(adjacency_.size(), -1);
            std::vector<bool> seen(adjacency_.size(), false);
            std::vector<int> pending = {s};
            seen[index(s)] = true;
            while (!pending.empty() && !seen[index(t)]) {
                int v = pending.back();
                pending.pop_back();
                for (int id : adjacency_[index(v)]) {
                    const Arc& arc = arcs_[index(id)];
                    if (arc.residual > 0 && !seen[index(arc.to)]) {
                        seen[index(arc.to)] = true;
                        arrived_by[index(arc.to)] = id;
                        pending.push_back(arc.to);
                    }
                }
            }
            if (!seen[index(t)]) return total;
            total += push_along(arrived_by, s, t);
        }
    }

    // Lo alcanzable desde s en la red residual. Al terminar el flujo, es el lado de s de un corte mínimo.
    std::vector<bool> reachable_from(int s) const {
        std::vector<bool> seen(adjacency_.size(), false);
        std::vector<int> pending = {s};
        seen[index(s)] = true;
        while (!pending.empty()) {
            int v = pending.back();
            pending.pop_back();
            for (int id : adjacency_[index(v)]) {
                const Arc& arc = arcs_[index(id)];
                if (arc.residual > 0 && !seen[index(arc.to)]) {
                    seen[index(arc.to)] = true;
                    pending.push_back(arc.to);
                }
            }
        }
        return seen;
    }

    // El flujo de una arista es lo que su inversa puede devolver.
    long long flow_on(int id) const { return arcs_[index(id ^ 1)].residual; }

    // Cuánta capacidad le queda a una arista.
    long long left_on(int id) const { return arcs_[index(id)].residual; }

    // Quita una arista de la red, junto con lo que deshacerla permitiría.
    void remove_edge(int id) {
        arcs_[index(id)].residual = 0;
        arcs_[index(id ^ 1)].residual = 0;
    }

    int augmentations() const { return augmentations_; }

    // Los caminos de aumento usados, si se pidió guardarlos.
    void keep_paths() { keep_paths_ = true; }
    const std::vector<std::vector<int>>& paths() const { return paths_; }

private:
    struct Arc {
        int to;
        long long residual;
    };

    static std::size_t index(int i) { return static_cast<std::size_t>(i); }

    long long push_along(const std::vector<int>& arrived_by, int s, int t) {
        long long amount = std::numeric_limits<long long>::max();
        std::vector<int> path = {t};
        for (int v = t; v != s; v = arcs_[index(arrived_by[index(v)] ^ 1)].to) {
            amount = std::min(amount, arcs_[index(arrived_by[index(v)])].residual);
            path.push_back(arcs_[index(arrived_by[index(v)] ^ 1)].to);
        }
        for (int v = t; v != s; v = arcs_[index(arrived_by[index(v)] ^ 1)].to) {
            int id = arrived_by[index(v)];
            arcs_[index(id)].residual -= amount;
            arcs_[index(id ^ 1)].residual += amount;
        }
        augmentations_++;
        if (keep_paths_) {
            std::reverse(path.begin(), path.end());
            paths_.push_back(path);
        }
        return amount;
    }

    std::vector<std::vector<int>> adjacency_;
    std::vector<Arc> arcs_;
    int augmentations_ = 0;
    bool keep_paths_ = false;
    std::vector<std::vector<int>> paths_;
};