basicgridgeometry.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_DISCRETIZATION_BASIC_GRID_GEOMETRY_HH
#define DUMUX_DISCRETIZATION_BASIC_GRID_GEOMETRY_HH
 
#include <memory>
#include <utility>
#include <type_traits>
 
#include <dune/grid/common/mcmgmapper.hh>
 
#include <dumux/common/entitymap.hh>
#include <dumux/common/indextraits.hh>
#include <dumux/geometry/boundingboxtree.hh>
#include <dumux/geometry/geometricentityset.hh>
 
namespace Dumux {
 
template<class GV, class EM, class VM>
class BasicGridGeometry
{
    using ElementMap = EntityMap<GV, 0>;
    using ElementSet = GridViewGeometricEntitySet<GV, 0, EM>;
    using BoundingBoxTree = Dumux::BoundingBoxTree<ElementSet>;
 
    static constexpr int dim = GV::dimension;
    static constexpr int dimWorld = GV::dimensionworld;
 
    using GridIndexType = typename IndexTraits<GV>::GridIndex;
    using Element = typename GV::template Codim<0>::Entity;
 
public:
    using Grid = typename GV::Grid;
    using GridView = GV;
    using GlobalCoordinate = typename Element::Geometry::GlobalCoordinate;
    using ElementMapper = EM;
    using VertexMapper = VM;
 
    BasicGridGeometry(const GridView& gridView)
    : gridView_(gridView)
    , elementMapper_(makeElementMapper_(gridView))
    , vertexMapper_(makeVertexMapper_(gridView))
    , bBoxMin_(std::numeric_limits<double>::max())
    , bBoxMax_(-std::numeric_limits<double>::max())
    {
        computeGlobalBoundingBox_();
        update_();
    }
 
    void update(const GridView& gridView)
    {
        gridView_ = gridView;
        update_();
    }
 
    void update(GridView&& gridView)
    {
        gridView_ = std::move(gridView);
        update_();
    }
 
    const GridView& gridView() const
    { return gridView_; }
 
    const VertexMapper &vertexMapper() const
    { return vertexMapper_; }
 
    const ElementMapper &elementMapper() const
    { return elementMapper_; }
 
    VertexMapper &vertexMapper()
    { return vertexMapper_; }
 
    ElementMapper &elementMapper()
    { return elementMapper_; }
 
    const BoundingBoxTree& boundingBoxTree() const
    { return *boundingBoxTree_; }
 
    const ElementMap& elementMap() const
    { return *elementMap_; }
 
    Element element(GridIndexType eIdx) const
    { return elementMap()[eIdx]; }
 
    const GlobalCoordinate &bBoxMin() const
    { return bBoxMin_; }
 
    const GlobalCoordinate &bBoxMax() const
    { return bBoxMax_; }
 
private:
 
    ElementMapper makeElementMapper_(const GridView& gridView) const
    {
        if constexpr (std::is_constructible<ElementMapper, GridView, Dune::MCMGLayout>())
            return ElementMapper(gridView, Dune::mcmgElementLayout());
        else
            return ElementMapper(gridView);
    }
 
    VertexMapper makeVertexMapper_(const GridView& gridView) const
    {
        if constexpr (std::is_constructible<VertexMapper, GridView, Dune::MCMGLayout>())
            return VertexMapper(gridView, Dune::mcmgVertexLayout());
        else
            return VertexMapper(gridView);
    }
 
    void computeGlobalBoundingBox_()
    {
        // calculate the bounding box of the local partition of the grid view
        for (const auto& vertex : vertices(gridView_))
        {
            const auto pos = vertex.geometry().corner(0);
            for (int i = 0; i < dimWorld; ++i)
            {
                using std::min;
                using std::max;
                bBoxMin_[i] = min(bBoxMin_[i], pos[i]);
                bBoxMax_[i] = max(bBoxMax_[i], pos[i]);
            }
        }
 
        // communicate to get the bounding box of the whole domain
        if (gridView_.comm().size() > 1)
        {
            for (int i = 0; i < dimWorld; ++i)
            {
                bBoxMin_[i] = gridView_.comm().min(bBoxMin_[i]);
                bBoxMax_[i] = gridView_.comm().max(bBoxMax_[i]);
            }
        }
    }
 
    void update_()
    {
        // Update the mappers
        elementMapper_.update(gridView_);
        vertexMapper_.update(gridView_);
 
        // Compute the bounding box of the entire domain, for e.g. setting boundary conditions
        computeGlobalBoundingBox_();
 
        // update element map and bounding box tree
        // always building these comes at a memory overhead but improved
        // performance and thread-safe element level access (e.g. during assembly)
        // for all simulation that use these features
        elementMap_ = std::make_shared<ElementMap>(gridView_.grid(), elementMapper_);
        boundingBoxTree_ = std::make_unique<BoundingBoxTree>(
            std::make_shared<ElementSet>(gridView_, elementMapper(), elementMap_)
        );
    }
 
    GridView gridView_;
 
    ElementMapper elementMapper_;
    VertexMapper vertexMapper_;
 
    std::unique_ptr<const BoundingBoxTree> boundingBoxTree_;
 
    std::shared_ptr<const ElementMap> elementMap_;
 
    GlobalCoordinate bBoxMin_;
    GlobalCoordinate bBoxMax_;
};
 
} // end namespace Dumux
 
#endif