DUNE-FEM (unstable)

referenceelementgeometry.hh
1// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2// vi: set et ts=4 sw=2 sts=2:
3// SPDX-FileCopyrightInfo: Copyright © DUNE Project contributors, see file LICENSE.md in module root
4// SPDX-License-Identifier: LicenseRef-GPL-2.0-only-with-DUNE-exception
5#ifndef DUNE_GEOMETRY_TEST_REFERENCEELEMENTGEOMETRY_HH
6#define DUNE_GEOMETRY_TEST_REFERENCEELEMENTGEOMETRY_HH
7
8#include <type_traits>
9
10#include <dune/geometry/referenceelements.hh>
11
12namespace Dune {
13namespace Impl {
14
15// Placeholder type for a trivial identity matrix without any functionality
16struct IdentityMatrix
17{
18 // multiply Id * A
19 template <class A>
20 friend const A& operator* (IdentityMatrix, const A& a) { return a; }
21
22 // multiply A * Id
23 template <class A>
24 friend const A& operator* (const A& a, IdentityMatrix) { return a; }
25
26 // multiply Id * Id
27 friend IdentityMatrix operator* (IdentityMatrix, IdentityMatrix) { return {}; }
28
29 friend std::ostream& operator<< (std::ostream& out, IdentityMatrix)
30 {
31 return out << "I";
32 }
33
34 // cast into FieldMatrix
35 template <class K, int n>
36 operator FieldMatrix<K,n,n> () const
37 {
38 FieldMatrix<K,n,n> I;
39 for (int i = 0; i < n; ++i)
40 I[i][i] = K(1);
41 return I;
42 }
43};
44
45
61template <class RefElem>
62class ReferenceElementGeometry
63 : public RefElem::template Codim<0>::Geometry
64{
65 using Base = typename RefElem::template Codim<0>::Geometry;
66
67public:
68 using LocalCoordinate = typename RefElem::Coordinate;
69 using GlobalCoordinate = typename RefElem::Coordinate;
70 using Jacobian = Impl::IdentityMatrix;
71 using JacobianTransposed = Impl::IdentityMatrix;
72 using JacobianInverse = Impl::IdentityMatrix;
73 using JacobianInverseTransposed = Impl::IdentityMatrix;
74
75public:
76 constexpr explicit ReferenceElementGeometry (const RefElem& refElem)
77 : Base{refElem.template geometry<0>(0)}
78 {}
79
81 constexpr const LocalCoordinate& local (const GlobalCoordinate& global) const noexcept
82 {
83 return global;
84 }
85
87 constexpr const GlobalCoordinate& global (const LocalCoordinate& local) const noexcept
88 {
89 return local;
90 }
91
93 constexpr Jacobian jacobian (const LocalCoordinate& local) const noexcept
94 {
95 return Jacobian{};
96 }
97
99 constexpr JacobianTransposed jacobianTransposed (const LocalCoordinate& local) const noexcept
100 {
101 return JacobianTransposed{};
102 }
103
105 constexpr JacobianInverse jacobianInverse (const LocalCoordinate& local) const noexcept
106 {
107 return JacobianInverse{};
108 }
109
111 constexpr JacobianInverseTransposed jacobianInverseTransposed (const LocalCoordinate& local) const noexcept
112 {
113 return JacobianInverseTransposed{};
114 }
115};
116
117
134template <class Geometry>
135class LocalDerivativeGeometry
136 : public ReferenceElementGeometry<
137 typename Dune::ReferenceElements<typename Geometry::ctype, Geometry::mydimension>::ReferenceElement >
138{
140 using ReferenceElement = typename ReferenceElements::ReferenceElement;
141 using Base = ReferenceElementGeometry<ReferenceElement>;
142
143public:
144 using LocalCoordinate = typename Geometry::LocalCoordinate;
145 using Jacobian = typename Geometry::Jacobian;
146 using JacobianTransposed = typename Geometry::JacobianTransposed;
147 using JacobianInverse = typename Geometry::JacobianInverse;
148 using JacobianInverseTransposed = typename Geometry::JacobianInverseTransposed;
149
150public:
151 explicit LocalDerivativeGeometry (const Geometry& geometry) noexcept
152 : Base(referenceElement(geometry))
153 , geometry_(geometry)
154 {}
155
157 Jacobian jacobian (const LocalCoordinate& local) const
158 {
159 return geometry_.jacobian(local);
160 }
161
163 JacobianTransposed jacobianTransposed (const LocalCoordinate& local) const
164 {
165 return geometry_.jacobianTransposed(local);
166 }
167
169 JacobianInverse jacobianInverse (const LocalCoordinate& local) const
170 {
171 return geometry_.jacobianInverse(local);
172 }
173
175 JacobianInverseTransposed jacobianInverseTransposed (const LocalCoordinate& local) const
176 {
177 return geometry_.jacobianInverseTransposed(local);
178 }
179
180private:
181 Geometry geometry_;
182};
183
184} // end namespace Impl
185} // end namespace Dune
186
187#endif // DUNE_GEOMETRY_TEST_REFERENCEELEMENTGEOMETRY_HH
Implementation::JacobianTransposed JacobianTransposed
type of jacobian transposed
Definition: geometry.hh:131
Implementation::JacobianInverseTransposed JacobianInverseTransposed
type of jacobian inverse transposed
Definition: geometry.hh:120
Std::detected_or_t< JacobianInverseDefault, JacobianInverseOfImplementation, Implementation > JacobianInverse
type of jacobian inverse
Definition: geometry.hh:178
FieldVector< ctype, mydim > LocalCoordinate
type of local coordinates
Definition: geometry.hh:103
Std::detected_or_t< JacobianDefault, JacobianOfImplementation, Implementation > Jacobian
type of jacobian
Definition: geometry.hh:189
unspecified value type referenceElement(T &&... t)
Returns a reference element for the objects t....
Dune namespace.
Definition: alignedallocator.hh:13
Class providing access to the singletons of the reference elements.
Definition: referenceelements.hh:128
typename Container::ReferenceElement ReferenceElement
The reference element type.
Definition: referenceelements.hh:146
Creative Commons License   |  Legal Statements / Impressum  |  Hosted by TU Dresden  |  generated with Hugo v0.111.3 (Nov 12, 23:30, 2024)