HeavyPathDecomposition.h

Go to the documentation of this file.

 
#pragma once
 
#include <ogdf/basic/Math.h>
#include <ogdf/graphalg/steiner_tree/EdgeWeightedGraphCopy.h>
 
#include <vector>
 
namespace ogdf {
namespace steiner_tree {
 
template<typename T>
class HeavyPathDecomposition {
private:
    const EdgeWeightedGraphCopy<T>& tree; 
    List<node> terminals; 
    const NodeArray<bool> isTerminal; 
 
    std::vector<std::vector<node>> chains; 
    std::vector<std::vector<node>> chainsOfTerminals; 
    NodeArray<int> chainOfNode; 
    NodeArray<int> positionOnChain; 
    NodeArray<int> weightOfSubtree; 
    NodeArray<int> nodeLevel; 
    NodeArray<T> distanceToRoot; 
    NodeArray<node> closestSteinerAncestor; 
    std::vector<node> fatherOfChain; 
 
    std::vector<std::vector<T>> longestDistToSteinerAncestorOnChain; 
    std::vector<std::vector<T>> longestDistToSteinerAncestorSegTree; 
 
    void buildMaxSegmentTree(std::vector<T>& segmentTree, const int nodeIndex, const int left,
            const int right, const std::vector<T>& baseArray) const {
        if (left == right) { // leaf
            segmentTree[nodeIndex] = baseArray[left];
            return;
        }
 
        int middle = (left + right) >> 1;
        int leftNodeIndex = nodeIndex + nodeIndex + 1;
        int rightNodeIndex = leftNodeIndex + 1;
 
        buildMaxSegmentTree(segmentTree, leftNodeIndex, left, middle, baseArray);
        buildMaxSegmentTree(segmentTree, rightNodeIndex, middle + 1, right, baseArray);
 
        segmentTree[nodeIndex] = max(segmentTree[leftNodeIndex], segmentTree[rightNodeIndex]);
    }
 
    T getMaxSegmentTree(const std::vector<T>& segmentTree, const int nodeIndex, const int left,
            const int right, const int queryLeft, const int queryRight) const {
        if (queryLeft > queryRight || left > queryRight || queryLeft > right) {
            return 0;
        }
 
        if (queryLeft <= left && right <= queryRight) { // included
            return segmentTree[nodeIndex];
        }
 
        int middleIndex = (left + right) >> 1;
        int leftNodeIndex = nodeIndex + nodeIndex + 1;
        int rightNodeIndex = leftNodeIndex + 1;
 
        T maxValue(0);
        if (queryLeft <= middleIndex) {
            Math::updateMax(maxValue,
                    getMaxSegmentTree(segmentTree, leftNodeIndex, left, middleIndex, queryLeft,
                            queryRight));
        }
        if (queryRight > middleIndex) {
            Math::updateMax(maxValue,
                    getMaxSegmentTree(segmentTree, rightNodeIndex, middleIndex + 1, right,
                            queryLeft, queryRight));
        }
 
        return maxValue;
    }
 
    T distanceToAncestor(node v, node ancestor) const {
        if (ancestor == nullptr) {
            return distanceToRoot[v];
        }
        return distanceToRoot[v] - distanceToRoot[ancestor];
    }
 
    void computeLongestDistToSteinerAncestorOnChain() {
        longestDistToSteinerAncestorOnChain.resize((int)chains.size());
        for (int i = 0; i < (int)chains.size(); ++i) {
            longestDistToSteinerAncestorOnChain[i].resize((int)chains[i].size());
            for (int j = 0; j < (int)chains[i].size(); ++j) {
                longestDistToSteinerAncestorOnChain[i][j] =
                        distanceToAncestor(chains[i][j], closestSteinerAncestor[chains[i][j]]);
                if (j > 0) {
                    Math::updateMax(longestDistToSteinerAncestorOnChain[i][j],
                            longestDistToSteinerAncestorOnChain[i][j - 1]);
                }
            }
        }
    }
 
    void computeLongestDistToSteinerAncestorSegTree() {
        longestDistToSteinerAncestorSegTree.resize((int)chains.size());
        for (int i = 0; i < (int)chains.size(); ++i) {
            longestDistToSteinerAncestorSegTree[i].resize(4 * (int)chains[i].size());
 
            std::vector<T> distanceToSteinerAncestor_byPosition;
            distanceToSteinerAncestor_byPosition.resize(chains[i].size());
            for (int j = 0; j < (int)chains[i].size(); ++j) {
                distanceToSteinerAncestor_byPosition[j] =
                        distanceToAncestor(chains[i][j], closestSteinerAncestor[chains[i][j]]);
            }
 
            buildMaxSegmentTree(longestDistToSteinerAncestorSegTree[i], 0, 0,
                    (int)chains[i].size() - 1, distanceToSteinerAncestor_byPosition);
        }
    }
 
    void dfsHeavyPathDecomposition(node v, node closestSteinerUpNode) {
        weightOfSubtree[v] = 1;
        node heaviestSon = nullptr;
        closestSteinerAncestor[v] = closestSteinerUpNode;
 
        for (adjEntry adj : v->adjEntries) {
            edge e = adj->theEdge();
            node son = adj->twinNode();
 
            if (weightOfSubtree[son] != 0) { // the parent
                continue;
            }
 
            nodeLevel[son] = nodeLevel[v] + 1;
            distanceToRoot[son] = distanceToRoot[v] + tree.weight(e);
 
            dfsHeavyPathDecomposition(son, v);
 
            fatherOfChain[chainOfNode[son]] = v;
            weightOfSubtree[v] += weightOfSubtree[son];
            if (heaviestSon == nullptr || weightOfSubtree[heaviestSon] < weightOfSubtree[son]) {
                heaviestSon = son;
            }
        }
 
        if (heaviestSon == nullptr) { // it's leaf => new chain
            OGDF_ASSERT((v->degree() <= 1));
 
            std::vector<node> new_chain;
            chains.push_back(new_chain);
            chainsOfTerminals.push_back(new_chain);
 
            fatherOfChain.push_back(nullptr);
            chainOfNode[v] = (int)chains.size() - 1;
        } else {
            chainOfNode[v] = chainOfNode[heaviestSon];
        }
 
        chains[chainOfNode[v]].push_back(v);
        chainsOfTerminals[chainOfNode[v]].push_back(v);
        positionOnChain[v] = (int)chains[chainOfNode[v]].size() - 1;
    }
 
    node binarySearchUpmostTerminal(node v, const std::vector<node>& chainOfTerminals) const {
        int left = 0, middle, right = (int)chainOfTerminals.size() - 1;
        while (left <= right) {
            middle = (left + right) >> 1;
 
            if (nodeLevel[chainOfTerminals[middle]] >= nodeLevel[v]) {
                right = middle - 1;
            } else {
                left = middle + 1;
            }
        }
 
        if (left == (int)chainOfTerminals.size()) {
            return nullptr;
        }
        return chainOfTerminals[left];
    }
 
    void computeBottleneckOnBranch(node x, node ancestor, T& longestPathDistance,
            T& fromLowestToAncestor) const {
        node upmostTerminal = x;
        while (closestSteinerAncestor[chains[chainOfNode[x]][0]] != nullptr
                && nodeLevel[closestSteinerAncestor[chains[chainOfNode[x]][0]]]
                        >= nodeLevel[ancestor]) {
            Math::updateMax(longestPathDistance,
                    longestDistToSteinerAncestorOnChain[chainOfNode[x]][positionOnChain[x]]);
 
            if (!chainsOfTerminals[chainOfNode[x]].empty()
                    && nodeLevel[chainsOfTerminals[chainOfNode[x]][0]] <= nodeLevel[x]) {
                upmostTerminal = chainsOfTerminals[chainOfNode[x]][0];
            }
            x = fatherOfChain[chainOfNode[x]];
        }
 
        // search the upmost terminal on the current chain that has level >= level[ancestor]
        node upmostTerminalLastChain =
                binarySearchUpmostTerminal(ancestor, chainsOfTerminals[chainOfNode[x]]);
        if (nodeLevel[upmostTerminalLastChain] > nodeLevel[x]) {
            upmostTerminalLastChain = nullptr;
        }
        if (upmostTerminalLastChain != nullptr) {
            upmostTerminal = upmostTerminalLastChain;
        }
 
        if (upmostTerminalLastChain != nullptr) {
            Math::updateMax(longestPathDistance,
                    getMaxSegmentTree(longestDistToSteinerAncestorSegTree[chainOfNode[x]], 0, 0,
                            static_cast<int>(chains[chainOfNode[x]].size()) - 1,
                            positionOnChain[upmostTerminalLastChain] + 1, positionOnChain[x]));
        }
 
        fromLowestToAncestor = distanceToAncestor(upmostTerminal, ancestor);
    }
 
public:
    HeavyPathDecomposition(const EdgeWeightedGraphCopy<T>& treeEWGraphCopy)
        : tree {treeEWGraphCopy}
        , chainOfNode {tree, -1}
        , positionOnChain {tree, -1}
        , weightOfSubtree {tree, 0}
        , nodeLevel {tree, 0}
        , distanceToRoot {tree, 0}
        , closestSteinerAncestor {tree, nullptr} {
        OGDF_ASSERT(!tree.empty());
        node root = tree.firstNode();
 
        dfsHeavyPathDecomposition(root, nullptr);
        fatherOfChain[chainOfNode[root]] = nullptr;
 
        // reverse the obtained chains
        int numberOfChains = static_cast<int>(chains.size());
        for (int i = 0; i < numberOfChains; ++i) {
            reverse(chains[i].begin(), chains[i].end());
            reverse(chainsOfTerminals[i].begin(), chainsOfTerminals[i].end());
            for (node v : chains[i]) {
                positionOnChain[v] = (int)chains[i].size() - 1 - positionOnChain[v];
            }
        }
 
        computeLongestDistToSteinerAncestorOnChain();
        computeLongestDistToSteinerAncestorSegTree();
    }
 
    node lowestCommonAncestor(node x, node y) const {
        while (chainOfNode[x] != chainOfNode[y]) {
            int xlevelOfFatherOfChain = fatherOfChain[chainOfNode[x]] != nullptr
                    ? nodeLevel[fatherOfChain[chainOfNode[x]]]
                    : -1;
            int ylevelOfFatherOfChain = fatherOfChain[chainOfNode[y]] != nullptr
                    ? nodeLevel[fatherOfChain[chainOfNode[y]]]
                    : -1;
 
            if (xlevelOfFatherOfChain >= ylevelOfFatherOfChain) {
                x = fatherOfChain[chainOfNode[x]];
            } else {
                y = fatherOfChain[chainOfNode[y]];
            }
        }
 
        if (nodeLevel[x] <= nodeLevel[y]) {
            return x;
        }
        return y;
    }
 
    T getBottleneckSteinerDistance(node x, node y) const {
        T xLongestPathDistance = 0, yLongestPathDistance = 0;
        T xFromLowestToLCA = 0, yFromLowestToLCA = 0;
        node xyLowestCommonAncestor = lowestCommonAncestor(x, y);
 
        computeBottleneckOnBranch(x, xyLowestCommonAncestor, xLongestPathDistance, xFromLowestToLCA);
        computeBottleneckOnBranch(y, xyLowestCommonAncestor, yLongestPathDistance, yFromLowestToLCA);
 
        T maxValue = max(xLongestPathDistance, yLongestPathDistance);
        Math::updateMax(maxValue, xFromLowestToLCA + yFromLowestToLCA);
 
        return maxValue;
    }
};
 
}
}