12 #ifndef BOUNDARYGRID_HH_
13 #define BOUNDARYGRID_HH_
15 #include <opm/common/utility/platform_dependent/disable_warnings.h>
17 #include <dune/common/version.hh>
18 #include <dune/common/fvector.hh>
19 #include <dune/common/fmatrix.hh>
20 #include <dune/geometry/referenceelements.hh>
21 #if ! DUNE_VERSION_NEWER(DUNE_GEOMETRY, 2, 5)
22 #include <dune/geometry/genericgeometry/matrixhelper.hh>
24 #include <dune/grid/common/mcmgmapper.hh>
26 #include <opm/common/utility/platform_dependent/reenable_warnings.h>
32 namespace Elasticity {
51 int k1,
int k2,
bool dc=
false);
84 return hypot(v2.
c[0]-
c[0],v2.
c[1]-
c[1]) < 1.e-8;
93 v[0].
i =
v[1].
i =
v[2].
i =
v[3].
i = 0;
94 v[0].
c =
v[1].
c =
v[2].
c =
v[3].
c = 0.f;
105 std::vector<double>
evalBasis(
double xi,
double eta)
const;
115 os <<
"Nodes: " << q.
v[0].
i <<
"," << q.
v[1].
i <<
"," << q.
v[2].
i <<
"," << q.
v[3].
i << std::endl;
116 os <<
"(" << q.
v[0].
c <<
")(" << q.
v[1].
c <<
")(" << q.
v[2].
c <<
")(" << q.
v[3].
c <<
")";
123 void add(
const Quad& quad);
125 void addToColumn(
size_t col,
const Quad& quad)
127 if (col >= colGrids.size())
128 colGrids.resize(col+1);
129 colGrids[col].push_back(quad);
148 const Quad& getQuad(
int col,
int index)
const
150 return colGrids[col][index];
159 size_t colSize(
int i)
const
161 return colGrids[i].size();
182 bool find(Vertex& res,
const Vertex& coord)
const;
195 static void extract(Vertex& res,
196 const Dune::FieldVector<double,3>& coord,
int dir);
226 return (
coord[0] >= q.
bb[0][0]-eps &&
238 std::vector<std::vector<Quad> > colGrids;
249 for (
size_t i=0;i<g.
size();++i)
250 os << g[i] << std::endl;
264 const double* A,
const double* B,
265 std::vector<double>& X,
266 std::vector<double>& Y)
const;
275 bool cubicSolve(
double eps,
double A,
double B,
double C,
276 double D, std::vector<double>& X)
const;
285 double epsZero,
double epsOut)
const;
289 template<
int mydim,
int cdim,
class Gr
idImp>
295 template<
int cdim,
class Gr
idImp>
300 enum { dimension = 2};
303 enum { mydimension = 2};
306 enum { coorddimension = cdim };
309 enum {dimensionworld = 2 };
321 typedef Dune::FieldMatrix<ctype,coorddimension,mydimension>
Jacobian;
332 Dune::LeafMultipleCodimMultipleGeomTypeMapper<GridImp,
333 Dune::MCMGVertexLayout> mapper(gv);
334 typename GridImp::LeafGridView::template Codim<3>::Iterator start=gv.leafGridView().template begin<3>();
335 const typename GridImp::LeafGridView::template Codim<3>::Iterator itend = gv.leafGridView().template end<3>();
336 for (
int i=0;i<4;++i) {
337 typename GridImp::LeafGridView::template Codim<3>::Iterator it=start;
338 for (; it != itend; ++it) {
339 if (mapper.index(*it) == q.
v[i].
i)
345 m_volume = (c[1][0]-c[0][0])*(c[2][1]-c[0][1]);
352 for (
int i=0;i<4;++i)
354 m_volume = (c[1][0]-c[0][0])*(c[2][1]-c[0][1]);
358 Dune::GeometryType
type()
const
360 Dune::GeometryType t;
361 t.makeCube(mydimension);
382 return Global(local);
400 const int pat[4][2] = {{ 0, 0 },
405 for (
int i = 0; i < 4; ++i) {
408 for (
int j = 0; j < 2; ++j) {
409 factor *= uvw[pat[i][j]][j];
411 corner_contrib *= factor;
412 xyz += corner_contrib;
421 const ctype epsilon = 1e-10;
422 const Dune::ReferenceElement< ctype , 2 > & refElement =
423 Dune::ReferenceElements< ctype, 2 >::general(type());
431 #if DUNE_VERSION_NEWER(DUNE_GEOMETRY, 2, 5)
432 Dune::Impl::FieldMatrixHelper<double>::template xTRightInvA<2, 2>(JT, z, dx );
434 using namespace Dune::GenericGeometry;
435 MatrixHelper<DuneCoordTraits<double> >::template xTRightInvA<2, 2>(JT, z, dx );
438 }
while (dx.two_norm2() > epsilon*epsilon);
444 const Dune::FieldMatrix<ctype, mydimension, coorddimension>
451 const int pat[4][3] = {{ 0, 0 },
455 Dune::FieldMatrix<ctype, mydimension, coorddimension> Jt(0.0);
456 for (
int i = 0; i < 4; ++i) {
457 for (
int deriv = 0; deriv < 2; ++deriv) {
460 for (
int j = 0; j < 2; ++j) {
461 factor *= (j != deriv) ? uvw[pat[i][j]][j]
462 : (pat[i][j] == 0 ? -1.0 : 1.0);
465 corner_contrib *= factor;
466 Jt[deriv] += corner_contrib;
474 const Dune::FieldMatrix<ctype, coorddimension, mydimension>
477 Dune::FieldMatrix<ctype, coorddimension, mydimension> Jti = jacobianTransposed(local);
486 Dune::FieldMatrix<ctype, coorddimension, mydimension> Jt = jacobianTransposed(local);
487 #if DUNE_VERSION_NEWER(DUNE_GEOMETRY, 2, 5)
488 return Dune::Impl::FieldMatrixHelper<double>::template sqrtDetAAT<2, 2>(Jt);
490 using namespace Dune::GenericGeometry;
491 return MatrixHelper<DuneCoordTraits<double> >::template sqrtDetAAT<2, 2>(Jt);
496 GlobalCoordinate c[4];
504 BoundaryGrid::Vertex minXminY(std::vector<BoundaryGrid::Vertex>& in);
508 BoundaryGrid::Vertex maxXminY(std::vector<BoundaryGrid::Vertex>& in);
512 BoundaryGrid::Vertex maxXmaxY(std::vector<BoundaryGrid::Vertex>& in);
516 BoundaryGrid::Vertex minXmaxY(std::vector<BoundaryGrid::Vertex>& in);
size_t totalNodes() const
Return the total number of nodes on the grid when known.
Definition: boundarygrid.hh:166
const Dune::FieldMatrix< ctype, coorddimension, mydimension > jacobianInverseTransposed(const LocalCoordinate &local) const
Returns the inverse, transposed Jacobian.
Definition: boundarygrid.hh:475
int dir
Compare using this direction, if -1 compare using index values.
Definition: boundarygrid.hh:213
Dune::GeometryType type() const
Returns entity type (a 2D cube)
Definition: boundarygrid.hh:358
bool isFixed(int node) const
Check if a given node is marked as fixed.
Definition: boundarygrid.hh:174
Dune::FieldVector< double, 3 > GlobalCoordinate
A coordinate on the underlying 3D grid.
Definition: boundarygrid.hh:38
FaceCoord c
Coordinates of the vertex (2D)
Definition: boundarygrid.hh:73
FaceCoord bb[2]
Bounding box.
Definition: boundarygrid.hh:110
Dune::FieldMatrix< ctype, mydimension, coorddimension > JacobianTransposed
Type of transposed Jacobian matrix.
Definition: boundarygrid.hh:324
HexGeometry(const BoundaryGrid::Quad &q, const GridImp &gv, int dir)
Construct integration element extracted from a 3D grid.
Definition: boundarygrid.hh:330
Dune::FieldMatrix< ctype, coorddimension, mydimension > Jacobian
Type of Jacobian matrix.
Definition: boundarygrid.hh:321
bool fixed
Whether or not this node is fixed.
Definition: boundarygrid.hh:79
std::vector< Quad > grid
Our quadrilateral elements.
Definition: boundarygrid.hh:237
FaceCoord pos(double xi, double eta) const
Return the physical coordinates corresponding to the given local coordinates.
Definition: boundarygrid.cpp:329
virtual ~BoundaryGrid()
Default (empty) destructor.
Definition: boundarygrid.hh:57
bool operator==(const Vertex &v2)
Comparison operator.
Definition: boundarygrid.hh:82
Dune::FieldVector< ctype, coorddimension > GlobalCoordinate
Range type.
Definition: boundarygrid.hh:318
bool find(Vertex &res, const Vertex &coord) const
Locate the position of a vertex on the grid.
Definition: boundarygrid.cpp:117
Quad * q
The quad containing the vertex (if search has been done)
Definition: boundarygrid.hh:76
GlobalCoordinate global(const LocalCoordinate &local) const
Map from local coordinates to global coordinates.
Definition: boundarygrid.hh:394
VertexLess(int comp)
Default constructor.
Definition: boundarygrid.hh:202
A class describing a 2D vertex.
Definition: boundarygrid.hh:64
A class describing a linear, quadrilateral element.
Definition: boundarygrid.hh:88
BoundedPredicate(const FaceCoord &coord_)
Default constructor.
Definition: boundarygrid.hh:220
A class describing a quad grid.
Definition: boundarygrid.hh:35
FaceCoord coord
The coordinates to check.
Definition: boundarygrid.hh:233
bool cubicSolve(double eps, double A, double B, double C, double D, std::vector< double > &X) const
Solves the cubic equation A*x^3+B*x^2+C*x+D=0.
Definition: boundarygrid.cpp:270
int i
Index of the vertex.
Definition: boundarygrid.hh:70
void add(const Quad &quad)
Add a quad to the grid.
Definition: boundarygrid.cpp:85
int Q4inv(FaceCoord &res, const Quad &q, const FaceCoord &coord, double epsZero, double epsOut) const
Find the local coordinates of a given point within a given quad.
Definition: boundarygrid.cpp:150
friend std::ostream & operator<<(std::ostream &os, const Quad &q)
Print to a stream.
Definition: boundarygrid.hh:113
static BoundaryGrid uniform(const FaceCoord &min, const FaceCoord &max, int k1, int k2, bool dc=false)
Establish an uniform quad grid.
Definition: boundarygrid.cpp:24
ctype integrationElement(const LocalCoordinate &local) const
Returns the integration element (|J'*J|)^(1/2)
Definition: boundarygrid.hh:484
friend std::ostream & operator<<(std::ostream &os, const BoundaryGrid &g)
Print to a stream.
Definition: boundarygrid.hh:247
size_t size() const
Obtain the number of quads in the grid.
Definition: boundarygrid.hh:154
const Dune::FieldMatrix< ctype, mydimension, coorddimension > jacobianTransposed(const LocalCoordinate &local) const
Return the transposed jacobian.
Definition: boundarygrid.hh:445
std::vector< bool > fixNodes
Whether or not a given node is marked as fixed.
Definition: boundarygrid.hh:241
GlobalCoordinate center() const
Returns center of quadrilateral.
Definition: boundarygrid.hh:378
const Quad & operator[](int index) const
Obtain a const reference to a quad.
Definition: boundarygrid.hh:143
Dune::FieldVector< ctype, mydimension > LocalCoordinate
Domain type.
Definition: boundarygrid.hh:315
size_t nodes
Total number of nodes on grid.
Definition: boundarygrid.hh:244
Predicate for sorting vertices.
Definition: boundarygrid.hh:199
Quad & operator[](int index)
Obtain a reference to a quad.
Definition: boundarygrid.hh:135
Predicate for checking if a vertex falls within a quads bounding box.
Definition: boundarygrid.hh:217
Vertex()
Default constructor.
Definition: boundarygrid.hh:67
bool operator()(const Quad &q)
The comparison operator.
Definition: boundarygrid.hh:223
Geometry class for general hexagons.
Definition: boundarygrid.hh:290
LocalCoordinate local(const GlobalCoordinate &y) const
Map from global coordinates to local coordinates.
Definition: boundarygrid.hh:419
bool bilinearSolve(double epsilon, double order, const double *A, const double *B, std::vector< double > &X, std::vector< double > &Y) const
Solves a bi-linear set of equations in x and y.
Definition: boundarygrid.cpp:213
BoundaryGrid()
Default (empty) constructor.
Definition: boundarygrid.hh:54
Dune::FieldVector< double, 2 > FaceCoord
A coordinate on the quad grid.
Definition: boundarygrid.hh:40
double ctype
Coordinate element type.
Definition: boundarygrid.hh:312
Vertex v[4]
Vertices.
Definition: boundarygrid.hh:108
GlobalCoordinate corner(int cor) const
Returns coordinates to requested corner.
Definition: boundarygrid.hh:387
Quad()
Default constructor.
Definition: boundarygrid.hh:91
ctype volume() const
Returns volume (area) of quadrilateral.
Definition: boundarygrid.hh:372
int corners() const
Returns number of corners.
Definition: boundarygrid.hh:366
HexGeometry(const BoundaryGrid::Quad &q)
Construct integration element.
Definition: boundarygrid.hh:350
bool operator()(const Vertex &q1, const Vertex &q2)
The comparison operator.
Definition: boundarygrid.hh:205
std::vector< double > evalBasis(double xi, double eta) const
Evaluate the basis functions.
Definition: boundarygrid.cpp:340
static void extract(FaceCoord &res, const GlobalCoordinate &coord, int dir)
Helper function for extracting given 2D coordinates from a 3D coordinate.
Definition: boundarygrid.cpp:70