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>
44 #include <dune/grid/common/geometry.hh>
46 #include <dune/geometry/type.hh>
48 #include <opm/grid/cpgrid/EntityRep.hpp>
49 #include <opm/grid/utility/platform_dependent/reenable_warnings.h>
51 #include <opm/grid/utility/ErrorMacros.hpp>
66 template <
int mydim,
int cdim>
78 static_assert(cdim == 3,
"");
81 enum { dimension = 3 };
83 enum { mydimension = 0};
85 enum { coorddimension = cdim };
87 enum { dimensionworld = 3 };
98 typedef FieldMatrix< ctype, coorddimension, mydimension >
Jacobian;
140 return Dune::GeometryTypes::cube(mydimension);
152 static_cast<void>(cor);
170 FieldMatrix<ctype, mydimension, coorddimension>
175 return FieldMatrix<ctype, mydimension, coorddimension>();
179 FieldMatrix<ctype, coorddimension, mydimension>
183 return FieldMatrix<ctype, coorddimension, mydimension>();
193 GlobalCoordinate pos_;
203 static_assert(cdim == 3,
"");
206 enum { dimension = 3 };
208 enum { mydimension = 3 };
210 enum { coorddimension = cdim };
212 enum { dimensionworld = 3 };
223 typedef FieldMatrix< ctype, coorddimension, mydimension >
Jacobian;
229 typedef Dune::Impl::FieldMatrixHelper< double > MatrixHelperType;
242 const int* corner_indices)
243 : pos_(pos), vol_(vol), allcorners_(allcorners.data()), cor_idx_(corner_indices)
245 assert(allcorners_ && corner_indices);
259 : pos_(pos), vol_(vol)
265 : pos_(0.0), vol_(0.0), allcorners_(0), cor_idx_(0)
275 static_assert(mydimension == 3,
"");
276 static_assert(coorddimension == 3,
"");
279 uvw[0] -= local_coord;
281 const int pat[8][3] = { { 0, 0, 0 },
290 for (
int i = 0; i < 8; ++i) {
293 for (
int j = 0; j < 3; ++j) {
294 factor *= uvw[pat[i][j]][j];
296 corner_contrib *= factor;
297 xyz += corner_contrib;
306 static_assert(mydimension == 3,
"");
307 static_assert(coorddimension == 3,
"");
310 const ctype epsilon = 1e-12;
311 auto refElement = Dune::ReferenceElements<ctype, 3>::cube();
319 MatrixHelperType::template xTRightInvA<3, 3>(JT, z, dx );
321 }
while (dx.two_norm2() > epsilon*epsilon);
332 return MatrixHelperType::template sqrtDetAAT<3, 3>(Jt);
339 return Dune::GeometryTypes::cube(mydimension);
352 assert(allcorners_ && cor_idx_);
353 return allcorners_[cor_idx_[cor]].center();
372 const JacobianTransposed
375 static_assert(mydimension == 3,
"");
376 static_assert(coorddimension == 3,
"");
380 uvw[0] -= local_coord;
382 const int pat[8][3] = { { 0, 0, 0 },
391 for (
int i = 0; i < 8; ++i) {
392 for (
int deriv = 0; deriv < 3; ++deriv) {
395 for (
int j = 0; j < 3; ++j) {
396 factor *= (j != deriv) ? uvw[pat[i][j]][j]
397 : (pat[i][j] == 0 ? -1.0 : 1.0);
400 corner_contrib *= factor;
401 Jt[deriv] += corner_contrib;
408 const JacobianInverseTransposed
423 GlobalCoordinate pos_;
438 static_assert(cdim == 3,
"");
441 enum { dimension = 3 };
443 enum { mydimension = 2 };
445 enum { coorddimension = cdim };
447 enum { dimensionworld = 3 };
458 typedef FieldMatrix< ctype, coorddimension, mydimension >
Jacobian;
469 : pos_(pos), vol_(vol)
475 : pos_(0.0), vol_(0.0)
482 OPM_THROW(std::runtime_error,
"Geometry::global() meaningless on singular geometry.");
488 OPM_THROW(std::runtime_error,
"Geometry::local() meaningless on singular geometry.");
501 return Dune::GeometryTypes::none(mydimension);
533 const FieldMatrix<ctype, mydimension, coorddimension>&
536 OPM_THROW(std::runtime_error,
"Meaningless to call jacobianTransposed() on singular geometries.");
540 const FieldMatrix<ctype, coorddimension, mydimension>&
543 OPM_THROW(std::runtime_error,
"Meaningless to call jacobianInverseTransposed() on singular geometries.");
553 GlobalCoordinate pos_;
563 template<
int mydim,
int cdim >
564 auto referenceElement(
const cpgrid::Geometry<mydim,cdim>& geo) -> decltype(referenceElement<double,mydim>(geo.type()))
566 return referenceElement<double,mydim>(geo.type());
A class design to hold a variable with a value for each entity of the given codimension,...
Definition: EntityRep.hpp:264
FieldVector< ctype, coorddimension > GlobalCoordinate
Range type of.
Definition: Geometry.hpp:95
bool affine() const
The mapping implemented by this geometry is constant, therefore affine.
Definition: Geometry.hpp:187
FieldMatrix< ctype, coorddimension, mydimension > Jacobian
Type of Jacobian matrix.
Definition: Geometry.hpp:98
FieldVector< ctype, mydimension > LocalCoordinate
Domain type of.
Definition: Geometry.hpp:93
const GlobalCoordinate & center() const
Returns the centroid of the geometry.
Definition: Geometry.hpp:164
double integrationElement(const LocalCoordinate &) const
Returns 1 for the vertex geometry.
Definition: Geometry.hpp:132
FieldMatrix< ctype, coorddimension, mydimension > JacobianInverseTransposed
Type of the inverse of the transposed Jacobian matrix.
Definition: Geometry.hpp:102
Geometry(const GlobalCoordinate &pos)
Construct from vertex position.
Definition: Geometry.hpp:107
FieldMatrix< ctype, coorddimension, mydimension > jacobianInverseTransposed(const LocalCoordinate &) const
This method is meaningless for singular geometries.
Definition: Geometry.hpp:180
GeometryType type() const
Using the cube type for vertices.
Definition: Geometry.hpp:138
ctype volume() const
Volume of vertex is arbitrarily set to 1.
Definition: Geometry.hpp:158
FieldMatrix< ctype, mydimension, coorddimension > jacobianTransposed(const LocalCoordinate &) const
This method is meaningless for singular geometries.
Definition: Geometry.hpp:171
FieldMatrix< ctype, mydimension, coorddimension > JacobianTransposed
Type of transposed Jacobian matrix.
Definition: Geometry.hpp:100
GlobalCoordinate corner(int cor) const
Returns the single corner: the vertex itself.
Definition: Geometry.hpp:150
int corners() const
A vertex is defined by a single corner.
Definition: Geometry.hpp:144
Geometry()
Default constructor, giving a non-valid geometry.
Definition: Geometry.hpp:113
LocalCoordinate local(const GlobalCoordinate &) const
Meaningless for the vertex geometry.
Definition: Geometry.hpp:125
const GlobalCoordinate & global(const LocalCoordinate &) const
Returns the position of the vertex.
Definition: Geometry.hpp:119
double ctype
Coordinate element type.
Definition: Geometry.hpp:90
const GlobalCoordinate & global(const LocalCoordinate &) const
This method is meaningless for singular geometries.
Definition: Geometry.hpp:480
int corners() const
The number of corners of this convex polytope.
Definition: Geometry.hpp:506
LocalCoordinate local(const GlobalCoordinate &) const
This method is meaningless for singular geometries.
Definition: Geometry.hpp:486
FieldMatrix< ctype, coorddimension, mydimension > Jacobian
Type of Jacobian matrix.
Definition: Geometry.hpp:458
bool affine() const
Since integrationElement() is constant, returns true.
Definition: Geometry.hpp:547
const FieldMatrix< ctype, coorddimension, mydimension > & jacobianInverseTransposed(const LocalCoordinate &) const
This method is meaningless for singular geometries.
Definition: Geometry.hpp:541
FieldVector< ctype, mydimension > LocalCoordinate
Domain type of.
Definition: Geometry.hpp:453
FieldMatrix< ctype, coorddimension, mydimension > JacobianInverseTransposed
Type of the inverse of the transposed Jacobian matrix.
Definition: Geometry.hpp:462
ctype volume() const
Volume (area, actually) of intersection.
Definition: Geometry.hpp:521
GeometryType type() const
We use the singular type (None) for intersections.
Definition: Geometry.hpp:499
FieldMatrix< ctype, mydimension, coorddimension > JacobianTransposed
Type of transposed Jacobian matrix.
Definition: Geometry.hpp:460
double integrationElement(const LocalCoordinate &) const
For the singular geometry, we return a constant integration element equal to the volume.
Definition: Geometry.hpp:493
Geometry()
Default constructor, giving a non-valid geometry.
Definition: Geometry.hpp:474
const FieldMatrix< ctype, mydimension, coorddimension > & jacobianTransposed(const LocalCoordinate &) const
This method is meaningless for singular geometries.
Definition: Geometry.hpp:534
Geometry(const GlobalCoordinate &pos, ctype vol)
Construct from centroid and volume (1- and 0-moments).
Definition: Geometry.hpp:467
double ctype
Coordinate element type.
Definition: Geometry.hpp:450
const GlobalCoordinate & center() const
Returns the centroid of the geometry.
Definition: Geometry.hpp:527
FieldVector< ctype, coorddimension > GlobalCoordinate
Range type of.
Definition: Geometry.hpp:455
GlobalCoordinate corner(int) const
This method is meaningless for singular geometries.
Definition: Geometry.hpp:512
const JacobianInverseTransposed jacobianInverseTransposed(const LocalCoordinate &local_coord) const
Inverse of Jacobian transposed.
Definition: Geometry.hpp:409
FieldMatrix< ctype, coorddimension, mydimension > JacobianInverseTransposed
Type of the inverse of the transposed Jacobian matrix.
Definition: Geometry.hpp:227
bool affine() const
The mapping implemented by this geometry is not generally affine.
Definition: Geometry.hpp:417
double ctype
Coordinate element type.
Definition: Geometry.hpp:215
Geometry(const GlobalCoordinate &pos, ctype vol, const EntityVariable< cpgrid::Geometry< 0, 3 >, 3 > &allcorners, const int *corner_indices)
Construct from centroid, volume (1- and 0-moments) and corners.
Definition: Geometry.hpp:239
FieldVector< ctype, coorddimension > GlobalCoordinate
Range type of.
Definition: Geometry.hpp:220
GlobalCoordinate corner(int cor) const
The 8 corners of the hexahedral base cell.
Definition: Geometry.hpp:350
FieldMatrix< ctype, mydimension, coorddimension > JacobianTransposed
Type of transposed Jacobian matrix.
Definition: Geometry.hpp:225
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:337
double integrationElement(const LocalCoordinate &local_coord) const
Equal to \sqrt{\det{J^T J}} where J is the Jacobian.
Definition: Geometry.hpp:329
const JacobianTransposed jacobianTransposed(const LocalCoordinate &local_coord) const
Jacobian transposed.
Definition: Geometry.hpp:373
GlobalCoordinate global(const LocalCoordinate &local_coord) const
Provide a trilinear mapping.
Definition: Geometry.hpp:273
Geometry()
Default constructor, giving a non-valid geometry.
Definition: Geometry.hpp:264
LocalCoordinate local(const GlobalCoordinate &y) const
Mapping from the cell to the reference domain.
Definition: Geometry.hpp:304
const GlobalCoordinate & center() const
Returns the centroid of the geometry.
Definition: Geometry.hpp:363
int corners() const
The number of corners of this convex polytope.
Definition: Geometry.hpp:344
ctype volume() const
Cell volume.
Definition: Geometry.hpp:357
Geometry(const GlobalCoordinate &pos, ctype vol)
Construct from centroid and volume (1- and 0-moments).
Definition: Geometry.hpp:257
FieldVector< ctype, mydimension > LocalCoordinate
Domain type of.
Definition: Geometry.hpp:218
FieldMatrix< ctype, coorddimension, mydimension > Jacobian
Type of Jacobian matrix.
Definition: Geometry.hpp:223
This class encapsulates geometry for both vertices, intersections and cells.
Definition: Geometry.hpp:68
Copyright 2019 Equinor AS.
Definition: CartesianIndexMapper.hpp:10
T volume(const Point< T, 3 > *c)
Computes the volume of a 3D simplex (embedded i 3D space).
Definition: Volumes.hpp:137