35 #ifndef OPM_VOLUMES_HEADER
36 #define OPM_VOLUMES_HEADER
41 #include <opm/grid/utility/platform_dependent/disable_warnings.h>
43 #include <dune/common/version.hh>
44 #include <dune/common/math.hh>
45 #include <dune/common/fvector.hh>
47 #include <opm/grid/utility/platform_dependent/reenable_warnings.h>
58 FieldVector<T, 3>
cross(
const FieldVector<T, 3>& a,
const FieldVector<T, 3>& b)
60 FieldVector<T, 3> res;
61 res[0] = a[1]*b[2] - a[2]*b[1];
62 res[1] = a[2]*b[0] - a[0]*b[2];
63 res[2] = a[0]*b[1] - a[1]*b[0];
72 template <
class Vector>
73 typename Vector::field_type
inner(
const Vector& a,
const Vector& b)
75 return std::inner_product(a.begin(), a.end(), b.begin(),
typename Vector::field_type());
81 template<
typename T,
template <
typename,
int>
class Point>
84 return a[0][0] * a[1][1] - a[1][0] * a[0][1];
91 template<
typename T,
template <
typename,
int>
class Point>
95 a[0][0] * (a[1][1] * a[2][2] - a[2][1] * a[1][2]) -
96 a[0][1] * (a[1][0] * a[2][2] - a[2][0] * a[1][2]) +
97 a[0][2] * (a[1][0] * a[2][1] - a[2][0] * a[1][1]);
103 template<
typename T,
template <
typename,
int>
class Point,
int Dim>
106 Point<T, Dim> tmp[Dim];
107 for (
int i = 0; i < Dim; ++i) {
108 tmp[i] = a[i+1] - a[i];
110 return determinantOf(tmp) / double(Factorial<Dim>::factorial);
117 template <
typename T,
template <
typename,
int>
class Point>
118 inline T
area(
const Point<T, 2>* c)
124 template <
typename T,
template <
typename,
int>
class Point>
125 inline T
area(
const Point<T, 3>* c)
128 Point<T, 3> d0 = c[1] - c[0];
129 Point<T, 3> d1 = c[2] - c[0];
130 Point<T, 3> crossprod =
cross(d0,d1);
131 return 0.5 * crossprod.two_norm();
136 template <
typename T,
template <
typename,
int>
class Point>
143 template <
typename T,
template <
typename,
int>
class Point>
147 Point<T, 3> d0 = c[1] - c[0];
148 Point<T, 3> d1 = c[2] - c[0];
149 Point<T, 3> crossprod =
cross(d0, d1);
150 if (
inner(crossprod, normal) > 0) {
151 return 0.5 * crossprod.two_norm();
153 return -0.5 * crossprod.two_norm();
162 #endif // OPM_VOLUMES_HEADER
T determinantOf(const Point< T, 2 > *a)
Calculates the determinant of a 2 x 2 matrix, represented in memory as an array of two-dimensional po...
Definition: Volumes.hpp:82
Vector::field_type inner(const Vector &a, const Vector &b)
Definition: Volumes.hpp:73
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
T volume(const Point< T, 3 > *c)
Computes the volume of a 3D simplex (embedded i 3D space).
Definition: Volumes.hpp:137
T area(const Point< T, 2 > *c)
Computes the area of a 2-dimensional triangle.
Definition: Volumes.hpp:118
FieldVector< T, 3 > cross(const FieldVector< T, 3 > &a, const FieldVector< T, 3 > &b)
Definition: Volumes.hpp:58
T signed_area(const Point< T, 3 > *c, const Point< T, 3 > &normal)
Computes the signed area of a triangle embedded in 3D space.
Definition: Volumes.hpp:144