mesh.h 7.25 KB
Newer Older
1
2
3
4
5
6
/*!
 * @file mesh.h
 * @brief Class for mesh description and partitioning
 * @author J. Wolters
 */

7
8
9
#ifndef MESH_H
#define MESH_H

10
#include "blaze/math/CompressedMatrix.h"
11
12
#include "common/globalconstants.h"
#include "common/typedef.h"
13

14
#include <vector>
vavrines's avatar
vavrines committed
15

16
17
18
class Mesh
{
  protected:
19
20
21
22
23
    const unsigned _dim;             /*! @brief: spatial dimension of the mesh, i.e. 1D,2D,3D */
    const unsigned _numCells;        /*! @brief: number of cells in the mesh */
    const unsigned _numNodes;        /*! @brief: number of nodes in the mesh (for node centered view)*/
    const unsigned _numNodesPerCell; /*! @brief: number of nodes per cell */
    const unsigned _numBoundaries;   /*! @brief: number of boundary cells in the mesh */
24
25
    const unsigned _ghostCellID; /*! @brief: Id of the ghost cell. (we use only one ghost cell). equal to _numCells and therefore has the ID of the
                                    last cell + 1 */
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47

    std::vector<std::pair<double, double>> _bounds;    // ???

    std::vector<Vector> _nodes;                /*! @brief: nodes in the mesh. dimension:_numNodes<_dim> */
    std::vector<std::vector<unsigned>> _cells; /*! @brief: cells in the mesh. dimension:_numCells<_numNodesPerCell>  */

    /*! @brief: boundary cells in the mesh. Pair defines boundary type of the boundary nodes of the cell. numBoundaries<(1,numBoundaryNodes)>*/
    std::vector<std::pair<BOUNDARY_TYPE, std::vector<unsigned>>> _boundaries;
    std::vector<double> _cellAreas;                    /*! @brief: cell areas of the mesh. dimension: numCells*/
    std::vector<Vector> _cellMidPoints;                /*! @brief: cell midpoints of the mesh. dimension: numCells<dim>*/
    std::vector<std::vector<unsigned>> _cellNeighbors; /*! @brief: neighbors of each cell. dimension: numCells<numNodesPerCell>*/

    /*! @brief: outward facing normals of each side of each cell. dimension: numCells<numNodesPerCell<dim>>, all
                normals are facing away from the cell center, and scaled with the edge length */
    std::vector<std::vector<Vector>> _cellNormals;
    /*! @brief: Tags each cell with its boundary type. None means no boundary. dimension: numCells */
    std::vector<BOUNDARY_TYPE> _cellBoundaryTypes;
    std::vector<unsigned> _colors;                /*! @brief: Color of each cell (for MPI mesh partitioning). dimension: numCells */
    blaze::CompressedMatrix<bool> _nodeNeighbors; /*! @brief: neighborshood relationship of nodes for (par-)metis */

    void ComputeCellAreas();     /*! @brief: Computes only the areas of the mesh cells. Write to _cellAreas. */
    void ComputeCellMidpoints(); /*! @brief: Compute only the midpoints of the cells. Write to _cellMidPoints*/
48
49
    void ComputeConnectivity();  /*! @brief: Computes _cellNeighbors and _nodeNeighbors, i.e. neighborship relation in mesh*/
    void ComputePartitioning();  /*! @brief: Computes local partitioning for openMP */
50
51
52
53
54
55
56
57

    /*! @brief: Computes outward facing normal of two neighboring nodes nodeA and nodeB with common cellCellcenter.
     *          Normals are scaled with their respective edge length
     *  @param: nodeA: first node
     *  @param: nodeB: neighboring node to nodeA
     *  @param: cellCenter: Center of the cell that has nodeA and nodeB as nodes.
     *  @return: outward facing normal */
    Vector ComputeOutwardFacingNormal( const Vector& nodeA, const Vector& nodeB, const Vector& cellCenter );
58
    void ComputeBounds(); /*! @brief: Computes the spatial bounds of a 2D domain. */
59
60

  public:
61
62
63
64
    Mesh() = delete;    //  no default constructor

    /*! @brief: Constructor of mesh. Needs nodes, cells, and boundary descriptions as specified above.
     *          See LoadSU2MeshFromFile in io.cpp for setup information*/
65
66
    Mesh( std::vector<Vector> nodes,
          std::vector<std::vector<unsigned>> cells,
67
          std::vector<std::pair<BOUNDARY_TYPE, std::vector<unsigned>>> boundaries );
68
69
    ~Mesh();

70
71
72
73
74
    inline unsigned GetDim() const { return _dim; }
    inline unsigned GetNumCells() const { return _numCells; }
    inline unsigned GetNumNodes() const { return _numNodes; }
    inline unsigned GetNumNodesPerCell() const { return _numNodesPerCell; }

75
76
    /*! @brief: Returns all node coordinates
     *  @return: dimension: numNodes x dim */
77
    const std::vector<Vector>& GetNodes() const;
78

79
80
    /*! @brief  Returns the mid point coordinates of each cell
     *  @return dimension: numCells x dim */
Jonas Kusch's avatar
Jonas Kusch committed
81
    const std::vector<Vector>& GetCellMidPoints() const;
82

83
84
    /*! @brief Returns all node IDs that construct up each cell
     *  @return dimension: numCells x numNodes */
85
    const std::vector<std::vector<unsigned>>& GetCells() const;
86

87
88
    /*! @brief Returns the cell area of each cell
     *  @return dimension: numCells */
89
    const std::vector<double>& GetCellAreas() const;
90

91
92
    /*! @brief Return the color/ID of the mesh partition
     *  @return dimension: numCells */
Jannick Wolters's avatar
Jannick Wolters committed
93
    const std::vector<unsigned>& GetPartitionIDs() const;
94

95
96
    /*! @brief Returns the neighbor cell IDs for every cell
     *  @return dimension: numCells x numNodes */
Jonas Kusch's avatar
Jonas Kusch committed
97
    const std::vector<std::vector<unsigned>>& GetNeighbours() const;
98

99
100
    /*! @brief Returns the edge length scaled normal vectors of each cell
     *  @return dimension: numCells x numNodes x dim */
Jonas Kusch's avatar
Jonas Kusch committed
101
    const std::vector<std::vector<Vector>>& GetNormals() const;
102

103
104
    /*! @brief Returns the boundary enum for each cell. BOUNDARY_TYPE::NONE is the default.
     *  @return dimension: numCells */
105
    const std::vector<BOUNDARY_TYPE>& GetBoundaryTypes() const;
106

107
108
    /*! @brief Returns the minimal and maximal coordinates of all nodes for each dimension
     *  @return dimension: dim */
109
110
    const std::vector<std::pair<double, double>> GetBounds() const;

111
112
    /*! @brief Returns distance of a specified cells center to the coordinate systems origin
     *  @return dimension: scalar */
steffen.schotthoefer's avatar
steffen.schotthoefer committed
113
    double GetDistanceToOrigin( unsigned idx_cell ) const;
114

115
    /*! @brief ComputeSlopes calculates the slope in every cell into x and y direction using the divergence theorem.
116
117
118
119
120
     *  @param nq is number of quadrature points
     *  @param psiDerX is slope in x direction (gets computed. Slope is stored here)
     *  @param psiDerY is slope in y direction (gets computed. Slope is stored here)
     *  @param psi is solution for which slope is computed */
    // Not used
121
    void ComputeSlopes( unsigned nq, VectorVector& psiDerX, VectorVector& psiDerY, const VectorVector& psi ) const;
vavrines's avatar
vavrines committed
122

Steffen Schotthöfer's avatar
Steffen Schotthöfer committed
123
124
125
126
127
    /*! @brief:Structured mesh slope reconstruction with flux limiters.
     *  @param nq is number of quadrature points
     *  @param psiDerX is slope in x direction (gets computed. Slope is stored here)
     *  @param psiDerY is slope in y direction (gets computed. Slope is stored here)
     *  @param psi is solution for which slope is computed */
vavrines's avatar
vavrines committed
128
    void ReconstructSlopesS( unsigned nq, VectorVector& psiDerX, VectorVector& psiDerY, const VectorVector& psi ) const;
Steffen Schotthöfer's avatar
Steffen Schotthöfer committed
129
130
131
132
133
    /*! @brief: Use gauss theorem and limiters. For unstructured mesh *
     *  @param nq is number of quadrature points
     *  @param psiDerX is slope in x direction (gets computed. Slope is stored here)
     *  @param psiDerY is slope in y direction (gets computed. Slope is stored here)
     *  @param psi is solution for which slope is computed */
vavrines's avatar
vavrines committed
134
    void ReconstructSlopesU( unsigned nq, VectorVector& psiDerX, VectorVector& psiDerY, const VectorVector& psi ) const;
135
136
137
};

#endif    // MESH_H