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#include <dune/grid/common/mcmgmapper.hh>
22#include <dune/grid/common/defaultgridview.hh>
24#include <opm/common/utility/platform_dependent/reenable_warnings.h>
49 int k1,
int k2,
bool dc=
false);
82 return hypot(v2.
c[0]-
c[0],v2.
c[1]-
c[1]) < 1.e-8;
91 v[0].
i =
v[1].
i =
v[2].
i =
v[3].
i = 0;
92 v[0].
c =
v[1].
c =
v[2].
c =
v[3].
c = 0.f;
103 std::vector<double>
evalBasis(
double xi,
double eta)
const;
113 os <<
"Nodes: " << q.
v[0].
i <<
"," << q.
v[1].
i <<
"," << q.
v[2].
i <<
"," << q.
v[3].
i << std::endl;
114 os <<
"(" << q.
v[0].
c <<
")(" << q.
v[1].
c <<
")(" << q.
v[2].
c <<
")(" << q.
v[3].
c <<
")";
121 void add(
const Quad& quad);
123 void addToColumn(
size_t col,
const Quad& quad)
125 if (col >= colGrids.size())
126 colGrids.resize(col+1);
127 colGrids[col].push_back(quad);
146 const Quad& getQuad(
int col,
int index)
const
148 return colGrids[col][index];
157 size_t colSize(
int i)
const
159 return colGrids[i].size();
180 bool find(Vertex& res,
const Vertex& coord)
const;
193 static void extract(Vertex& res,
194 const Dune::FieldVector<double,3>& coord,
int dir);
224 return (
coord[0] >= q.
bb[0][0]-eps &&
236 std::vector<std::vector<Quad> > colGrids;
247 for (
size_t i=0;i<g.
size();++i)
248 os << g[i] << std::endl;
262 const double* A,
const double* B,
263 std::vector<double>& X,
264 std::vector<double>& Y)
const;
273 bool cubicSolve(
double eps,
double A,
double B,
double C,
274 double D, std::vector<double>& X)
const;
283 double epsZero,
double epsOut)
const;
287 template<
int mydim,
int cdim,
class Gr
idImp>
293 template<
int cdim,
class Gr
idImp>
298 enum { dimension = 2};
301 enum { mydimension = 2};
304 enum { coorddimension = cdim };
307 enum {dimensionworld = 2 };
319 typedef Dune::FieldMatrix<ctype,coorddimension,mydimension>
Jacobian;
330 using LeafGridView = Dune::GridView<Dune::DefaultLeafGridViewTraits<GridImp>>;
331 Dune::MultipleCodimMultipleGeomTypeMapper<LeafGridView> mapper(gv.leafGridView(), Dune::mcmgVertexLayout());
332 typename LeafGridView::template Codim<3>::Iterator start=gv.leafGridView().template begin<3>();
333 const typename LeafGridView::template Codim<3>::Iterator itend = gv.leafGridView().template end<3>();
334 for (
int i=0;i<4;++i) {
335 typename LeafGridView::template Codim<3>::Iterator it=start;
336 for (; it != itend; ++it) {
337 if (mapper.index(*it) == q.
v[i].
i)
343 m_volume = (c[1][0]-c[0][0])*(c[2][1]-c[0][1]);
350 for (
int i=0;i<4;++i)
352 m_volume = (c[1][0]-c[0][0])*(c[2][1]-c[0][1]);
356 Dune::GeometryType
type()
const
358 Dune::GeometryType t;
359 t = Dune::GeometryTypes::cube(mydimension);
380 return Global(local);
398 const int pat[4][2] = {{ 0, 0 },
403 for (
int i = 0; i < 4; ++i) {
406 for (
int j = 0; j < 2; ++j) {
407 factor *= uvw[pat[i][j]][j];
409 corner_contrib *= factor;
410 xyz += corner_contrib;
419 const ctype epsilon = 1e-10;
420 auto refElement = Dune::ReferenceElements<ctype,2>::cube();
428 Dune::Impl::FieldMatrixHelper<double>::template xTRightInvA<2, 2>(JT, z, dx );
430 }
while (dx.two_norm2() > epsilon*epsilon);
436 const Dune::FieldMatrix<ctype, mydimension, coorddimension>
443 const int pat[4][3] = {{ 0, 0 },
447 Dune::FieldMatrix<ctype, mydimension, coorddimension> Jt(0.0);
448 for (
int i = 0; i < 4; ++i) {
449 for (
int deriv = 0; deriv < 2; ++deriv) {
452 for (
int j = 0; j < 2; ++j) {
453 factor *= (j != deriv) ? uvw[pat[i][j]][j]
454 : (pat[i][j] == 0 ? -1.0 : 1.0);
457 corner_contrib *= factor;
458 Jt[deriv] += corner_contrib;
466 const Dune::FieldMatrix<ctype, coorddimension, mydimension>
469 Dune::FieldMatrix<ctype, coorddimension, mydimension> Jti = jacobianTransposed(local);
478 Dune::FieldMatrix<ctype, coorddimension, mydimension> Jt = jacobianTransposed(local);
479 return Dune::Impl::FieldMatrixHelper<double>::template sqrtDetAAT<2, 2>(Jt);
484 GlobalCoordinate c[4];
492BoundaryGrid::Vertex minXminY(
const std::vector<BoundaryGrid::Vertex>& in);
496BoundaryGrid::Vertex maxXminY(
const std::vector<BoundaryGrid::Vertex>& in);
500BoundaryGrid::Vertex maxXmaxY(
const std::vector<BoundaryGrid::Vertex>& in);
504BoundaryGrid::Vertex minXmaxY(
const std::vector<BoundaryGrid::Vertex>& in);
A class describing a linear, quadrilateral element.
Definition boundarygrid.hh:86
Quad()
Default constructor.
Definition boundarygrid.hh:89
Vertex v[4]
Vertices.
Definition boundarygrid.hh:106
FaceCoord bb[2]
Bounding box.
Definition boundarygrid.hh:108
friend std::ostream & operator<<(std::ostream &os, const Quad &q)
Print to a stream.
Definition boundarygrid.hh:111
std::vector< double > evalBasis(double xi, double eta) const
Evaluate the basis functions.
Definition boundarygrid.cpp:340
FaceCoord pos(double xi, double eta) const
Return the physical coordinates corresponding to the given local coordinates.
Definition boundarygrid.cpp:329
A class describing a 2D vertex.
Definition boundarygrid.hh:62
bool operator==(const Vertex &v2)
Comparison operator.
Definition boundarygrid.hh:80
Quad * q
The quad containing the vertex (if search has been done)
Definition boundarygrid.hh:74
bool fixed
Whether or not this node is fixed.
Definition boundarygrid.hh:77
FaceCoord c
Coordinates of the vertex (2D)
Definition boundarygrid.hh:71
int i
Index of the vertex.
Definition boundarygrid.hh:68
Vertex()
Default constructor.
Definition boundarygrid.hh:65
A class describing a quad grid.
Definition boundarygrid.hh:33
std::vector< Quad > grid
Our quadrilateral elements.
Definition boundarygrid.hh:235
std::vector< bool > fixNodes
Whether or not a given node is marked as fixed.
Definition boundarygrid.hh:239
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
Quad & operator[](int index)
Obtain a reference to a quad.
Definition boundarygrid.hh:133
virtual ~BoundaryGrid()
Default (empty) destructor.
Definition boundarygrid.hh:55
Dune::FieldVector< double, 3 > GlobalCoordinate
A coordinate on the underlying 3D grid.
Definition boundarygrid.hh:36
size_t totalNodes() const
Return the total number of nodes on the grid when known.
Definition boundarygrid.hh:164
size_t nodes
Total number of nodes on grid.
Definition boundarygrid.hh:242
BoundaryGrid()
Default (empty) constructor.
Definition boundarygrid.hh:52
const Quad & operator[](int index) const
Obtain a const reference to a quad.
Definition boundarygrid.hh:141
friend std::ostream & operator<<(std::ostream &os, const BoundaryGrid &g)
Print to a stream.
Definition boundarygrid.hh:245
void add(const Quad &quad)
Add a quad to the grid.
Definition boundarygrid.cpp:85
bool find(Vertex &res, const Vertex &coord) const
Locate the position of a vertex on the grid.
Definition boundarygrid.cpp:117
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
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
Dune::FieldVector< double, 2 > FaceCoord
A coordinate on the quad grid.
Definition boundarygrid.hh:38
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
bool isFixed(int node) const
Check if a given node is marked as fixed.
Definition boundarygrid.hh:172
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
size_t size() const
Obtain the number of quads in the grid.
Definition boundarygrid.hh:152
Dune::FieldMatrix< ctype, coorddimension, mydimension > Jacobian
Type of Jacobian matrix.
Definition boundarygrid.hh:319
HexGeometry(const BoundaryGrid::Quad &q, const GridImp &gv, int dir)
Construct integration element extracted from a 3D grid.
Definition boundarygrid.hh:328
int corners() const
Returns number of corners.
Definition boundarygrid.hh:364
GlobalCoordinate center() const
Returns center of quadrilateral.
Definition boundarygrid.hh:376
LocalCoordinate local(const GlobalCoordinate &y) const
Map from global coordinates to local coordinates.
Definition boundarygrid.hh:417
const Dune::FieldMatrix< ctype, coorddimension, mydimension > jacobianInverseTransposed(const LocalCoordinate &local) const
Returns the inverse, transposed Jacobian.
Definition boundarygrid.hh:467
Dune::FieldVector< ctype, mydimension > LocalCoordinate
Domain type.
Definition boundarygrid.hh:313
Dune::GeometryType type() const
Returns entity type (a 2D cube)
Definition boundarygrid.hh:356
ctype volume() const
Returns volume (area) of quadrilateral.
Definition boundarygrid.hh:370
double ctype
Coordinate element type.
Definition boundarygrid.hh:310
Dune::FieldVector< ctype, coorddimension > GlobalCoordinate
Range type.
Definition boundarygrid.hh:316
Dune::FieldMatrix< ctype, mydimension, coorddimension > JacobianTransposed
Type of transposed Jacobian matrix.
Definition boundarygrid.hh:322
GlobalCoordinate corner(int cor) const
Returns coordinates to requested corner.
Definition boundarygrid.hh:385
ctype integrationElement(const LocalCoordinate &local) const
Returns the integration element (|J'*J|)^(1/2)
Definition boundarygrid.hh:476
const Dune::FieldMatrix< ctype, mydimension, coorddimension > jacobianTransposed(const LocalCoordinate &local) const
Return the transposed jacobian.
Definition boundarygrid.hh:437
HexGeometry(const BoundaryGrid::Quad &q)
Construct integration element.
Definition boundarygrid.hh:348
GlobalCoordinate global(const LocalCoordinate &local) const
Map from local coordinates to global coordinates.
Definition boundarygrid.hh:392
Geometry class for general hexagons.
Definition boundarygrid.hh:289
Inverting small matrices.
Definition ImplicitAssembly.hpp:43
Predicate for checking if a vertex falls within a quads bounding box.
Definition boundarygrid.hh:215
bool operator()(const Quad &q)
The comparison operator.
Definition boundarygrid.hh:221
FaceCoord coord
The coordinates to check.
Definition boundarygrid.hh:231
BoundedPredicate(const FaceCoord &coord_)
Default constructor.
Definition boundarygrid.hh:218
Predicate for sorting vertices.
Definition boundarygrid.hh:197
VertexLess(int comp)
Default constructor.
Definition boundarygrid.hh:200
int dir
Compare using this direction, if -1 compare using index values.
Definition boundarygrid.hh:211
bool operator()(const Vertex &q1, const Vertex &q2)
The comparison operator.
Definition boundarygrid.hh:203