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
56namespace Dune
57{
58 namespace cpgrid
59 {
60
69 template <int mydim, int cdim>
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
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
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
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);
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
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
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
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
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
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
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
const GlobalCoordinate & center() const
Returns the centroid of the geometry.
Definition: Geometry.hpp:169
Geometry()
Default constructor, giving a non-valid geometry.
Definition: Geometry.hpp:118
const GlobalCoordinate & global(const LocalCoordinate &) const
Returns the position of the vertex.
Definition: Geometry.hpp:124
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
double ctype
Coordinate element type.
Definition: Geometry.hpp:93
JacobianInverse jacobianInverse(const LocalCoordinate &) const
The inverse of the jacobian.
Definition: Geometry.hpp:742
const FieldMatrix< ctype, mydimension, coorddimension > & jacobianTransposed(const LocalCoordinate &) const
This method is meaningless for singular geometries.
Definition: Geometry.hpp:721
const GlobalCoordinate & center() const
Returns the centroid of the geometry.
Definition: Geometry.hpp:714
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
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
const FieldMatrix< ctype, coorddimension, mydimension > & jacobianInverseTransposed(const LocalCoordinate &) const
This method is meaningless for singular geometries.
Definition: Geometry.hpp:728
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 GlobalCoordinate & global(const LocalCoordinate &) const
This method is meaningless for singular geometries.
Definition: Geometry.hpp:667
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
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
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
int corners() const
The number of corners of this convex polytope.
Definition: Geometry.hpp:365
ctype volume() const
Cell volume.
Definition: Geometry.hpp:378
const GlobalCoordinate & center() const
Returns the centroid of the geometry.
Definition: Geometry.hpp:388
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
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
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