distributedboundingboxtree.hh

Go to the documentation of this file.

// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
// vi: set et ts=4 sw=4 sts=4:
//
// SPDX-FileCopyrightText: Copyright © DuMux Project contributors, see AUTHORS.md in root folder
// SPDX-License-Identifier: GPL-3.0-or-later
//
#ifndef DUMUX_GEOMETRY_DISTRIBUTED_BOUNDINGBOXTREE_HH
#define DUMUX_GEOMETRY_DISTRIBUTED_BOUNDINGBOXTREE_HH
 
#include <cassert>
#include <vector>
#include <array>
#include <span>
#include <memory>
#include <algorithm>
#include <type_traits>
#include <utility>
 
#include <dune/common/fvector.hh>
#include <dune/common/exceptions.hh>
 
#include <dumux/geometry/boundingboxgeometry.hh>
#include <dumux/geometry/boundingboxtree.hh>
#include <dumux/geometry/geometricentityset.hh>
 
namespace Dumux {

template<class GeometricEntitySet>
  requires requires(const GeometricEntitySet& set) { { set.comm() }; }
class DistributedBoundingBoxTree
{
    static constexpr int dimworld = GeometricEntitySet::dimensionworld;
    using ctype = typename GeometricEntitySet::ctype;
 
    // internal types for the replicated process tree (the tree over the per-process boxes,
    // each leaf being a process partition box carrying its owning rank)
    using ProcessBoxGeometry = BoundingBoxGeometry<ctype, dimworld>;
    using ProcessEntitySet = ProcessGeometricEntitySet<ctype, dimworld>;
    using ProcessEntity = typename ProcessEntitySet::Entity;
 
public:
    using LocalTree = BoundingBoxTree<GeometricEntitySet>;
    using EntitySet = GeometricEntitySet;
    using Communication = std::decay_t<decltype(std::declval<const GeometricEntitySet>().comm())>;
    using ProcessTree = BoundingBoxTree<ProcessEntitySet>;

    DistributedBoundingBoxTree() = default;

    explicit DistributedBoundingBoxTree(std::shared_ptr<const GeometricEntitySet> set)
    { build(set); }

    void build(std::shared_ptr<const GeometricEntitySet> set)
    {
        entitySet_ = set;
        comm_ = set->comm();
        const int numProc = comm_.size();
 
        // build the process-local tree (guarding against empty partitions for
        // which the sequential tree is not well-defined)
        const bool hasLocalEntities = set->size() > 0;
        if (hasLocalEntities)
            localTree_ = std::make_shared<const LocalTree>(set);
        else
            localTree_.reset();
 
        // determine the bounding box of the local partition (the root box of
        // the local tree) or a sentinel for empty partitions
        std::array<ctype, 2*dimworld> localBox{};
        if (hasLocalEntities)
            std::copy_n(localBoundingBox().data(), 2*dimworld, localBox.begin());
 
        // communicate the per-process bounding boxes to all processes
        std::vector<ctype> processBoxes(2*dimworld*numProc, 0.0);
        comm_.allgather(localBox.data(), 2*dimworld, processBoxes.data());
 
        // communicate which processes actually hold entities (valid boxes)
        const int localValid = hasLocalEntities ? 1 : 0;
        std::vector<int> processValid(numProc, 0);
        comm_.allgather(&localValid, 1, processValid.data());
 
        // build the replicated process tree over the non-empty process boxes;
        // each entity carries the rank it belongs to
        std::vector<ProcessEntity> processEntities;
        processEntities.reserve(numProc);
        for (int p = 0; p < numProc; ++p)
        {
            if (!processValid[p])
                continue;
 
            const ctype* b = processBoxes.data() + std::size_t(p)*2*dimworld;
            Dune::FieldVector<ctype, dimworld> lower, upper;
            for (int d = 0; d < dimworld; ++d)
            {
                lower[d] = b[d];
                upper[d] = b[d + dimworld];
            }
            // the entity's consecutive index is its position in the set
            processEntities.emplace_back(ProcessBoxGeometry{lower, upper}, processEntities.size(), p);
        }
 
        // the global entity set is assumed non-empty (at least one process holds
        // entities), just like for the sequential bounding box tree
        assert(!processEntities.empty());
        processTree_ = std::make_shared<const ProcessTree>(
            std::make_shared<const ProcessEntitySet>(std::move(processEntities)));
    }

    const EntitySet& entitySet() const
    { return *entitySet_; }

    const Communication& comm() const
    { return comm_; }

    int rank() const
    { return comm_.rank(); }

    bool hasLocalEntities() const
    { return static_cast<bool>(localTree_); }

    const LocalTree& localTree() const
    {
        if (!localTree_)
            DUNE_THROW(Dune::InvalidStateException, "No local bounding box tree on a process without entities");
        return *localTree_;
    }

    const ProcessTree& processTree() const
    { return *processTree_; }

    int processForLeaf(std::size_t leafIdx) const
    { return processTree_->entitySet().rank(leafIdx); }

    std::span<const ctype, 2*dimworld> localBoundingBox() const
    {
        if (!localTree_)
            DUNE_THROW(Dune::InvalidStateException, "No local bounding box on a process without entities");
        const auto* root = localTree_->getBoundingBoxCoordinates(localTree_->numBoundingBoxes() - 1);
        return std::span<const ctype, 2*dimworld>{root, 2*dimworld};
    }
 
private:
    std::shared_ptr<const EntitySet> entitySet_;

    std::shared_ptr<const LocalTree> localTree_;

    std::shared_ptr<const ProcessTree> processTree_;

    Communication comm_{};
};
 
} // end namespace Dumux
 
#endif