MeshGraph.h

Go to the documentation of this file.

////////////////////////////////////////////////////////////////////////////////
//
//  File: MeshGraph.h
//
//  For more information, please see: http://www.nektar.info/
//
//  The MIT License
//
//  Copyright (c) 2006 Division of Applied Mathematics, Brown University (USA),
//  Department of Aeronautics, Imperial College London (UK), and Scientific
//  Computing and Imaging Institute, University of Utah (USA).
//
//  Permission is hereby granted, free of charge, to any person obtaining a
//  copy of this software and associated documentation files (the "Software"),
//  to deal in the Software without restriction, including without limitation
//  the rights to use, copy, modify, merge, publish, distribute, sublicense,
//  and/or sell copies of the Software, and to permit persons to whom the
//  Software is furnished to do so, subject to the following conditions:
//
//  The above copyright notice and this permission notice shall be included
//  in all copies or substantial portions of the Software.
//
//  THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
//  OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
//  FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
//  THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
//  LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
//  FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
//  DEALINGS IN THE SOFTWARE.
//
//  Description:
//
////////////////////////////////////////////////////////////////////////////////
 
#ifndef NEKTAR_SPATIALDOMAINS_MESHGRAPH_H
#define NEKTAR_SPATIALDOMAINS_MESHGRAPH_H
 
#include <unordered_map>
 
#include <LibUtilities/BasicUtils/DomainRange.h>
#include <LibUtilities/BasicUtils/FieldIO.h>
#include <LibUtilities/BasicUtils/SessionReader.h>
 
#include <SpatialDomains/HexGeom.h>
#include <SpatialDomains/MeshEntities.hpp>
#include <SpatialDomains/PrismGeom.h>
#include <SpatialDomains/PyrGeom.h>
#include <SpatialDomains/QuadGeom.h>
#include <SpatialDomains/SegGeom.h>
#include <SpatialDomains/TetGeom.h>
#include <SpatialDomains/TriGeom.h>
 
#include <SpatialDomains/Curve.hpp>
#include <SpatialDomains/SpatialDomainsDeclspec.h>
 
class TiXmlDocument;
 
namespace Nektar
{
namespace SpatialDomains
{
typedef std::map<int, std::pair<LibUtilities::ShapeType, std::vector<int>>>
    CompositeDescriptor;
 
enum ExpansionType
{
    eNoExpansionType,
    eModified,
    eModifiedQuadPlus1,
    eModifiedQuadPlus2,
    eModifiedGLLRadau10,
    eOrthogonal,
    eGLL_Lagrange,
    eGLL_Lagrange_SEM,
    eGauss_Lagrange,
    eGauss_Lagrange_SEM,
    eFourier,
    eFourierSingleMode,
    eFourierHalfModeRe,
    eFourierHalfModeIm,
    eChebyshev,
    eFourierChebyshev,
    eChebyshevFourier,
    eFourierModified,
    eExpansionTypeSize
};
 
// Keep this consistent with the enums in ExpansionType.
// This is used in the BC file to specify the expansion type.
const std::string kExpansionTypeStr[] = {"NOTYPE",
                                         "MODIFIED",
                                         "MODIFIEDQUADPLUS1",
                                         "MODIFIEDQUADPLUS2",
                                         "MODIFIEDGLLRADAU10",
                                         "ORTHOGONAL",
                                         "GLL_LAGRANGE",
                                         "GLL_LAGRANGE_SEM",
                                         "GAUSS_LAGRANGE",
                                         "GAUSS_LAGRANGE_SEM",
                                         "FOURIER",
                                         "FOURIERSINGLEMODE",
                                         "FOURIERHALFMODERE",
                                         "FOURIERHALFMODEIM",
                                         "CHEBYSHEV",
                                         "FOURIER-CHEBYSHEV",
                                         "CHEBYSHEV-FOURIER",
                                         "FOURIER-MODIFIED"};
 
typedef std::map<int, std::vector<unsigned int>> CompositeOrdering;
typedef std::map<int, std::vector<unsigned int>> BndRegionOrdering;
 
// set restriction on domain range for post-processing.
// struct DomainRange
// {
//     bool m_doXrange;
//     NekDouble m_xmin;
//     NekDouble m_xmax;
//     bool m_doYrange;
//     NekDouble m_ymin;
//     NekDouble m_ymax;
//     bool m_doZrange;
//     NekDouble m_zmin;
//     NekDouble m_zmax;
 
//     bool m_checkShape;
//     LibUtilities::ShapeType m_shapeType;
// };
 
// typedef std::shared_ptr<DomainRange> DomainRangeShPtr;
// static DomainRangeShPtr NullDomainRangeShPtr;
 
struct Composite
{
    std::vector<std::shared_ptr<Geometry>> m_geomVec;
};
 
typedef std::shared_ptr<Composite> CompositeSharedPtr;
typedef std::map<int, CompositeSharedPtr> CompositeMap;
 
struct ExpansionInfo;
 
typedef std::shared_ptr<ExpansionInfo> ExpansionInfoShPtr;
typedef std::map<int, ExpansionInfoShPtr> ExpansionInfoMap;
 
typedef std::shared_ptr<ExpansionInfoMap> ExpansionInfoMapShPtr;
typedef std::map<std::string, ExpansionInfoMapShPtr> ExpansionInfoMapShPtrMap;
 
struct ExpansionInfo
{
    ExpansionInfo(GeometrySharedPtr geomShPtr,
                  const LibUtilities::BasisKeyVector basiskeyvec)
        : m_geomShPtr(geomShPtr), m_basisKeyVector(basiskeyvec)
    {
    }
 
    ExpansionInfo(ExpansionInfoShPtr ExpInfo)
        : m_geomShPtr(ExpInfo->m_geomShPtr),
          m_basisKeyVector(ExpInfo->m_basisKeyVector)
    {
    }
 
    GeometrySharedPtr m_geomShPtr;
    LibUtilities::BasisKeyVector m_basisKeyVector;
};
 
typedef std::map<std::string, std::string> GeomInfoMap;
typedef std::shared_ptr<std::vector<std::pair<GeometrySharedPtr, int>>>
    GeometryLinkSharedPtr;
 
// Forward declaration
class RefRegion;
 
typedef std::map<std::string, std::string> MeshMetaDataMap;
 
class MeshGraph;
typedef std::shared_ptr<MeshGraph> MeshGraphSharedPtr;
 
class Movement;
typedef std::shared_ptr<Movement> MovementSharedPtr;
 
/// Base class for a spectral/hp element mesh.
class MeshGraph
{
public:
    SPATIAL_DOMAINS_EXPORT MeshGraph();
    SPATIAL_DOMAINS_EXPORT virtual ~MeshGraph();
 
    SPATIAL_DOMAINS_EXPORT static MeshGraphSharedPtr Read(
        const LibUtilities::SessionReaderSharedPtr pSession,
        LibUtilities::DomainRangeShPtr rng = LibUtilities::NullDomainRangeShPtr,
        bool fillGraph                     = true,
        SpatialDomains::MeshGraphSharedPtr partitionedGraph = nullptr);
 
    SPATIAL_DOMAINS_EXPORT void WriteGeometry(
        std::string &outfilename, bool defaultExp = false,
        const LibUtilities::FieldMetaDataMap &metadata =
            LibUtilities::NullFieldMetaDataMap);
 
    void Empty(int dim, int space)
    {
        m_meshDimension  = dim;
        m_spaceDimension = space;
    }
 
    /*transfers the minial data structure to full meshgraph*/
    SPATIAL_DOMAINS_EXPORT void FillGraph();
 
    SPATIAL_DOMAINS_EXPORT void FillBoundingBoxTree();
 
    SPATIAL_DOMAINS_EXPORT std::vector<int> GetElementsContainingPoint(
        PointGeomSharedPtr p);
 
    ////////////////////
    SPATIAL_DOMAINS_EXPORT void ReadExpansionInfo();
 
    /// Read refinement info.
    SPATIAL_DOMAINS_EXPORT void ReadRefinementInfo();
 
    /* ---- Helper functions ---- */
    /// Dimension of the mesh (can be a 1D curve in 3D space).
    int GetMeshDimension()
    {
        return m_meshDimension;
    }
 
    /// Dimension of the space (can be a 1D curve in 3D space).
    int GetSpaceDimension()
    {
        return m_spaceDimension;
    }
 
    /* Range definitions for postprorcessing */
    SPATIAL_DOMAINS_EXPORT void SetDomainRange(
        NekDouble xmin, NekDouble xmax,
        NekDouble ymin = NekConstants::kNekUnsetDouble,
        NekDouble ymax = NekConstants::kNekUnsetDouble,
        NekDouble zmin = NekConstants::kNekUnsetDouble,
        NekDouble zmax = NekConstants::kNekUnsetDouble);
 
    /// Check if goemetry is in range definition if activated
    SPATIAL_DOMAINS_EXPORT bool CheckRange(Geometry2D &geom);
 
    /// Check if goemetry is in range definition if activated
    SPATIAL_DOMAINS_EXPORT bool CheckRange(Geometry3D &geom);
 
    /* ---- Composites and Domain ---- */
    CompositeSharedPtr GetComposite(int whichComposite)
    {
        ASSERTL0(m_meshComposites.find(whichComposite) !=
                     m_meshComposites.end(),
                 "Composite not found.");
        return m_meshComposites.find(whichComposite)->second;
    }
 
    SPATIAL_DOMAINS_EXPORT GeometrySharedPtr
    GetCompositeItem(int whichComposite, int whichItem);
 
    SPATIAL_DOMAINS_EXPORT void GetCompositeList(
        const std::string &compositeStr, CompositeMap &compositeVector) const;
 
    std::map<int, CompositeSharedPtr> &GetComposites()
    {
        return m_meshComposites;
    }
 
    std::map<int, std::string> &GetCompositesLabels()
    {
        return m_compositesLabels;
    }
 
    std::map<int, std::map<int, CompositeSharedPtr>> &GetDomain()
    {
        return m_domain;
    }
 
    std::map<int, CompositeSharedPtr> &GetDomain(int domain)
    {
        ASSERTL1(m_domain.count(domain),
                 "Request for domain which does not exist");
        return m_domain[domain];
    }
 
    SPATIAL_DOMAINS_EXPORT const ExpansionInfoMap &GetExpansionInfo(
        const std::string variable = "DefaultVar");
 
    SPATIAL_DOMAINS_EXPORT ExpansionInfoShPtr GetExpansionInfo(
        GeometrySharedPtr geom, const std::string variable = "DefaultVar");
 
    /// Sets expansions given field definitions
    SPATIAL_DOMAINS_EXPORT void SetExpansionInfo(
        std::vector<LibUtilities::FieldDefinitionsSharedPtr> &fielddef);
 
    /// Sets expansions given field definition, quadrature points.
    SPATIAL_DOMAINS_EXPORT void SetExpansionInfo(
        std::vector<LibUtilities::FieldDefinitionsSharedPtr> &fielddef,
        std::vector<std::vector<LibUtilities::PointsType>> &pointstype);
 
    /// Sets expansions to have equispaced points
    SPATIAL_DOMAINS_EXPORT void SetExpansionInfoToEvenlySpacedPoints(
        int npoints = 0);
 
    /// Reset expansion to have specified polynomial order \a nmodes
    SPATIAL_DOMAINS_EXPORT void SetExpansionInfoToNumModes(int nmodes);
 
    /// Reset expansion to have specified point order \a
    /// npts
    SPATIAL_DOMAINS_EXPORT void SetExpansionInfoToPointOrder(int npts);
    /// This function sets the expansion #exp in map with
    /// entry #variable
 
    /// Set refinement info.
    SPATIAL_DOMAINS_EXPORT void SetRefinementInfo(
        ExpansionInfoMapShPtr &expansionMap);
 
    /// Perform the p-refinement in the selected elements
    SPATIAL_DOMAINS_EXPORT void PRefinementElmts(
        ExpansionInfoMapShPtr &expansionMap, RefRegion *&region,
        GeometrySharedPtr geomVecIter);
 
    inline void SetExpansionInfo(const std::string variable,
                                 ExpansionInfoMapShPtr &exp);
 
    inline void SetSession(LibUtilities::SessionReaderSharedPtr pSession);
 
    /// Sets the basis key for all expansions of the given shape.
    SPATIAL_DOMAINS_EXPORT void SetBasisKey(LibUtilities::ShapeType shape,
                                            LibUtilities::BasisKeyVector &keys,
                                            std::string var = "DefaultVar");
 
    SPATIAL_DOMAINS_EXPORT void ResetExpansionInfoToBasisKey(
        ExpansionInfoMapShPtr &expansionMap, LibUtilities::ShapeType shape,
        LibUtilities::BasisKeyVector &keys);
 
    inline bool SameExpansionInfo(const std::string var1,
                                  const std::string var2);
 
    inline bool ExpansionInfoDefined(const std::string var);
 
    inline bool CheckForGeomInfo(std::string parameter);
 
    inline const std::string GetGeomInfo(std::string parameter);
 
    SPATIAL_DOMAINS_EXPORT static LibUtilities::BasisKeyVector
    DefineBasisKeyFromExpansionType(GeometrySharedPtr in, ExpansionType type,
                                    const int order);
 
    SPATIAL_DOMAINS_EXPORT LibUtilities::BasisKeyVector
    DefineBasisKeyFromExpansionTypeHomo(
        GeometrySharedPtr in, ExpansionType type_x, ExpansionType type_y,
        ExpansionType type_z, const int nummodes_x, const int nummodes_y,
        const int nummodes_z);
 
    /* ---- Manipulation of mesh ---- */
    int GetNvertices()
    {
        return m_vertSet.size();
    }
 
    PointGeomSharedPtr GetVertex(int id)
    {
        return m_vertSet[id];
    }
 
    SPATIAL_DOMAINS_EXPORT SegGeomSharedPtr GetSegGeom(int id)
    {
        return m_segGeoms[id];
    }
 
    SPATIAL_DOMAINS_EXPORT CurveMap &GetCurvedEdges()
    {
        return m_curvedEdges;
    }
    SPATIAL_DOMAINS_EXPORT CurveMap &GetCurvedFaces()
    {
        return m_curvedFaces;
    }
 
    SPATIAL_DOMAINS_EXPORT std::map<int, PointGeomSharedPtr> &GetAllPointGeoms()
    {
        return m_vertSet;
    }
    SPATIAL_DOMAINS_EXPORT std::map<int, SegGeomSharedPtr> &GetAllSegGeoms()
    {
        return m_segGeoms;
    }
    SPATIAL_DOMAINS_EXPORT TriGeomMap &GetAllTriGeoms()
    {
        return m_triGeoms;
    }
    SPATIAL_DOMAINS_EXPORT QuadGeomMap &GetAllQuadGeoms()
    {
        return m_quadGeoms;
    }
    SPATIAL_DOMAINS_EXPORT TetGeomMap &GetAllTetGeoms()
    {
        return m_tetGeoms;
    }
    SPATIAL_DOMAINS_EXPORT PyrGeomMap &GetAllPyrGeoms()
    {
        return m_pyrGeoms;
    }
    SPATIAL_DOMAINS_EXPORT PrismGeomMap &GetAllPrismGeoms()
    {
        return m_prismGeoms;
    }
    SPATIAL_DOMAINS_EXPORT HexGeomMap &GetAllHexGeoms()
    {
        return m_hexGeoms;
    }
    SPATIAL_DOMAINS_EXPORT std::unordered_map<int, GeometryLinkSharedPtr>
        &GetAllFaceToElMap()
    {
        return m_faceToElMap;
    }
 
    SPATIAL_DOMAINS_EXPORT int GetNumElements();
 
    Geometry2DSharedPtr GetGeometry2D(int gID)
    {
        auto it1 = m_triGeoms.find(gID);
        if (it1 != m_triGeoms.end())
            return it1->second;
 
        auto it2 = m_quadGeoms.find(gID);
        if (it2 != m_quadGeoms.end())
            return it2->second;
 
        return Geometry2DSharedPtr();
    };
 
    SPATIAL_DOMAINS_EXPORT GeometryLinkSharedPtr
    GetElementsFromEdge(Geometry1DSharedPtr edge);
 
    SPATIAL_DOMAINS_EXPORT GeometryLinkSharedPtr
    GetElementsFromFace(Geometry2DSharedPtr face);
 
    void SetPartition(SpatialDomains::MeshGraphSharedPtr graph);
 
    CompositeOrdering &GetCompositeOrdering()
    {
        return m_compOrder;
    }
 
    void SetCompositeOrdering(CompositeOrdering p_compOrder)
    {
        m_compOrder = p_compOrder;
    }
 
    BndRegionOrdering &GetBndRegionOrdering()
    {
        return m_bndRegOrder;
    }
 
    void SetBndRegionOrdering(BndRegionOrdering p_bndRegOrder)
    {
        m_bndRegOrder = p_bndRegOrder;
    }
 
    /*an inital read which loads a very light weight data structure*/
    SPATIAL_DOMAINS_EXPORT void ReadGeometry(LibUtilities::DomainRangeShPtr rng,
                                             bool fillGraph);
    SPATIAL_DOMAINS_EXPORT void PartitionMesh(
        LibUtilities::SessionReaderSharedPtr session);
 
    SPATIAL_DOMAINS_EXPORT std::map<int, MeshEntity> CreateMeshEntities();
    SPATIAL_DOMAINS_EXPORT CompositeDescriptor CreateCompositeDescriptor();
 
    SPATIAL_DOMAINS_EXPORT inline MovementSharedPtr &GetMovement()
    {
        return m_movement;
    }
 
protected:
    SPATIAL_DOMAINS_EXPORT virtual void v_WriteGeometry(
        std::string &outfilename, bool defaultExp = false,
        const LibUtilities::FieldMetaDataMap &metadata =
            LibUtilities::NullFieldMetaDataMap) = 0;
    SPATIAL_DOMAINS_EXPORT virtual void v_ReadGeometry(
        LibUtilities::DomainRangeShPtr rng, bool fillGraph) = 0;
    SPATIAL_DOMAINS_EXPORT virtual void v_PartitionMesh(
        LibUtilities::SessionReaderSharedPtr session) = 0;
 
    void PopulateFaceToElMap(Geometry3DSharedPtr element, int kNfaces);
    ExpansionInfoMapShPtr SetUpExpansionInfoMap();
    std::string GetCompositeString(CompositeSharedPtr comp);
 
    LibUtilities::SessionReaderSharedPtr m_session;
    PointGeomMap m_vertSet;
 
    CurveMap m_curvedEdges;
    CurveMap m_curvedFaces;
 
    SegGeomMap m_segGeoms;
 
    TriGeomMap m_triGeoms;
    QuadGeomMap m_quadGeoms;
    TetGeomMap m_tetGeoms;
    PyrGeomMap m_pyrGeoms;
    PrismGeomMap m_prismGeoms;
    HexGeomMap m_hexGeoms;
 
    int m_meshDimension;
    int m_spaceDimension;
    int m_partition;
    bool m_meshPartitioned = false;
    bool m_useExpansionType;
 
    // Refinement attributes (class members)
    /// Link the refinement id with the composites
    std::map<int, CompositeMap> m_refComposite;
    // std::map<int, LibUtilities::BasisKeyVector> m_refBasis;
    /// Link the refinement id with the surface region data
    std::map<int, RefRegion *> m_refRegion;
    bool m_refFlag = false;
 
    CompositeMap m_meshComposites;
    std::map<int, std::string> m_compositesLabels;
    std::map<int, CompositeMap> m_domain;
    LibUtilities::DomainRangeShPtr m_domainRange;
 
    ExpansionInfoMapShPtrMap m_expansionMapShPtrMap;
 
    std::unordered_map<int, GeometryLinkSharedPtr> m_faceToElMap;
 
    TiXmlElement *m_xmlGeom;
 
    CompositeOrdering m_compOrder;
    BndRegionOrdering m_bndRegOrder;
 
    struct GeomRTree;
    std::unique_ptr<GeomRTree> m_boundingBoxTree;
    MovementSharedPtr m_movement = nullptr;
};
typedef std::shared_ptr<MeshGraph> MeshGraphSharedPtr;
typedef LibUtilities::NekFactory<std::string, MeshGraph> MeshGraphFactory;
 
SPATIAL_DOMAINS_EXPORT MeshGraphFactory &GetMeshGraphFactory();
 
/**
 *
 */
inline void MeshGraph::WriteGeometry(
    std::string &outfilename, bool defaultExp,
    const LibUtilities::FieldMetaDataMap &metadata)
{
    v_WriteGeometry(outfilename, defaultExp, metadata);
}
 
/**
 *
 */
inline void MeshGraph::ReadGeometry(LibUtilities::DomainRangeShPtr rng,
                                    bool fillGraph)
{
    v_ReadGeometry(rng, fillGraph);
}
 
/**
 *
 */
inline void MeshGraph::PartitionMesh(
    LibUtilities::SessionReaderSharedPtr session)
{
    v_PartitionMesh(session);
}
 
/**
 *
 */
void MeshGraph::SetExpansionInfo(const std::string variable,
                                 ExpansionInfoMapShPtr &exp)
{
    if (m_expansionMapShPtrMap.count(variable) != 0)
    {
        ASSERTL0(
            false,
            (std::string("ExpansionInfo field is already set for variable ") +
             variable)
                .c_str());
    }
    else
    {
        m_expansionMapShPtrMap[variable] = exp;
    }
}
 
/**
 *
 */
inline bool MeshGraph::SameExpansionInfo(const std::string var1,
                                         const std::string var2)
{
    ExpansionInfoMapShPtr expVec1 = m_expansionMapShPtrMap.find(var1)->second;
    ExpansionInfoMapShPtr expVec2 = m_expansionMapShPtrMap.find(var2)->second;
 
    if (expVec1.get() == expVec2.get())
    {
        return true;
    }
 
    return false;
}
 
/**
 *
 */
inline bool MeshGraph::ExpansionInfoDefined(const std::string var)
{
    return m_expansionMapShPtrMap.count(var);
}
 
} // namespace SpatialDomains
} // namespace Nektar
#endif