dune-fem 2.8.0
Loading...
Searching...
No Matches
transformed.hh
Go to the documentation of this file.
1#ifndef DUNE_FEM_SPACE_BASISFUNCTIONSET_TRANSFORMED_HH
2#define DUNE_FEM_SPACE_BASISFUNCTIONSET_TRANSFORMED_HH
3
4// C++ includes
5#include <cassert>
6#include <cstddef>
7
8// dune-geometry includes
9#include <dune/geometry/referenceelements.hh>
10#include <dune/geometry/type.hh>
11
12// dune-fem includes
17
18
19namespace Dune
20{
21
22 namespace Fem
23 {
24
25 // TransformedBasisFunctionSet
26 // ---------------------------
27
42 template< class Entity, class ShapeFunctionSet, class Transformation >
44 {
46
47 public:
49 typedef Entity EntityType;
52
53 protected:
54 typedef typename EntityType::Geometry GeometryType;
55 typedef typename GeometryType::ctype ctype;
56 typedef typename GeometryType::JacobianTransposed JacobianTransposed;
57
58 public:
61
70
72
74 typedef std::decay_t< decltype( Dune::ReferenceElements< ctype, GeometryType::coorddimension >::general( std::declval< const Dune::GeometryType & >() ) ) > ReferenceElementType;
75
78 : entity_( nullptr )
79 {}
80
83 : entity_( &entity ),
84 shapeFunctionSet_( shapeFunctionSet )
85 {}
86
87
88 // Basis Function Set Interface Methods
89 // ------------------------------------
90
92 int order () const { return shapeFunctionSet().order(); }
93
95 std::size_t size () const { return shapeFunctionSet().size(); }
96
98 auto referenceElement () const
99 -> decltype( Dune::ReferenceElements< ctype, GeometryType::coorddimension >::general( std::declval< const Dune::GeometryType & >() ) )
100 {
101 return Dune::ReferenceElements< ctype, GeometryType::coorddimension >::general( type() );
102 }
103
107 template< class QuadratureType, class Vector, class DofVector >
108 void axpy ( const QuadratureType &quad, const Vector &values, DofVector &dofs ) const
109 {
110 // call axpy method for each entry of the given vector, e.g. rangeVector or jacobianVector
111 const unsigned int nop = quad.nop();
112 for( unsigned int qp = 0; qp < nop; ++qp )
113 axpy( quad[ qp ], values[ qp ], dofs );
114 }
115
122 template< class QuadratureType, class VectorA, class VectorB, class DofVector >
123 void axpy ( const QuadratureType &quad, const VectorA &valuesA, const VectorB &valuesB, DofVector &dofs ) const
124 {
125 // call axpy method for each entry of the given vector, e.g. rangeVector or jacobianVector
126 const unsigned int nop = quad.nop();
127 for( unsigned int qp = 0; qp < nop; ++qp )
128 {
129 axpy( quad[ qp ], valuesA[ qp ], dofs );
130 axpy( quad[ qp ], valuesB[ qp ], dofs );
131 }
132 }
133
137 template< class Point, class DofVector >
138 void axpy ( const Point &x, const RangeType &valueFactor, DofVector &dofs ) const
139 {
140 RangeType factor = transformation( coordinate( x ) ).apply_t( valueFactor );
143 }
144
148 template< class Point, class DofVector >
149 void axpy ( const Point &x, const JacobianRangeType &jacobianFactor, DofVector &dofs ) const
150 {
151 typedef typename GeometryType::JacobianInverseTransposed GeometryJacobianInverseTransposedType;
152 const GeometryType &geo = geometry();
153 const GeometryJacobianInverseTransposedType &gjit = geo.jacobianInverseTransposed( coordinate( x ) );
154 JacobianRangeType tmpJacobianFactor( RangeFieldType( 0 ) );
155 for( int r = 0; r < FunctionSpaceType::dimRange; ++r )
156 gjit.mtv( jacobianFactor[ r ], tmpJacobianFactor[ r ] );
157
158 tmpJacobianFactor = transformation( coordinate( x ) ).apply_t( tmpJacobianFactor );
161 }
162
166 template< class Point, class DofVector >
167 void axpy( const Point &x, const HessianRangeType &hessianFactor, DofVector &dofs ) const
168 {
169 DUNE_THROW( NotImplemented, "hessian axpy for TransformedBasisFunctionSet not implemented." );
170 }
171
175 template< class Point, class DofVector >
176 void axpy ( const Point &x, const RangeType &valueFactor, const JacobianRangeType &jacobianFactor,
177 DofVector &dofs ) const
178 {
179 axpy( x, valueFactor, dofs );
180 axpy( x, jacobianFactor, dofs );
181 }
182
184 template< class QuadratureType, class DofVector, class RangeArray >
185 void evaluateAll ( const QuadratureType &quad, const DofVector &dofs, RangeArray &ranges ) const
186 {
187 // call axpy method for each entry of the given vector, e.g. rangeVector or jacobianVector
188 const unsigned int nop = quad.nop();
189 for( unsigned int qp = 0; qp < nop; ++qp )
190 evaluateAll( quad[ qp ], dofs, ranges[ qp ] );
191 }
192
194 template< class Point, class DofVector >
195 void evaluateAll ( const Point &x, const DofVector &dofs, RangeType &value ) const
196 {
197 value = 0;
200 value = transformation( coordinate( x ) ).apply( value );
201 }
202
204 template< class Point, class RangeArray >
205 void evaluateAll ( const Point &x, RangeArray &values ) const
206 {
207 assert( values.size() >= size() );
208 AssignFunctor< RangeArray > f( values );
210 auto transform = transformation( coordinate( x ) );
211 for( auto &e : values )
212 e = transform.apply( e );
213 }
214
216 template< class QuadratureType, class DofVector, class JacobianArray >
217 void jacobianAll ( const QuadratureType &quad, const DofVector &dofs, JacobianArray &jacobians ) const
218 {
219 // call axpy method for each entry of the given vector, e.g. rangeVector or jacobianVector
220 const unsigned int nop = quad.nop();
221 for( unsigned int qp = 0; qp < nop; ++qp )
222 jacobianAll( quad[ qp ], dofs, jacobians[ qp ] );
223 }
224
226 template< class Point, class DofVector >
227 void jacobianAll ( const Point &x, const DofVector &dofs, JacobianRangeType &jacobian ) const
228 {
229 JacobianRangeType localJacobian( RangeFieldType( 0 ) );
232 const GeometryType &geo = geometry();
233 JacobianTransformation< GeometryType >( geo, coordinate( x ) )( localJacobian, jacobian );
234 jacobian = transformation( coordinate( x ) ).apply( jacobian );
235 }
236
238 template< class Point, class JacobianRangeArray >
239 void jacobianAll ( const Point &x, JacobianRangeArray &jacobians ) const
240 {
241 assert( jacobians.size() >= size() );
242 const GeometryType &geo = geometry();
247 auto transform = transformation( coordinate( x ) );
248 for( auto &jacobian : jacobians )
249 jacobian = transform.apply( jacobian );
250 }
251
253 template< class Point, class DofVector >
254 void hessianAll ( const Point &x, const DofVector &dofs, HessianRangeType &hessian ) const
255 {
256 DUNE_THROW( NotImplemented, "hessianAll for TransformedBasisFunctionSet not implemented." );
257 }
258
260 template< class Point, class HessianRangeArray >
261 void hessianAll ( const Point &x, HessianRangeArray &hessians ) const
262 {
263 DUNE_THROW( NotImplemented, "hessianAll for TransformedBasisFunctionSet not implemented." );
264 }
265
267 const Entity &entity () const
268 {
269 assert( entity_ );
270 return *entity_;
271 }
272
274 bool valid () const { return bool(entity_); }
275
277 Dune::GeometryType type () const { return entity().type(); }
278
279
280 // Non-interface methods
281 // ---------------------
282
284 const ShapeFunctionSetType &shapeFunctionSet () const { return shapeFunctionSet_; }
285
286 Transformation transformation ( const DomainType& x ) const
287 {
288 return Transformation( geometry(), x );
289 }
290
291 protected:
292 GeometryType geometry () const { return entity().geometry(); }
293
294 private:
295 const EntityType *entity_;
296 ShapeFunctionSetType shapeFunctionSet_;
297 };
298
299 } // namespace Fem
300
301} // namespace Dune
302
303#endif // #ifndef DUNE_FEM_SPACE_BASISFUNCTIONSET_TRANSFORMED_HH
STL namespace.
Definition: bindguard.hh:11
void jacobianTransformation(const GeometryJacobianInverseTransposed &gjit, const FieldMatrix< K, ROWS, GeometryJacobianInverseTransposed::cols > &a, FieldMatrix< K, ROWS, GeometryJacobianInverseTransposed::rows > &b)
Definition: transformation.hh:21
IteratorRange< typename DF::DofIteratorType > dofs(DF &df)
Iterates over all DOFs.
Definition: rangegenerators.hh:76
static const Point & coordinate(const Point &x)
Definition: coordinate.hh:14
Definition: explicitfieldvector.hh:75
Definition: misc/functor.hh:31
Definition: space/basisfunctionset/functor.hh:108
Definition: space/basisfunctionset/functor.hh:132
Definition: transformation.hh:36
implementation of a basis function set for given entity
Definition: transformed.hh:44
void hessianAll(const Point &x, HessianRangeArray &hessians) const
Definition: transformed.hh:261
const ShapeFunctionSetType & shapeFunctionSet() const
return shape function set
Definition: transformed.hh:284
void jacobianAll(const Point &x, JacobianRangeArray &jacobians) const
Definition: transformed.hh:239
ShapeFunctionSet ShapeFunctionSetType
shape function set type
Definition: transformed.hh:51
void hessianAll(const Point &x, const DofVector &dofs, HessianRangeType &hessian) const
Definition: transformed.hh:254
GeometryType geometry() const
Definition: transformed.hh:292
void evaluateAll(const QuadratureType &quad, const DofVector &dofs, RangeArray &ranges) const
Definition: transformed.hh:185
EntityType::Geometry GeometryType
Definition: transformed.hh:54
int order() const
return order of basis function set
Definition: transformed.hh:92
void axpy(const QuadratureType &quad, const Vector &values, DofVector &dofs) const
evaluate all basis function and multiply with given values and add to dofs
Definition: transformed.hh:108
void evaluateAll(const Point &x, RangeArray &values) const
Definition: transformed.hh:205
GeometryType::JacobianTransposed JacobianTransposed
Definition: transformed.hh:56
std::size_t size() const
return size of basis function set
Definition: transformed.hh:95
void axpy(const Point &x, const JacobianRangeType &jacobianFactor, DofVector &dofs) const
evaluate all basis function and multiply with given values and add to dofs
Definition: transformed.hh:149
TransformedBasisFunctionSet(const EntityType &entity, const ShapeFunctionSet &shapeFunctionSet=ShapeFunctionSet())
constructor
Definition: transformed.hh:82
void axpy(const QuadratureType &quad, const VectorA &valuesA, const VectorB &valuesB, DofVector &dofs) const
evaluate all basis function and multiply with given values and add to dofs
Definition: transformed.hh:123
GeometryType::ctype ctype
Definition: transformed.hh:55
bool valid() const
return true if entity pointer is set
Definition: transformed.hh:274
Entity EntityType
entity type
Definition: transformed.hh:49
Dune::GeometryType type() const
return geometry type
Definition: transformed.hh:277
auto referenceElement() const -> decltype(Dune::ReferenceElements< ctype, GeometryType::coorddimension >::general(std::declval< const Dune::GeometryType & >()))
return reference element
Definition: transformed.hh:98
FunctionSpaceType::RangeType RangeType
range type
Definition: transformed.hh:65
void axpy(const Point &x, const RangeType &valueFactor, DofVector &dofs) const
evaluate all basis function and multiply with given values and add to dofs
Definition: transformed.hh:138
std::decay_t< decltype(Dune::ReferenceElements< ctype, GeometryType::coorddimension >::general(std::declval< const Dune::GeometryType & >())) > ReferenceElementType
type of reference element
Definition: transformed.hh:74
FunctionSpaceType::JacobianRangeType JacobianRangeType
jacobian range type
Definition: transformed.hh:67
ShapeFunctionSetType::FunctionSpaceType FunctionSpaceType
type of function space
Definition: transformed.hh:60
FunctionSpaceType::DomainType DomainType
domain type
Definition: transformed.hh:63
Transformation transformation(const DomainType &x) const
Definition: transformed.hh:286
const Entity & entity() const
return entity
Definition: transformed.hh:267
void evaluateAll(const Point &x, const DofVector &dofs, RangeType &value) const
Definition: transformed.hh:195
FunctionSpaceType::HessianRangeType HessianRangeType
hessian range type
Definition: transformed.hh:69
void axpy(const Point &x, const RangeType &valueFactor, const JacobianRangeType &jacobianFactor, DofVector &dofs) const
evaluate all basis function and multiply with given values and add to dofs
Definition: transformed.hh:176
TransformedBasisFunctionSet()
constructor
Definition: transformed.hh:77
void jacobianAll(const QuadratureType &quad, const DofVector &dofs, JacobianArray &jacobians) const
Definition: transformed.hh:217
void axpy(const Point &x, const HessianRangeType &hessianFactor, DofVector &dofs) const
evaluate all basis function and multiply with given values and add to dofs
Definition: transformed.hh:167
void jacobianAll(const Point &x, const DofVector &dofs, JacobianRangeType &jacobian) const
Definition: transformed.hh:227
FunctionSpaceType::RangeFieldType RangeFieldType
Definition: transformed.hh:71
A vector valued function space.
Definition: functionspace.hh:60
FunctionSpaceTraits::RangeType RangeType
Type of range vector (using type of range field) has a Dune::FieldVector type interface.
Definition: functionspaceinterface.hh:71
@ dimRange
dimension of range vector space
Definition: functionspaceinterface.hh:48
FunctionSpaceTraits::LinearMappingType JacobianRangeType
Intrinsic type used for the jacobian values has a Dune::FieldMatrix type interface.
Definition: functionspaceinterface.hh:75
FunctionSpaceTraits::DomainType DomainType
Type of domain vector (using type of domain field) has a Dune::FieldVector type interface.
Definition: functionspaceinterface.hh:67
FunctionSpaceTraits::RangeFieldType RangeFieldType
Intrinsic type used for values in the range field (usually a double)
Definition: functionspaceinterface.hh:63
Interface class for shape function sets.
Definition: shapefunctionset/shapefunctionset.hh:33
void evaluateEach(const Point &x, Functor functor) const
evalute each shape function
std::size_t size() const
return number of shape functions
int order() const
return order of shape functions
void jacobianEach(const Point &x, Functor functor) const
evalute jacobian of each shape function