mesh.h 7.35 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
16
17
#include "metis.h"
#include "parmetis.h"
#include "toolboxes/errormessages.h"
#include "toolboxes/reconstructor.h"

18
19
20
class Mesh
{
  protected:
Steffen Schotthöfer's avatar
Steffen Schotthöfer committed
21
22
23
24
25
26
27
    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 */
    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 */
28
    unsigned _numNodesPerBoundary;
Steffen Schotthöfer's avatar
Steffen Schotthöfer committed
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
    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*/
    void ComputeConnectivity();  /*! @brief: Computes _cellNeighbors and _nodeNeighbors, i.e. neighborship relation in mesh*/
    void ComputePartitioning();  /*! @brief: Computes local partitioning for openMP */

    /*! @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 );
    void ComputeBounds(); /*! @brief: Computes the spatial bounds of a 2D domain. */
61
62

  public:
Steffen Schotthöfer's avatar
Steffen Schotthöfer committed
63
64
65
66
    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*/
67
68
    Mesh( std::vector<Vector> nodes,
          std::vector<std::vector<unsigned>> cells,
Steffen Schotthöfer's avatar
Steffen Schotthöfer committed
69
          std::vector<std::pair<BOUNDARY_TYPE, std::vector<unsigned>>> boundaries );
70
71
    ~Mesh();

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

Steffen Schotthöfer's avatar
Steffen Schotthöfer committed
77
78
    /*! @brief: Returns all node coordinates
     *  @return: dimension: numNodes x dim */
79
    const std::vector<Vector>& GetNodes() const;
80

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

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

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

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

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

Steffen Schotthöfer's avatar
Steffen Schotthöfer committed
101
102
    /*! @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
103
    const std::vector<std::vector<Vector>>& GetNormals() const;
104

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

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

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

Steffen Schotthöfer's avatar
Steffen Schotthöfer committed
117
118
119
120
121
122
    /*! @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
123
    void ComputeSlopes( unsigned nq, VectorVector& psiDerX, VectorVector& psiDerY, const VectorVector& psi ) const;
vavrines's avatar
vavrines committed
124

Steffen Schotthöfer's avatar
Steffen Schotthöfer committed
125
126
127
128
129
    /*! @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
130
    void ReconstructSlopesS( unsigned nq, VectorVector& psiDerX, VectorVector& psiDerY, const VectorVector& psi ) const;
Steffen Schotthöfer's avatar
Steffen Schotthöfer committed
131
132
133
134
135
    /*! @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
136
    void ReconstructSlopesU( unsigned nq, VectorVector& psiDerX, VectorVector& psiDerY, const VectorVector& psi ) const;
137
138
139
};

#endif    // MESH_H