36 #ifndef OPM_GEOMETRY_HEADER 37 #define OPM_GEOMETRY_HEADER 40 #include <opm/grid/utility/platform_dependent/disable_warnings.h> 42 #include <dune/common/version.hh> 43 #include <dune/geometry/referenceelements.hh> 45 #if DUNE_VERSION_NEWER(DUNE_GEOMETRY, 2, 5 ) 46 #include <dune/geometry/type.hh> 48 #include <dune/geometry/genericgeometry/geometrytraits.hh> 49 #include <dune/geometry/genericgeometry/matrixhelper.hh> 52 #include <opm/grid/utility/platform_dependent/reenable_warnings.h> 54 #include <opm/grid/utility/ErrorMacros.hpp> 69 template <
int mydim,
int cdim>
81 static_assert(cdim == 3,
"");
84 enum { dimension = 3 };
86 enum { mydimension = 3 };
88 enum { coorddimension = cdim };
90 enum { dimensionworld = 3 };
101 typedef FieldMatrix< ctype, coorddimension, mydimension >
Jacobian;
107 #if DUNE_VERSION_NEWER(DUNE_GRID,2,5) 108 typedef Dune::Impl::FieldMatrixHelper< double > MatrixHelperType;
110 typedef Dune::GenericGeometry::MatrixHelper< Dune::GenericGeometry::DuneCoordTraits<double> > MatrixHelperType;
124 const int* corner_indices)
125 : pos_(pos), vol_(vol), allcorners_(allcorners), cor_idx_(corner_indices)
127 assert(allcorners && corner_indices);
141 : pos_(pos), vol_(vol)
147 : pos_(0.0), vol_(0.0), allcorners_(0), cor_idx_(0)
157 static_assert(mydimension == 3,
"");
158 static_assert(coorddimension == 3,
"");
161 uvw[0] -= local_coord;
163 const int pat[8][3] = { { 0, 0, 0 },
172 for (
int i = 0; i < 8; ++i) {
175 for (
int j = 0; j < 3; ++j) {
176 factor *= uvw[pat[i][j]][j];
178 corner_contrib *= factor;
179 xyz += corner_contrib;
188 static_assert(mydimension == 3,
"");
189 static_assert(coorddimension == 3,
"");
192 const ctype epsilon = 1e-12;
193 const ReferenceElement< ctype , 3 > & refElement =
194 ReferenceElements< ctype, 3 >::general(type());
202 MatrixHelperType::template xTRightInvA<3, 3>(JT, z, dx );
204 }
while (dx.two_norm2() > epsilon*epsilon);
215 return MatrixHelperType::template sqrtDetAAT<3, 3>(Jt);
223 t.makeCube(mydimension);
237 assert(allcorners_ && cor_idx_);
238 return allcorners_[cor_idx_[cor]];
257 const JacobianTransposed
260 static_assert(mydimension == 3,
"");
261 static_assert(coorddimension == 3,
"");
265 uvw[0] -= local_coord;
267 const int pat[8][3] = { { 0, 0, 0 },
276 for (
int i = 0; i < 8; ++i) {
277 for (
int deriv = 0; deriv < 3; ++deriv) {
280 for (
int j = 0; j < 3; ++j) {
281 factor *= (j != deriv) ? uvw[pat[i][j]][j]
282 : (pat[i][j] == 0 ? -1.0 : 1.0);
285 corner_contrib *= factor;
286 Jt[deriv] += corner_contrib;
293 const JacobianInverseTransposed
308 GlobalCoordinate pos_;
310 const GlobalCoordinate* allcorners_;
323 static_assert(cdim == 3,
"");
326 enum { dimension = 3 };
328 enum { mydimension = 2 };
330 enum { coorddimension = cdim };
332 enum { dimensionworld = 3 };
343 typedef FieldMatrix< ctype, coorddimension, mydimension >
Jacobian;
354 : pos_(pos), vol_(vol)
360 : pos_(0.0), vol_(0.0)
367 OPM_THROW(std::runtime_error,
"Geometry::global() meaningless on singular geometry.");
373 OPM_THROW(std::runtime_error,
"Geometry::local() meaningless on singular geometry.");
387 t.makeNone(mydimension);
420 const FieldMatrix<ctype, mydimension, coorddimension>&
423 OPM_THROW(std::runtime_error,
"Meaningless to call jacobianTransposed() on singular geometries.");
427 const FieldMatrix<ctype, coorddimension, mydimension>&
430 OPM_THROW(std::runtime_error,
"Meaningless to call jacobianInverseTransposed() on singular geometries.");
440 GlobalCoordinate pos_;
451 static_assert(cdim == 3,
"");
454 enum { dimension = 3 };
456 enum { mydimension = 0};
458 enum { coorddimension = cdim };
460 enum { dimensionworld = 3 };
471 typedef FieldMatrix< ctype, coorddimension, mydimension >
Jacobian;
514 t.makeCube(mydimension);
527 static_cast<void>(cor);
545 FieldMatrix<ctype, mydimension, coorddimension>
550 return FieldMatrix<ctype, mydimension, coorddimension>();
554 FieldMatrix<ctype, coorddimension, mydimension>
558 return FieldMatrix<ctype, coorddimension, mydimension>();
568 GlobalCoordinate pos_;
578 #endif // OPM_GEOMETRY_HEADER const GlobalCoordinate & global(const LocalCoordinate &) const
Returns the position of the vertex.
Definition: Geometry.hpp:492
FieldMatrix< ctype, coorddimension, mydimension > Jacobian
Type of Jacobian matrix.
Definition: Geometry.hpp:343
GeometryType type() const
Using the cube type for all entities now (cells and vertices), but we use the singular type for inter...
Definition: Geometry.hpp:220
ctype volume() const
Volume of vertex is arbitrarily set to 1.
Definition: Geometry.hpp:533
Geometry(const GlobalCoordinate &pos, ctype vol)
Construct from centroid and volume (1- and 0-moments).
Definition: Geometry.hpp:139
int corners() const
The number of corners of this convex polytope.
Definition: Geometry.hpp:229
Geometry(const GlobalCoordinate &pos, ctype vol, const GlobalCoordinate *allcorners, const int *corner_indices)
Construct from centroid, volume (1- and 0-moments) and corners.
Definition: Geometry.hpp:121
double ctype
Coordinate element type.
Definition: Geometry.hpp:463
LocalCoordinate local(const GlobalCoordinate &) const
This method is meaningless for singular geometries.
Definition: Geometry.hpp:371
FieldMatrix< ctype, coorddimension, mydimension > JacobianInverseTransposed
Type of the inverse of the transposed Jacobian matrix.
Definition: Geometry.hpp:347
Geometry(const GlobalCoordinate &pos, ctype vol)
Construct from centroid and volume (1- and 0-moments).
Definition: Geometry.hpp:352
Geometry()
Default constructor, giving a non-valid geometry.
Definition: Geometry.hpp:146
FieldMatrix< ctype, coorddimension, mydimension > Jacobian
Type of Jacobian matrix.
Definition: Geometry.hpp:471
const JacobianTransposed jacobianTransposed(const LocalCoordinate &local_coord) const
Jacobian transposed.
Definition: Geometry.hpp:258
FieldMatrix< ctype, coorddimension, mydimension > JacobianInverseTransposed
Type of the inverse of the transposed Jacobian matrix.
Definition: Geometry.hpp:475
Holds the implementation of the CpGrid as a pimple.
Definition: OpmParserIncludes.hpp:42
Geometry()
Default constructor, giving a non-valid geometry.
Definition: Geometry.hpp:486
FieldVector< ctype, coorddimension > GlobalCoordinate
Range type of.
Definition: Geometry.hpp:98
double ctype
Coordinate element type.
Definition: Geometry.hpp:93
This class encapsulates geometry for both vertices, intersections and cells.
Definition: CpGridData.hpp:92
double integrationElement(const LocalCoordinate &) const
Returns 1 for the vertex geometry.
Definition: Geometry.hpp:505
bool affine() const
The mapping implemented by this geometry is constant, therefore affine.
Definition: Geometry.hpp:562
GlobalCoordinate global(const LocalCoordinate &local_coord) const
Provide a trilinear mapping.
Definition: Geometry.hpp:155
double integrationElement(const LocalCoordinate &local_coord) const
Equal to {{J^T J}} where J is the Jacobian.
Definition: Geometry.hpp:212
ctype volume() const
Volume (area, actually) of intersection.
Definition: Geometry.hpp:408
int corners() const
A vertex is defined by a single corner.
Definition: Geometry.hpp:519
const FieldMatrix< ctype, mydimension, coorddimension > & jacobianTransposed(const LocalCoordinate &) const
This method is meaningless for singular geometries.
Definition: Geometry.hpp:421
FieldVector< ctype, mydimension > LocalCoordinate
Domain type of.
Definition: Geometry.hpp:466
const GlobalCoordinate & center() const
Returns the centroid of the geometry.
Definition: Geometry.hpp:248
Geometry()
Default constructor, giving a non-valid geometry.
Definition: Geometry.hpp:359
const JacobianInverseTransposed jacobianInverseTransposed(const LocalCoordinate &local_coord) const
Inverse of Jacobian transposed.
Definition: Geometry.hpp:294
int corners() const
The number of corners of this convex polytope.
Definition: Geometry.hpp:393
const FieldMatrix< ctype, coorddimension, mydimension > & jacobianInverseTransposed(const LocalCoordinate &) const
This method is meaningless for singular geometries.
Definition: Geometry.hpp:428
double integrationElement(const LocalCoordinate &) const
For the singular geometry, we return a constant integration element equal to the volume.
Definition: Geometry.hpp:378
FieldMatrix< ctype, mydimension, coorddimension > jacobianTransposed(const LocalCoordinate &) const
This method is meaningless for singular geometries.
Definition: Geometry.hpp:546
const GlobalCoordinate & center() const
Returns the centroid of the geometry.
Definition: Geometry.hpp:539
LocalCoordinate local(const GlobalCoordinate &y) const
Mapping from the cell to the reference domain.
Definition: Geometry.hpp:186
bool affine() const
The mapping implemented by this geometry is not generally affine.
Definition: Geometry.hpp:302
GeometryType type() const
Using the cube type for vertices.
Definition: Geometry.hpp:511
FieldVector< ctype, coorddimension > GlobalCoordinate
Range type of.
Definition: Geometry.hpp:468
GlobalCoordinate corner(int cor) const
Returns the single corner: the vertex itself.
Definition: Geometry.hpp:525
FieldVector< ctype, mydimension > LocalCoordinate
Domain type of.
Definition: Geometry.hpp:338
LocalCoordinate local(const GlobalCoordinate &) const
Meaningless for the vertex geometry.
Definition: Geometry.hpp:498
FieldMatrix< ctype, mydimension, coorddimension > JacobianTransposed
Type of transposed Jacobian matrix.
Definition: Geometry.hpp:345
FieldMatrix< ctype, coorddimension, mydimension > JacobianInverseTransposed
Type of the inverse of the transposed Jacobian matrix.
Definition: Geometry.hpp:105
const GlobalCoordinate & center() const
Returns the centroid of the geometry.
Definition: Geometry.hpp:414
FieldMatrix< ctype, coorddimension, mydimension > jacobianInverseTransposed(const LocalCoordinate &) const
This method is meaningless for singular geometries.
Definition: Geometry.hpp:555
ctype volume() const
Cell volume.
Definition: Geometry.hpp:242
Geometry(const GlobalCoordinate &pos)
Construct from vertex position.
Definition: Geometry.hpp:480
T volume(const Point< T, 3 > *c)
Computes the volume of a 3D simplex (embedded i 3D space).
Definition: Volumes.hpp:137
const GlobalCoordinate & global(const LocalCoordinate &) const
This method is meaningless for singular geometries.
Definition: Geometry.hpp:365
bool affine() const
Since integrationElement() is constant, returns true.
Definition: Geometry.hpp:434
GlobalCoordinate corner(int) const
This method is meaningless for singular geometries.
Definition: Geometry.hpp:399
FieldMatrix< ctype, mydimension, coorddimension > JacobianTransposed
Type of transposed Jacobian matrix.
Definition: Geometry.hpp:473
FieldVector< ctype, mydimension > LocalCoordinate
Domain type of.
Definition: Geometry.hpp:96
GeometryType type() const
We use the singular type (None) for intersections.
Definition: Geometry.hpp:384
FieldMatrix< ctype, mydimension, coorddimension > JacobianTransposed
Type of transposed Jacobian matrix.
Definition: Geometry.hpp:103
FieldVector< ctype, coorddimension > GlobalCoordinate
Range type of.
Definition: Geometry.hpp:340
FieldMatrix< ctype, coorddimension, mydimension > Jacobian
Type of Jacobian matrix.
Definition: Geometry.hpp:101
double ctype
Coordinate element type.
Definition: Geometry.hpp:335
GlobalCoordinate corner(int cor) const
The 8 corners of the hexahedral base cell.
Definition: Geometry.hpp:235