My Project
Loading...
Searching...
No Matches
WachspressCoord.hpp
1/*
2 Copyright 2012 SINTEF ICT, Applied Mathematics.
3
4 This file is part of the Open Porous Media project (OPM).
5
6 OPM is free software: you can redistribute it and/or modify
7 it under the terms of the GNU General Public License as published by
8 the Free Software Foundation, either version 3 of the License, or
9 (at your option) any later version.
10
11 OPM is distributed in the hope that it will be useful,
12 but WITHOUT ANY WARRANTY; without even the implied warranty of
13 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 GNU General Public License for more details.
15
16 You should have received a copy of the GNU General Public License
17 along with OPM. If not, see <http://www.gnu.org/licenses/>.
18*/
19
20#ifndef OPM_WACHSPRESSCOORD_HEADER_INCLUDED
21#define OPM_WACHSPRESSCOORD_HEADER_INCLUDED
22
23#include <opm/grid/utility/SparseTable.hpp>
24#include <vector>
25
26struct UnstructuredGrid;
27
28namespace Opm
29{
30
37 {
38 public:
41 explicit WachspressCoord(const UnstructuredGrid& grid);
42
46 int numCorners(const int cell) const;
47
55 void cartToBary(const int cell,
56 const double* x,
57 double* xb) const;
58
59 // A corner is here defined as a {cell, vertex} pair where the
60 // vertex is adjacent to the cell.
62 {
63 int corner_id; // Unique for each corner.
64 int vertex; // Shared between corners belonging to different cells.
65 double volume; // Defined as det(N) where N is the matrix of adjacent face normals.
66 };
67
71
74 const std::vector<int>& adjacentFaces() const;
75
76 private:
77 const UnstructuredGrid& grid_;
78 SparseTable<CornerInfo> corner_info_; // Corner info by cell.
79 std::vector<int> adj_faces_; // Set of adjacent faces, by corner id. Contains dim face indices per corner.
80 SparseTable<int> nonadj_faces_; // Set of nonadjacent faces, by corner id.
81 };
82
83} // namespace Opm
84
85#endif // OPM_WACHSPRESSCOORD_HEADER_INCLUDED
A SparseTable stores a table with rows of varying size as efficiently as possible.
Definition SparseTable.hpp:55
Class capable of computing Wachspress coordinates in 2d and 3d.
Definition WachspressCoord.hpp:37
void cartToBary(const int cell, const double *x, double *xb) const
Compute generalized barycentric coordinates for some point x with respect to the vertices of a grid c...
Definition WachspressCoord.cpp:202
const SparseTable< CornerInfo > & cornerInfo() const
The class stores some info for each corner, made accessible for user convenience.
Definition WachspressCoord.cpp:179
const std::vector< int > & adjacentFaces() const
The class stores adjacent faces for each corner, made accessible for user convenience.
Definition WachspressCoord.cpp:188
int numCorners(const int cell) const
Count of vertices adjacent to a call.
Definition WachspressCoord.cpp:171
Holds the implementation of the CpGrid as a pimple.
Definition CellQuadrature.cpp:68
Definition WachspressCoord.hpp:62
Data structure for an unstructured grid, unstructured meaning that any cell may have an arbitrary num...
Definition UnstructuredGrid.h:99