mesh.h 6.93 KB
Newer Older
1
2
3
#ifndef MESH_H
#define MESH_H

4
#include "blaze/math/CompressedMatrix.h"
5
6
#include "common/globalconstants.h"
#include "common/typedef.h"
7

8
9
10
#include <algorithm>
#include <mpi.h>
#include <omp.h>
11
#include <vector>
vavrines's avatar
vavrines committed
12

13
14
15
#include "toolboxes/errormessages.h"
#include "toolboxes/reconstructor.h"

16
17
18
class Mesh
{
  protected:
Steffen Schotthöfer's avatar
Steffen Schotthöfer committed
19
20
21
    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)*/
22
    const unsigned _numNodesPerCell; /*!< @brief number of nodes per cell */
Steffen Schotthöfer's avatar
Steffen Schotthöfer committed
23
    const unsigned _numBoundaries;   /*!< @brief number of boundary cells in the mesh */
24
    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
Steffen Schotthöfer's avatar
Steffen Schotthöfer committed
25
                                    last cell + 1 */
26
    unsigned _numNodesPerBoundary;
Steffen Schotthöfer's avatar
Steffen Schotthöfer committed
27
28
    std::vector<std::pair<double, double>> _bounds;    // ???

Steffen Schotthöfer's avatar
Steffen Schotthöfer committed
29
    std::vector<Vector> _nodes;                /*!< @brief nodes in the mesh. dimension:_numNodes<_dim> */
30
    std::vector<std::vector<unsigned>> _cells; /*!< @brief cells in the mesh. dimension:_numCells<_numNodesPerCell>  */
Steffen Schotthöfer's avatar
Steffen Schotthöfer committed
31

32
    /*! @brief boundary cells in the mesh. Pair defines boundary type of the boundary nodes of the cell. numBoundaries<(1,numBoundaryNodes)>*/
Steffen Schotthöfer's avatar
Steffen Schotthöfer committed
33
    std::vector<std::pair<BOUNDARY_TYPE, std::vector<unsigned>>> _boundaries;
Steffen Schotthöfer's avatar
Steffen Schotthöfer committed
34
35
    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>*/
36
    std::vector<std::vector<unsigned>> _cellNeighbors; /*!< @brief neighbors of each cell. dimension: numCells<numNodesPerCell>*/
Steffen Schotthöfer's avatar
Steffen Schotthöfer committed
37

38
    /*! @brief outward facing normals of each side of each cell. dimension: numCells<numNodesPerCell<dim>>, all
Steffen Schotthöfer's avatar
Steffen Schotthöfer committed
39
40
                normals are facing away from the cell center, and scaled with the edge length */
    std::vector<std::vector<Vector>> _cellNormals;
41
    /*! @brief Tags each cell with its boundary type. None means no boundary. dimension: numCells */
Steffen Schotthöfer's avatar
Steffen Schotthöfer committed
42
    std::vector<BOUNDARY_TYPE> _cellBoundaryTypes;
43
    blaze::CompressedMatrix<bool> _nodeNeighbors; /*!< @brief neighborshood relationship of nodes for (par-)metis */
Steffen Schotthöfer's avatar
Steffen Schotthöfer committed
44

Steffen Schotthöfer's avatar
Steffen Schotthöfer committed
45
    void ComputeCellAreas();     /*!< @brief Computes only the areas of the mesh cells. Write to _cellAreas. */
46
    void ComputeCellMidpoints(); /*!< @brief Compute only the midpoints of the cells. Write to _cellMidPoints*/
Steffen Schotthöfer's avatar
Steffen Schotthöfer committed
47
    void ComputeConnectivity();  /*!< @brief Computes _cellNeighbors and _nodeNeighbors, i.e. neighborship relation in mesh*/
Steffen Schotthöfer's avatar
Steffen Schotthöfer committed
48

49
    /*! @brief Computes outward facing normal of two neighboring nodes nodeA and nodeB with common cellCellcenter.
Steffen Schotthöfer's avatar
Steffen Schotthöfer committed
50
     *          Normals are scaled with their respective edge length
51
52
53
54
     *  @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 */
Steffen Schotthöfer's avatar
Steffen Schotthöfer committed
55
    Vector ComputeOutwardFacingNormal( const Vector& nodeA, const Vector& nodeB, const Vector& cellCenter );
56
    void ComputeBounds(); /*!< @brief Computes the spatial bounds of a 2D domain. */
57
58

  public:
Steffen Schotthöfer's avatar
Steffen Schotthöfer committed
59
60
    Mesh() = delete;    //  no default constructor

61
    /*! @brief Constructor of mesh. Needs nodes, cells, and boundary descriptions as specified above.
Steffen Schotthöfer's avatar
Steffen Schotthöfer committed
62
     *          See LoadSU2MeshFromFile in io.cpp for setup information*/
63
64
    Mesh( std::vector<Vector> nodes,
          std::vector<std::vector<unsigned>> cells,
Steffen Schotthöfer's avatar
Steffen Schotthöfer committed
65
          std::vector<std::pair<BOUNDARY_TYPE, std::vector<unsigned>>> boundaries );
66
67
    ~Mesh();

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

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

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

Steffen Schotthöfer's avatar
Steffen Schotthöfer committed
81
82
    /*! @brief Returns all node IDs that construct up each cell
     *  @return dimension: numCells x numNodes */
83
    const std::vector<std::vector<unsigned>>& GetCells() const;
84

Steffen Schotthöfer's avatar
Steffen Schotthöfer committed
85
86
    /*! @brief Returns the cell area of each cell
     *  @return dimension: numCells */
87
    const std::vector<double>& GetCellAreas() const;
88

Steffen Schotthöfer's avatar
Steffen Schotthöfer committed
89
90
    /*! @brief Returns the neighbor cell IDs for every cell
     *  @return dimension: numCells x numNodes */
Jonas Kusch's avatar
Jonas Kusch committed
91
    const std::vector<std::vector<unsigned>>& GetNeighbours() const;
92

Steffen Schotthöfer's avatar
Steffen Schotthöfer committed
93
94
    /*! @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
95
    const std::vector<std::vector<Vector>>& GetNormals() const;
96

Steffen Schotthöfer's avatar
Steffen Schotthöfer committed
97
98
    /*! @brief Returns the boundary enum for each cell. BOUNDARY_TYPE::NONE is the default.
     *  @return dimension: numCells */
99
    const std::vector<BOUNDARY_TYPE>& GetBoundaryTypes() const;
100

Steffen Schotthöfer's avatar
Steffen Schotthöfer committed
101
102
    /*! @brief Returns the minimal and maximal coordinates of all nodes for each dimension
     *  @return dimension: dim */
103
104
    const std::vector<std::pair<double, double>> GetBounds() const;

Steffen Schotthöfer's avatar
Steffen Schotthöfer committed
105
106
    /*! @brief Returns distance of a specified cells center to the coordinate systems origin
     *  @return dimension: scalar */
steffen.schotthoefer's avatar
steffen.schotthoefer committed
107
    double GetDistanceToOrigin( unsigned idx_cell ) const;
108

Steffen Schotthöfer's avatar
Steffen Schotthöfer committed
109
110
111
112
113
114
    /*! @brief ComputeSlopes calculates the slope in every cell into x and y direction using the divergence theorem.
     *  @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
115
    void ComputeSlopes( unsigned nq, VectorVector& psiDerX, VectorVector& psiDerY, const VectorVector& psi ) const;
vavrines's avatar
vavrines committed
116

117
    /*! @briefStructured mesh slope reconstruction with flux limiters.
Steffen Schotthöfer's avatar
Steffen Schotthöfer committed
118
119
120
121
     *  @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
122
    void ReconstructSlopesS( unsigned nq, VectorVector& psiDerX, VectorVector& psiDerY, const VectorVector& psi ) const;
123
    /*! @brief Use gauss theorem and limiters. For unstructured mesh *
Steffen Schotthöfer's avatar
Steffen Schotthöfer committed
124
125
126
127
     *  @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 ReconstructSlopesU( unsigned nq, VectorVector& psiDerX, VectorVector& psiDerY, const VectorVector& psi ) const;
129
130
131
};

#endif    // MESH_H