My Project
Geometry.hpp
1 //===========================================================================
2 //
3 // File: Geometry.hpp
4 //
5 // Created: Fri May 29 23:29:24 2009
6 //
7 // Author(s): Atgeirr F Rasmussen <atgeirr@sintef.no>
8 // B�rd Skaflestad <bard.skaflestad@sintef.no>
9 //
10 // $Date$
11 //
12 // $Revision$
13 //
14 //===========================================================================
15 
16 /*
17  Copyright 2009, 2010, 2011 SINTEF ICT, Applied Mathematics.
18  Copyright 2009, 2010, 2011, 2022 Equinor ASA.
19 
20  This file is part of the Open Porous Media project (OPM).
21 
22  OPM is free software: you can redistribute it and/or modify
23  it under the terms of the GNU General Public License as published by
24  the Free Software Foundation, either version 3 of the License, or
25  (at your option) any later version.
26 
27  OPM is distributed in the hope that it will be useful,
28  but WITHOUT ANY WARRANTY; without even the implied warranty of
29  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
30  GNU General Public License for more details.
31 
32  You should have received a copy of the GNU General Public License
33  along with OPM. If not, see <http://www.gnu.org/licenses/>.
34 */
35 
36 #ifndef OPM_GEOMETRY_HEADER
37 #define OPM_GEOMETRY_HEADER
38 
39 #include <cmath>
40 
41 // Warning suppression for Dune includes.
42 #include <opm/grid/utility/platform_dependent/disable_warnings.h>
43 
44 #include <dune/common/version.hh>
45 #include <dune/geometry/referenceelements.hh>
46 #include <dune/grid/common/geometry.hh>
47 
48 #include <dune/geometry/type.hh>
49 
50 #include <opm/grid/cpgrid/EntityRep.hpp>
51 #include <opm/grid/common/Volumes.hpp>
52 #include <opm/grid/utility/platform_dependent/reenable_warnings.h>
53 
54 #include <opm/grid/utility/ErrorMacros.hpp>
55 
56 namespace Dune
57 {
58  namespace cpgrid
59  {
60 
69  template <int mydim, int cdim>
70  class Geometry
71  {
72  };
73 
74 
75 
76 
78  template <int cdim> // GridImp arg never used
79  class Geometry<0, cdim>
80  {
81  static_assert(cdim == 3, "");
82  public:
84  enum { dimension = 3 };
86  enum { mydimension = 0};
88  enum { coorddimension = cdim };
90  enum { dimensionworld = 3 };
91 
93  typedef double ctype;
94 
96  typedef FieldVector<ctype, mydimension> LocalCoordinate;
98  typedef FieldVector<ctype, coorddimension> GlobalCoordinate;
99 
101  typedef FieldMatrix< ctype, coorddimension, mydimension > Jacobian;
103  typedef FieldMatrix< ctype, coorddimension, mydimension > JacobianInverse;
105  typedef FieldMatrix< ctype, mydimension, coorddimension > JacobianTransposed;
107  typedef FieldMatrix< ctype, coorddimension, mydimension > JacobianInverseTransposed;
108 
109 
113  : pos_(pos)
114  {
115  }
116 
119  : pos_(0.0)
120  {
121  }
122 
125  {
126  return pos_;
127  }
128 
131  {
132  // return 0 to make the geometry check happy.
133  return LocalCoordinate(0.0);
134  }
135 
137  double integrationElement(const LocalCoordinate&) const
138  {
139  return volume();
140  }
141 
143  GeometryType type() const
144  {
145  return Dune::GeometryTypes::cube(mydimension);
146  }
147 
149  int corners() const
150  {
151  return 1;
152  }
153 
155  GlobalCoordinate corner(int cor) const
156  {
157  static_cast<void>(cor);
158  assert(cor == 0);
159  return pos_;
160  }
161 
163  ctype volume() const
164  {
165  return 1.0;
166  }
167 
169  const GlobalCoordinate& center() const
170  {
171  return pos_;
172  }
173 
175  JacobianTransposed
176  jacobianTransposed(const LocalCoordinate& /* local */) const
177  {
178 
179  // Meaningless to call jacobianTransposed() on singular geometries. But we need to make DUNE happy.
180  return {};
181  }
182 
184  JacobianInverseTransposed
186  {
187  // Meaningless to call jacobianInverseTransposed() on singular geometries. But we need to make DUNE happy.
188  return {};
189  }
190 
192  Jacobian
193  jacobian(const LocalCoordinate& /*local*/) const
194  {
195  return {};
196  }
197 
199  JacobianInverse
200  jacobianInverse(const LocalCoordinate& /*local*/) const
201  {
202  return {};
203  }
204 
206  bool affine() const
207  {
208  return true;
209  }
210 
211  private:
212  GlobalCoordinate pos_;
213  };
214 
215 
216 
217 
219  template <int cdim>
220  class Geometry<3, cdim>
221  {
222  static_assert(cdim == 3, "");
223  public:
225  enum { dimension = 3 };
227  enum { mydimension = 3 };
229  enum { coorddimension = cdim };
231  enum { dimensionworld = 3 };
232 
234  typedef double ctype;
235 
237  typedef FieldVector<ctype, mydimension> LocalCoordinate;
239  typedef FieldVector<ctype, coorddimension> GlobalCoordinate;
240 
242  typedef FieldMatrix< ctype, coorddimension, mydimension > Jacobian;
244  typedef FieldMatrix< ctype, coorddimension, mydimension > JacobianInverse;
246  typedef FieldMatrix< ctype, mydimension, coorddimension > JacobianTransposed;
248  typedef FieldMatrix< ctype, coorddimension, mydimension > JacobianInverseTransposed;
249 
250  typedef Dune::Impl::FieldMatrixHelper< double > MatrixHelperType;
251 
261  ctype vol,
262  const EntityVariable<cpgrid::Geometry<0, 3>, 3>& allcorners,
263  const int* corner_indices)
264  : pos_(pos), vol_(vol), allcorners_(allcorners.data()), cor_idx_(corner_indices)
265  {
266  assert(allcorners_ && corner_indices);
267  }
268 
279  ctype vol)
280  : pos_(pos), vol_(vol)
281  {
282  }
283 
286  : pos_(0.0), vol_(0.0), allcorners_(0), cor_idx_(0)
287  {
288  }
289 
294  GlobalCoordinate global(const LocalCoordinate& local_coord) const
295  {
296  static_assert(mydimension == 3, "");
297  static_assert(coorddimension == 3, "");
298  // uvw = { (1-u, 1-v, 1-w), (u, v, w) }
299  LocalCoordinate uvw[2] = { LocalCoordinate(1.0), local_coord };
300  uvw[0] -= local_coord;
301  // Access pattern for uvw matching ordering of corners.
302  const int pat[8][3] = { { 0, 0, 0 },
303  { 1, 0, 0 },
304  { 0, 1, 0 },
305  { 1, 1, 0 },
306  { 0, 0, 1 },
307  { 1, 0, 1 },
308  { 0, 1, 1 },
309  { 1, 1, 1 } };
310  GlobalCoordinate xyz(0.0);
311  for (int i = 0; i < 8; ++i) {
312  GlobalCoordinate corner_contrib = corner(i);
313  double factor = 1.0;
314  for (int j = 0; j < 3; ++j) {
315  factor *= uvw[pat[i][j]][j];
316  }
317  corner_contrib *= factor;
318  xyz += corner_contrib;
319  }
320  return xyz;
321  }
322 
326  {
327  static_assert(mydimension == 3, "");
328  static_assert(coorddimension == 3, "");
329  // This code is modified from dune/grid/genericgeometry/mapping.hh
330  // \todo: Implement direct computation.
331  const ctype epsilon = 1e-12;
332  auto refElement = Dune::ReferenceElements<ctype, 3>::cube();
333  LocalCoordinate x = refElement.position(0,0);
334  LocalCoordinate dx;
335  do {
336  // DF^n dx^n = F^n, x^{n+1} -= dx^n
337  JacobianTransposed JT = jacobianTransposed(x);
338  GlobalCoordinate z = global(x);
339  z -= y;
340  MatrixHelperType::template xTRightInvA<3, 3>(JT, z, dx );
341  x -= dx;
342  } while (dx.two_norm2() > epsilon*epsilon);
343  return x;
344  }
345 
350  double integrationElement(const LocalCoordinate& local_coord) const
351  {
352  JacobianTransposed Jt = jacobianTransposed(local_coord);
353  return MatrixHelperType::template sqrtDetAAT<3, 3>(Jt);
354  }
355 
358  GeometryType type() const
359  {
360  return Dune::GeometryTypes::cube(mydimension);
361  }
362 
365  int corners() const
366  {
367  return 8;
368  }
369 
371  GlobalCoordinate corner(int cor) const
372  {
373  assert(allcorners_ && cor_idx_);
374  return allcorners_[cor_idx_[cor]].center();
375  }
376 
378  ctype volume() const
379  {
380  return vol_;
381  }
382 
383  void set_volume(ctype volume) {
384  vol_ = volume;
385  }
386 
388  const GlobalCoordinate& center() const
389  {
390  return pos_;
391  }
392 
397  const JacobianTransposed
398  jacobianTransposed(const LocalCoordinate& local_coord) const
399  {
400  static_assert(mydimension == 3, "");
401  static_assert(coorddimension == 3, "");
402 
403  // uvw = { (1-u, 1-v, 1-w), (u, v, w) }
404  LocalCoordinate uvw[2] = { LocalCoordinate(1.0), local_coord };
405  uvw[0] -= local_coord;
406  // Access pattern for uvw matching ordering of corners.
407  const int pat[8][3] = { { 0, 0, 0 },
408  { 1, 0, 0 },
409  { 0, 1, 0 },
410  { 1, 1, 0 },
411  { 0, 0, 1 },
412  { 1, 0, 1 },
413  { 0, 1, 1 },
414  { 1, 1, 1 } };
415  JacobianTransposed Jt(0.0);
416  for (int i = 0; i < 8; ++i) {
417  for (int deriv = 0; deriv < 3; ++deriv) {
418  // This part contributing to dg/du_{deriv}
419  double factor = 1.0;
420  for (int j = 0; j < 3; ++j) {
421  factor *= (j != deriv) ? uvw[pat[i][j]][j]
422  : (pat[i][j] == 0 ? -1.0 : 1.0);
423  }
424  GlobalCoordinate corner_contrib = corner(i);
425  corner_contrib *= factor;
426  Jt[deriv] += corner_contrib; // using FieldMatrix row access.
427  }
428  }
429  return Jt;
430  }
431 
433  const JacobianInverseTransposed
434  jacobianInverseTransposed(const LocalCoordinate& local_coord) const
435  {
436  JacobianInverseTransposed Jti = jacobianTransposed(local_coord);
437  Jti.invert();
438  return Jti;
439  }
440 
442  Jacobian
443  jacobian(const LocalCoordinate& local_coord) const
444  {
445  return jacobianTransposed(local_coord).transposed();
446  }
447 
449  JacobianInverse
450  jacobianInverse(const LocalCoordinate& local_coord) const
451  {
452  return jacobianInverseTransposed(local_coord).transposed();
453  }
454 
456  bool affine() const
457  {
458  return false;
459  }
460 
477  std::vector<Geometry<3, cdim>> refine(const std::array<int, 3>& cells_per_dim,
478  std::vector<EntityVariable<Geometry<0, 3>, 3>>& corner_storage,
479  std::vector<std::array<int, 8>>& indices_storage)
480  {
481  // The center of the parent in local coordinates.
482  const Geometry<3, cdim>::LocalCoordinate parent_center(this->local(this->center()));
483 
484  // Corners of the parent hexahedron in order, in local coordinates.
485  const Geometry<3, cdim>::LocalCoordinate parent_corners[8] = {
486  {0.0, 0.0, 0.0},
487  {1.0, 0.0, 0.0},
488  {0.0, 1.0, 0.0},
489  {1.0, 1.0, 0.0},
490  {0.0, 0.0, 1.0},
491  {1.0, 0.0, 1.0},
492  {0.0, 1.0, 1.0},
493  {1.0, 1.0, 1.0},
494  };
495 
496  // Indices of the corners of the 6 faces of the hexahedrons.
497  const int face_corner_indices[6][4] = {
498  {0, 1, 2, 3},
499  {0, 1, 4, 5},
500  {0, 2, 4, 6},
501  {1, 3, 5, 7},
502  {2, 3, 6, 7},
503  {4, 5, 6, 7},
504  };
505 
506  // To calculate a refined cell's volume, the hexahedron is
507  // divided in 24 tetrahedrons, each of which is defined by the
508  // center of the cell, the center of one face, and by one edge
509  // of that face. This struct defines that edge for each face,
510  // for each of the four possible tetrahedrons that are based on
511  // that face.
512  const int tetra_edge_indices[6][4][2] = {
513  {{0, 1}, {0, 2}, {1, 3}, {2, 3}},
514  {{0, 1}, {0, 4}, {1, 5}, {4, 5}},
515  {{0, 2}, {0, 4}, {2, 6}, {4, 6}},
516  {{1, 3}, {1, 5}, {3, 7}, {5, 7}},
517  {{2, 3}, {2, 6}, {3, 7}, {6, 7}},
518  {{4, 5}, {4, 6}, {5, 7}, {6, 7}},
519  };
520 
521 
522  std::vector<Geometry<3, cdim>> result;
523  result.reserve(cells_per_dim[0] * cells_per_dim[1] * cells_per_dim[2]);
524 
525  auto pcs = corner_storage.begin();
526  auto pis = indices_storage.begin();
527 
528  Geometry<3, cdim>::ctype total_volume = 0.0;
529  for (int k = 0; k < cells_per_dim[2]; k++) {
530  Geometry<3, cdim>::LocalCoordinate refined_corners[8];
531  Geometry<3, cdim>::LocalCoordinate refined_center(0.0);
532 
533  refined_center[2] = (parent_center[2] + k) / cells_per_dim[2];
534  for (int h = 0; h < 8; h++) {
535  refined_corners[h][2] = (parent_corners[h][2] + k) / cells_per_dim[2];
536  }
537  for (int j = 0; j < cells_per_dim[1]; j++) {
538  refined_center[1] = (parent_center[1] + j) / cells_per_dim[1];
539  for (int h = 0; h < 8; h++) {
540  refined_corners[h][1] = (parent_corners[h][1] + j) / cells_per_dim[1];
541  }
542  for (int i = 0; i < cells_per_dim[0]; i++) {
543  refined_center[0] = (parent_center[0] + i) / cells_per_dim[0];
544  for (int h = 0; h < 8; h++) {
545  refined_corners[h][0] = (parent_corners[h][0] + i) / cells_per_dim[0];
546  }
547 
548  auto& global_refined_corners = *pcs++;
549  global_refined_corners.reserve(8);
550  for (const auto& corner : refined_corners) {
551  global_refined_corners.push_back(Geometry<0, 3>(this->global(corner)));
552  }
553 
554  // The indices must match the order of the constant
555  // arrays containing unit corners, face indices, and
556  // tetrahedron edge indices. Do not reorder.
557  auto& indices = *pis++;
558  indices = {0, 1, 2, 3, 4, 5, 6, 7};
559 
560  // Get the center of the cell.
561  const Geometry<3, cdim>::GlobalCoordinate global_refined_center(
562  this->global(refined_center));
563 
564  // Calculate the centers of the 6 faces.
565  const auto& hex_corners = global_refined_corners.data();
566  Geometry<0, 3>::GlobalCoordinate face_centers[6];
567  for (int f = 0; f < 6; f++) {
568  face_centers[f] = hex_corners[face_corner_indices[f][0]].center();
569  face_centers[f] += hex_corners[face_corner_indices[f][1]].center();
570  face_centers[f] += hex_corners[face_corner_indices[f][2]].center();
571  face_centers[f] += hex_corners[face_corner_indices[f][3]].center();
572  face_centers[f] /= 4;
573  }
574 
575  // Calculate the volume of the cell by adding the 4 tetrahedrons at each face.
577  for (int f = 0; f < 6; f++) {
578  for (int e = 0; e < 4; e++) {
579  const Geometry<0, 3>::GlobalCoordinate tetra_corners[4]
580  = {hex_corners[tetra_edge_indices[f][e][0]].center(),
581  hex_corners[tetra_edge_indices[f][e][1]].center(),
582  face_centers[f],
583  global_refined_center};
584  volume += std::fabs(simplex_volume(tetra_corners));
585  }
586  }
587  total_volume += volume;
588 
589  result.push_back(Geometry<3, cdim>(
590  global_refined_center, volume, global_refined_corners, indices.data()));
591  }
592  }
593  }
594 
595  // Rescale all volumes if the sum of volumes does not match the parent.
596  if (std::fabs(total_volume - this->volume())
597  > std::numeric_limits<Geometry<3, cdim>::ctype>::epsilon()) {
598  Geometry<3, cdim>::ctype correction = this->volume() / total_volume;
599  for (auto& r : result) {
600  r.set_volume(r.volume() * correction);
601  }
602  }
603 
604  return result;
605  }
606 
607  private:
608  GlobalCoordinate pos_;
609  double vol_;
610  const cpgrid::Geometry<0, 3>* allcorners_; // For dimension 3 only
611  const int* cor_idx_; // For dimension 3 only
612  };
613 
614 
615 
616 
617 
620  template <int cdim> // GridImp arg never used
621  class Geometry<2, cdim>
622  {
623  static_assert(cdim == 3, "");
624  public:
626  enum { dimension = 3 };
628  enum { mydimension = 2 };
630  enum { coorddimension = cdim };
632  enum { dimensionworld = 3 };
633 
635  typedef double ctype;
636 
638  typedef FieldVector<ctype, mydimension> LocalCoordinate;
640  typedef FieldVector<ctype, coorddimension> GlobalCoordinate;
641 
643  typedef FieldMatrix< ctype, coorddimension, mydimension > Jacobian;
645  typedef FieldMatrix< ctype, coorddimension, mydimension > JacobianInverse;
647  typedef FieldMatrix< ctype, mydimension, coorddimension > JacobianTransposed;
649  typedef FieldMatrix< ctype, coorddimension, mydimension > JacobianInverseTransposed;
650 
655  ctype vol)
656  : pos_(pos), vol_(vol)
657  {
658  }
659 
662  : pos_(0.0), vol_(0.0)
663  {
664  }
665 
668  {
669  OPM_THROW(std::runtime_error, "Geometry::global() meaningless on singular geometry.");
670  }
671 
674  {
675  OPM_THROW(std::runtime_error, "Geometry::local() meaningless on singular geometry.");
676  }
677 
680  double integrationElement(const LocalCoordinate&) const
681  {
682  return vol_;
683  }
684 
686  GeometryType type() const
687  {
688  return Dune::GeometryTypes::none(mydimension);
689  }
690 
693  int corners() const
694  {
695  return 0;
696  }
697 
699  GlobalCoordinate corner(int /* cor */) const
700  {
701  // Meaningless call to cpgrid::Geometry::corner(int):
702  //"singular geometry has no corners.
703  // But the DUNE tests assume at least one corner.
704  return GlobalCoordinate( 0.0 );
705  }
706 
708  ctype volume() const
709  {
710  return vol_;
711  }
712 
714  const GlobalCoordinate& center() const
715  {
716  return pos_;
717  }
718 
720  const FieldMatrix<ctype, mydimension, coorddimension>&
721  jacobianTransposed(const LocalCoordinate& /* local */) const
722  {
723  OPM_THROW(std::runtime_error, "Meaningless to call jacobianTransposed() on singular geometries.");
724  }
725 
727  const FieldMatrix<ctype, coorddimension, mydimension>&
729  {
730  OPM_THROW(std::runtime_error, "Meaningless to call jacobianInverseTransposed() on singular geometries.");
731  }
732 
734  Jacobian
735  jacobian(const LocalCoordinate& /*local*/) const
736  {
737  return jacobianTransposed({}).transposed();
738  }
739 
741  JacobianInverse
742  jacobianInverse(const LocalCoordinate& /*local*/) const
743  {
744  return jacobianInverseTransposed({}).transposed();
745  }
746 
748  bool affine() const
749  {
750  return true;
751  }
752 
753  private:
754  GlobalCoordinate pos_;
755  ctype vol_;
756  };
757 
758 
759 
760 
761 
762  } // namespace cpgrid
763 
764  template< int mydim, int cdim >
765  auto referenceElement(const cpgrid::Geometry<mydim,cdim>& geo) -> decltype(referenceElement<double,mydim>(geo.type()))
766  {
767  return referenceElement<double,mydim>(geo.type());
768  }
769 
770 } // namespace Dune
771 
772 #endif // OPM_GEOMETRY_HEADER
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:98
bool affine() const
The mapping implemented by this geometry is constant, therefore affine.
Definition: Geometry.hpp:206
FieldMatrix< ctype, coorddimension, mydimension > Jacobian
Type of Jacobian matrix.
Definition: Geometry.hpp:101
JacobianInverse jacobianInverse(const LocalCoordinate &) const
This method is meaningless for singular geometries.
Definition: Geometry.hpp:200
FieldVector< ctype, mydimension > LocalCoordinate
Domain type of.
Definition: Geometry.hpp:96
Jacobian jacobian(const LocalCoordinate &) const
This method is meaningless for singular geometries.
Definition: Geometry.hpp:193
const GlobalCoordinate & center() const
Returns the centroid of the geometry.
Definition: Geometry.hpp:169
double integrationElement(const LocalCoordinate &) const
Returns 1 for the vertex geometry.
Definition: Geometry.hpp:137
FieldMatrix< ctype, coorddimension, mydimension > JacobianInverseTransposed
Type of the inverse of the transposed Jacobian matrix.
Definition: Geometry.hpp:107
Geometry(const GlobalCoordinate &pos)
Construct from vertex position.
Definition: Geometry.hpp:112
FieldMatrix< ctype, coorddimension, mydimension > JacobianInverse
Type of inverse of Jacobian matrix.
Definition: Geometry.hpp:103
JacobianTransposed jacobianTransposed(const LocalCoordinate &) const
This method is meaningless for singular geometries.
Definition: Geometry.hpp:176
GeometryType type() const
Using the cube type for vertices.
Definition: Geometry.hpp:143
ctype volume() const
Volume of vertex is arbitrarily set to 1.
Definition: Geometry.hpp:163
FieldMatrix< ctype, mydimension, coorddimension > JacobianTransposed
Type of transposed Jacobian matrix.
Definition: Geometry.hpp:105
GlobalCoordinate corner(int cor) const
Returns the single corner: the vertex itself.
Definition: Geometry.hpp:155
int corners() const
A vertex is defined by a single corner.
Definition: Geometry.hpp:149
Geometry()
Default constructor, giving a non-valid geometry.
Definition: Geometry.hpp:118
JacobianInverseTransposed jacobianInverseTransposed(const LocalCoordinate &) const
This method is meaningless for singular geometries.
Definition: Geometry.hpp:185
LocalCoordinate local(const GlobalCoordinate &) const
Meaningless for the vertex geometry.
Definition: Geometry.hpp:130
const GlobalCoordinate & global(const LocalCoordinate &) const
Returns the position of the vertex.
Definition: Geometry.hpp:124
double ctype
Coordinate element type.
Definition: Geometry.hpp:93
JacobianInverse jacobianInverse(const LocalCoordinate &) const
The inverse of the jacobian.
Definition: Geometry.hpp:742
const GlobalCoordinate & global(const LocalCoordinate &) const
This method is meaningless for singular geometries.
Definition: Geometry.hpp:667
int corners() const
The number of corners of this convex polytope.
Definition: Geometry.hpp:693
LocalCoordinate local(const GlobalCoordinate &) const
This method is meaningless for singular geometries.
Definition: Geometry.hpp:673
FieldMatrix< ctype, coorddimension, mydimension > Jacobian
Type of Jacobian matrix.
Definition: Geometry.hpp:643
bool affine() const
Since integrationElement() is constant, returns true.
Definition: Geometry.hpp:748
const FieldMatrix< ctype, coorddimension, mydimension > & jacobianInverseTransposed(const LocalCoordinate &) const
This method is meaningless for singular geometries.
Definition: Geometry.hpp:728
FieldVector< ctype, mydimension > LocalCoordinate
Domain type of.
Definition: Geometry.hpp:638
FieldMatrix< ctype, coorddimension, mydimension > JacobianInverseTransposed
Type of the inverse of the transposed Jacobian matrix.
Definition: Geometry.hpp:649
ctype volume() const
Volume (area, actually) of intersection.
Definition: Geometry.hpp:708
GeometryType type() const
We use the singular type (None) for intersections.
Definition: Geometry.hpp:686
FieldMatrix< ctype, mydimension, coorddimension > JacobianTransposed
Type of transposed Jacobian matrix.
Definition: Geometry.hpp:647
double integrationElement(const LocalCoordinate &) const
For the singular geometry, we return a constant integration element equal to the volume.
Definition: Geometry.hpp:680
Geometry()
Default constructor, giving a non-valid geometry.
Definition: Geometry.hpp:661
const FieldMatrix< ctype, mydimension, coorddimension > & jacobianTransposed(const LocalCoordinate &) const
This method is meaningless for singular geometries.
Definition: Geometry.hpp:721
FieldMatrix< ctype, coorddimension, mydimension > JacobianInverse
Type of Jacobian matrix.
Definition: Geometry.hpp:645
Geometry(const GlobalCoordinate &pos, ctype vol)
Construct from centroid and volume (1- and 0-moments).
Definition: Geometry.hpp:654
double ctype
Coordinate element type.
Definition: Geometry.hpp:635
Jacobian jacobian(const LocalCoordinate &) const
The jacobian.
Definition: Geometry.hpp:735
const GlobalCoordinate & center() const
Returns the centroid of the geometry.
Definition: Geometry.hpp:714
FieldVector< ctype, coorddimension > GlobalCoordinate
Range type of.
Definition: Geometry.hpp:640
GlobalCoordinate corner(int) const
This method is meaningless for singular geometries.
Definition: Geometry.hpp:699
Specialization for 3-dimensional geometries, i.e. cells.
Definition: Geometry.hpp:221
const JacobianInverseTransposed jacobianInverseTransposed(const LocalCoordinate &local_coord) const
Inverse of Jacobian transposed.
Definition: Geometry.hpp:434
FieldMatrix< ctype, coorddimension, mydimension > JacobianInverseTransposed
Type of the inverse of the transposed Jacobian matrix.
Definition: Geometry.hpp:248
bool affine() const
The mapping implemented by this geometry is not generally affine.
Definition: Geometry.hpp:456
double ctype
Coordinate element type.
Definition: Geometry.hpp:234
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:260
FieldVector< ctype, coorddimension > GlobalCoordinate
Range type of.
Definition: Geometry.hpp:239
GlobalCoordinate corner(int cor) const
The 8 corners of the hexahedral base cell.
Definition: Geometry.hpp:371
FieldMatrix< ctype, mydimension, coorddimension > JacobianTransposed
Type of transposed Jacobian matrix.
Definition: Geometry.hpp:246
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:358
std::vector< Geometry< 3, cdim > > refine(const std::array< int, 3 > &cells_per_dim, std::vector< EntityVariable< Geometry< 0, 3 >, 3 >> &corner_storage, std::vector< std::array< int, 8 >> &indices_storage)
Refine a single cell with regular intervals.
Definition: Geometry.hpp:477
double integrationElement(const LocalCoordinate &local_coord) const
Equal to \sqrt{\det{J^T J}} where J is the Jacobian.
Definition: Geometry.hpp:350
const JacobianTransposed jacobianTransposed(const LocalCoordinate &local_coord) const
Jacobian transposed.
Definition: Geometry.hpp:398
GlobalCoordinate global(const LocalCoordinate &local_coord) const
Provide a trilinear mapping.
Definition: Geometry.hpp:294
Geometry()
Default constructor, giving a non-valid geometry.
Definition: Geometry.hpp:285
FieldMatrix< ctype, coorddimension, mydimension > JacobianInverse
Type of inverse of Jacobian matrix.
Definition: Geometry.hpp:244
LocalCoordinate local(const GlobalCoordinate &y) const
Mapping from the cell to the reference domain.
Definition: Geometry.hpp:325
const GlobalCoordinate & center() const
Returns the centroid of the geometry.
Definition: Geometry.hpp:388
int corners() const
The number of corners of this convex polytope.
Definition: Geometry.hpp:365
ctype volume() const
Cell volume.
Definition: Geometry.hpp:378
Jacobian jacobian(const LocalCoordinate &local_coord) const
The jacobian.
Definition: Geometry.hpp:443
Geometry(const GlobalCoordinate &pos, ctype vol)
Construct from centroid and volume (1- and 0-moments).
Definition: Geometry.hpp:278
FieldVector< ctype, mydimension > LocalCoordinate
Domain type of.
Definition: Geometry.hpp:237
JacobianInverse jacobianInverse(const LocalCoordinate &local_coord) const
The inverse of the jacobian.
Definition: Geometry.hpp:450
FieldMatrix< ctype, coorddimension, mydimension > Jacobian
Type of Jacobian matrix.
Definition: Geometry.hpp:242
This class encapsulates geometry for both vertices, intersections and cells.
Definition: Geometry.hpp:71
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
T simplex_volume(const Point< T, Dim > *a)
Computes the volume of a simplex consisting of (Dim+1) vertices embedded in Euclidean space of dimens...
Definition: Volumes.hpp:104