BoxMesh.h

Go to the documentation of this file.

#ifndef included_AMP_BoxMesh
#define included_AMP_BoxMesh
 
#include "AMP/geometry/Geometry.h"
#include "AMP/mesh/Mesh.h"
#include "AMP/mesh/MeshID.h"
#include "AMP/mesh/MeshIterator.h"
#include "AMP/utils/ArraySize.h"
#include "AMP/utils/Utilities.h"
 
namespace AMP::LinearAlgebra {
class Vector;
}
 
#include <array>
#include <map>
#include <memory>
#include <vector>
 
 
namespace AMP::Mesh {
 
class structuredMeshElement;
class structuredMeshIterator;
 
 
class BoxMesh : public AMP::Mesh::Mesh
{
public:
    class Box
    {
    public:
        constexpr explicit Box(
            int ifirst, int ilast, int jfirst = 0, int jlast = 0, int kfirst = 0, int klast = 0 );
        constexpr Box();                  
        constexpr ArraySize size() const; 
        int first[3];                     
        int last[3];                      
        std::string print() const;
 
    private:
    };
 
    class MeshElementIndex
    {
    public:
        constexpr MeshElementIndex() = default;
        constexpr explicit MeshElementIndex(
            GeomType type, uint8_t side, int x, int y = 0, int z = 0 );
        constexpr void reset();
        constexpr void reset( GeomType type, uint8_t side, int x, int y = 0, int z = 0 );
        constexpr bool isNull() const { return d_side == 255; }
        constexpr bool operator==( const MeshElementIndex &rhs ) const; 
        constexpr bool operator!=( const MeshElementIndex &rhs ) const; 
        constexpr bool operator>( const MeshElementIndex &rhs ) const;  
        constexpr bool operator>=( const MeshElementIndex &rhs ) const; 
        constexpr bool operator<( const MeshElementIndex &rhs ) const;  
        constexpr bool operator<=( const MeshElementIndex &rhs ) const; 
        constexpr const auto &index() const { return d_index; }
        constexpr int index( int d ) const { return d_index[d]; }
        constexpr int &index( int d ) { return d_index[d]; }
        constexpr int operator[]( int d ) const { return d_index[d]; }
        constexpr int &operator[]( int d ) { return d_index[d]; }
        constexpr GeomType type() const { return static_cast<GeomType>( d_type ); }
        constexpr uint8_t side() const { return d_side; }
        static constexpr size_t numElements( const MeshElementIndex &first,
                                             const MeshElementIndex &last );
        std::string print() const;
 
    private:
        uint8_t d_type             = 0;           
        uint8_t d_side             = 255;         
        std::array<int, 3> d_index = { 0, 0, 0 }; 
        friend class BoxMesh;
        friend class structuredMeshElement;
    };
 
    class MeshElementIndexIterator final
    {
    public: // iterator_traits
        using iterator_category = std::random_access_iterator_tag;
        using value_type        = MeshElementIndex;
        using difference_type   = ptrdiff_t;
        using pointer           = const MeshElementIndex *;
        using reference         = const MeshElementIndex &;
 
    public:
        MeshElementIndexIterator() = default;
        MeshElementIndexIterator( const MeshElementIndex &first,
                                  const MeshElementIndex &last,
                                  const AMP::Mesh::BoxMesh *mesh,
                                  size_t pos = 0 );
        MeshElementIndexIterator &operator++();
        MeshElementIndexIterator &operator--();
        MeshElementIndexIterator &operator+=( int N );
        MeshElementIndexIterator &operator+=( const MeshElementIndexIterator &it );
        MeshElementIndexIterator &operator[]( int );
        MeshElementIndexIterator begin() const;
        MeshElementIndexIterator end() const;
        MeshElementIndex operator*() const;
        bool operator==( const MeshElementIndexIterator &rhs ) const;
        bool operator!=( const MeshElementIndexIterator &rhs ) const;
        bool empty() const { return d_size == 0; }
        void set( uint32_t i ) { d_pos = i; }
        size_t size() const { return d_size; }
        size_t position() const { return d_pos; }
        inline const auto &first() const { return d_first; }
        inline const auto &last() const { return d_last; }
 
    private:
        // Data members
        bool d_checkBoundary             = false;
        std::array<bool, 3> d_isPeriodic = { false, false, false };
        std::array<int, 3> d_globalSize  = { 0, 0, 0 };
        uint32_t d_pos                   = { 0 };
        uint32_t d_size                  = { 0 };
        MeshElementIndex d_first;
        MeshElementIndex d_last;
    };
 
 
public:
    static std::shared_ptr<BoxMesh> generate( std::shared_ptr<const MeshParameters> params );
 
 
    std::unique_ptr<Mesh> clone() const override = 0;
 
 
    static size_t estimateMeshSize( std::shared_ptr<const MeshParameters> params );
 
    static ArraySize estimateLogicalMeshSize( std::shared_ptr<const MeshParameters> params );
 
 
    static size_t maxProcs( std::shared_ptr<const MeshParameters> params );
 
 
    virtual ~BoxMesh();
 
 
    /* Return the number of local element of the given type
     * \param type   Geometric type
     */
    size_t numLocalElements( const GeomType type ) const override final;
 
 
    /* Return the global number of elements of the given type
     * Note: depending on the mesh this routine may require global communication across the mesh.
     * \param type   Geometric type
     */
    size_t numGlobalElements( const GeomType type ) const override final;
 
 
    /* Return the number of ghost elements of the given type on the current processor
     * \param type   Geometric type
     */
    size_t numGhostElements( const GeomType type, const int gcw ) const override final;
 
 
    MeshIterator getIterator( const GeomType type, const int gcw = 0 ) const override final;
 
 
    virtual MeshIterator getSurfaceIterator( const GeomType type,
                                             const int gcw = 0 ) const override final;
 
 
    std::vector<int> getBoundaryIDs() const override final;
 
 
    virtual MeshIterator getBoundaryIDIterator( const GeomType type,
                                                const int id,
                                                const int gcw = 0 ) const override final;
 
    std::vector<int> getBlockIDs() const override final;
 
 
    virtual MeshIterator
    getBlockIDIterator( const GeomType type, const int id, const int gcw = 0 ) const override final;
 
 
    MeshElementPtr getElement( const MeshElementID &id ) const override final;
 
 
    virtual ElementListPtr getElementParents( const MeshElement &elem,
                                              const GeomType type ) const override final;
 
 
    inline Box getGlobalBox( int gcw = 0 ) const;
 
 
    inline Box getLocalBox( int gcw = 0 ) const;
 
 
    inline std::vector<bool> periodic() const;
 
    inline ArraySize size() const;
 
    inline ArraySize localSize() const;
 
    inline std::vector<size_t> numBlocks() const;
 
 
    bool operator==( const Mesh &mesh ) const override;
 
 
    structuredMeshElement getElement( const MeshElementIndex &index ) const;
 
    virtual void coord( const MeshElementIndex &index, double *pos ) const = 0;
 
    virtual std::array<AMP::Array<double>, 3> localCoord() const = 0;
 
    virtual std::array<AMP::Array<double>, 3> globalCoord() const = 0;
 
    virtual bool isOnSurface( const MeshElementIndex &index ) const;
 
    virtual bool isOnBoundary( const MeshElementIndex &index, int id ) const;
 
 
public: // BoxMesh specific functionality
    virtual AMP::Geometry::Point physicalToLogical( const AMP::Geometry::Point &x ) const = 0;
 
    MeshElementIndex getElementFromLogical( const AMP::Geometry::Point &x, GeomType type ) const;
 
    MeshElementIndex getElementFromPhysical( const AMP::Geometry::Point &x, GeomType type ) const;
 
    inline int getRank( const MeshElementIndex &id ) const;
 
    inline MeshElementID convert( const MeshElementIndex &id ) const;
 
    inline MeshElementIndex convert( const MeshElementID &id ) const;
 
 
    std::shared_ptr<AMP::LinearAlgebra::Vector> createVector( const std::string &name,
                                                              int gcw = 0 );
 
    inline void fixPeriodic( MeshElementIndex &id ) const;
 
public: // Convenience typedef
    typedef AMP::Utilities::stackVector<std::pair<MeshElementIndex, MeshElementIndex>, 32>
        ElementBlocks;
 
 
public: // Advanced functions
    // Get the surface id for a given surface
    int getSurfaceID( int surface ) const;
 
    // Get the surface set for a given surface/type
    ElementBlocks getSurface( int surface, GeomType type ) const;
 
    // Helper function to return the indices of the local block owned by the given processor
    inline std::array<int, 6> getLocalBlock( int rank ) const;
 
    // Helper functions to identify the iterator blocks
    ElementBlocks
    getIteratorRange( std::array<int, 6> range, const GeomType type, const int gcw ) const;
    ElementBlocks intersect( const ElementBlocks &v1, const ElementBlocks &v2 ) const;
 
    // Helper function to create an iterator from an ElementBlocks list
    MeshIterator createIterator( const ElementBlocks &list ) const;
 
 
protected:
    // Constructor
    explicit BoxMesh( std::shared_ptr<const MeshParameters> );
    explicit BoxMesh( const BoxMesh & );
    BoxMesh &operator=( const BoxMesh & ) = delete;
 
    // Function to create the load balancing
    static void loadBalance( std::array<int, 3> size,
                             int N_procs,
                             std::vector<int> *startIndex,
                             std::vector<int> minSize = {} );
 
    // Function to initialize the mesh data once the logical mesh info has been created
    void initialize( const std::array<int, 3> &size,
                     const std::array<int, 6> &ids,
                     const std::vector<int> &minSize = {} );
 
    // Function to finalize the mesh data once the coordinates have been set
    void finalize( const std::string &name, const std::vector<double> &displacement );
 
    // Function to finalize the mesh data once the coordinates have been set
    std::vector<double> getDisplacement( std::shared_ptr<const AMP::Database> db );
 
    // Function to finalize the mesh data once the coordinates have been set
    virtual void createBoundingBox();
 
    // Helper function to fill the node data for a uniform cartesian mesh
    static void fillCartesianNodes( int dim,
                                    const int *globalSize,
                                    const double *range,
                                    const std::vector<MeshElementIndex> &index,
                                    std::vector<double> *coord );
 
    // Helper function to check if an element is on a given side
    bool onSide( const MeshElementIndex &index, int d, int s ) const;
 
    // Get the surface (including mapped surfaces)
    ElementBlocks getSurface2( int surface, GeomType type ) const;
 
    // Helper function to map points on the boundary
    template<uint8_t NDIM>
    std::vector<MeshElementIndex> createMap( const std::vector<MeshElementIndex> ) const;
    void createMaps();
 
protected: // Write/read restart data
    void writeRestart( int64_t ) const override;
    BoxMesh( int64_t, AMP::IO::RestartManager * );
 
 
protected: // Friend functions to access protected functions
    friend class structuredMeshElement;
    friend class structuredMeshIterator;
    typedef std::vector<MeshElementIndex> IndexList;
    typedef std::shared_ptr<IndexList> ListPtr;
    typedef std::tuple<IndexList, IndexList> SurfaceMapStruct;
 
 
protected:                                     // Internal data
    const int d_rank, d_size;                  // Cached values for the rank and size
    int d_blockID;                             // Block id for the mesh
    const std::array<int, 3> d_globalSize;     // The size of the logical domain in each direction
    const std::array<int, 3> d_numBlocks;      // The number of local box in each direction
    const std::vector<int> d_startIndex[3];    // The first index for each block
    const std::vector<int> d_endIndex[3];      // The end index (last=1) for each block
    const std::array<int, 6> d_localIndex;     // Local index range (cached for performance)
    const std::array<int, 3> d_indexSize;      // Local index size (local box + 2) (cached)
    const std::array<int, 6> d_surfaceId;      // ID of each surface (if any, -1 if not)
    mutable std::vector<ListPtr> d_surface[4]; // List of surface elements
    mutable std::vector<ListPtr> d_bnd[4][6];  // List of boundary elements
    SurfaceMapStruct d_surfaceMaps[3];         // Surface maps
 
protected:
    BoxMesh();
};
 
 
// Function to write an index to std::ostream
std::ostream &operator<<( std::ostream &out, const BoxMesh::MeshElementIndex &x );
 
 
} // namespace AMP::Mesh
 
#include "AMP/mesh/structured/BoxMesh.inline.h"
 
 
#endif