DUNE-FEM (unstable)

transformed.hh
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
13#include <dune/fem/common/coordinate.hh>
14#include <dune/fem/space/basisfunctionset/default.hh>
15#include <dune/fem/space/basisfunctionset/functor.hh>
16#include <dune/fem/space/basisfunctionset/transformation.hh>
17#include <dune/fem/space/common/functionspace.hh>
18
19
20namespace Dune
21{
22
23 namespace Fem
24 {
25
26 // TransformedBasisFunctionSet
27 // ---------------------------
28
43 template< class Entity, class ShapeFunctionSet, class Transformation >
45 : public BasisFunctionSetStorage< Entity, ShapeFunctionSet >
46 {
47 typedef BasisFunctionSetStorage< Entity, ShapeFunctionSet > BaseType;
49
50 public:
52 typedef typename BaseType :: EntityType EntityType;
53
55 typedef typename BaseType :: Geometry Geometry;
56
58 typedef typename BaseType :: ShapeFunctionSetType ShapeFunctionSetType;
59
60 protected:
61 typedef Geometry GeometryType;
62 typedef typename GeometryType::ctype ctype;
63 typedef typename GeometryType::JacobianTransposed JacobianTransposed;
64
65 public:
68
77
78 typedef typename FunctionSpaceType::RangeFieldType RangeFieldType;
79
81 typedef typename BaseType :: ReferenceElementType ReferenceElementType;
82
85
87 explicit TransformedBasisFunctionSet ( const EntityType &entity, const ShapeFunctionSet &shapeFunctionSet = ShapeFunctionSet() )
88 : BaseType( entity, shapeFunctionSet )
89 {
90 }
91
93 TransformedBasisFunctionSet &operator= ( const TransformedBasisFunctionSet &other ) = default;
94
95 using BaseType :: entity;
96 using BaseType :: geometry;
97 using BaseType :: type;
98 using BaseType :: shapeFunctionSet;
99 using BaseType :: order;
100 using BaseType :: size;
101 using BaseType :: referenceElement;
102
106 template< class QuadratureType, class Vector, class DofVector >
107 void axpy ( const QuadratureType &quad, const Vector &values, DofVector &dofs ) const
108 {
109 // call axpy method for each entry of the given vector, e.g. rangeVector or jacobianVector
110 const unsigned int nop = quad.nop();
111 for( unsigned int qp = 0; qp < nop; ++qp )
112 axpy( quad[ qp ], values[ qp ], dofs );
113 }
114
121 template< class QuadratureType, class VectorA, class VectorB, class DofVector >
122 void axpy ( const QuadratureType &quad, const VectorA &valuesA, const VectorB &valuesB, DofVector &dofs ) const
123 {
124 // call axpy method for each entry of the given vector, e.g. rangeVector or jacobianVector
125 const unsigned int nop = quad.nop();
126 for( unsigned int qp = 0; qp < nop; ++qp )
127 {
128 axpy( quad[ qp ], valuesA[ qp ], dofs );
129 axpy( quad[ qp ], valuesB[ qp ], dofs );
130 }
131 }
132
136 template< class Point, class DofVector >
137 void axpy ( const Point &x, const RangeType &valueFactor, DofVector &dofs ) const
138 {
139 RangeType factor = transformation( coordinate( x ) ).apply_t( valueFactor );
140 FunctionalAxpyFunctor< RangeType, DofVector > f( factor, dofs );
141 shapeFunctionSet().evaluateEach( x, f );
142 }
143
147 template< class Point, class DofVector >
148 void axpy ( const Point &x, const JacobianRangeType &jacobianFactor, DofVector &dofs ) const
149 {
150 typedef typename GeometryType::JacobianInverseTransposed GeometryJacobianInverseTransposedType;
151 const Geometry &geo = geometry();
152 const GeometryJacobianInverseTransposedType &gjit = geo.jacobianInverseTransposed( coordinate( x ) );
153 JacobianRangeType tmpJacobianFactor( RangeFieldType( 0 ) );
154 for( int r = 0; r < FunctionSpaceType::dimRange; ++r )
155 gjit.mtv( jacobianFactor[ r ], tmpJacobianFactor[ r ] );
156
157 tmpJacobianFactor = transformation( coordinate( x ) ).apply_t( tmpJacobianFactor );
158 FunctionalAxpyFunctor< JacobianRangeType, DofVector > f( tmpJacobianFactor, dofs );
159 shapeFunctionSet().jacobianEach( x, f );
160 }
161
165 template< class Point, class DofVector >
166 void axpy( const Point &x, const HessianRangeType &hessianFactor, DofVector &dofs ) const
167 {
168 DUNE_THROW( NotImplemented, "hessian axpy for TransformedBasisFunctionSet not implemented." );
169 }
170
174 template< class Point, class DofVector >
175 void axpy ( const Point &x, const RangeType &valueFactor, const JacobianRangeType &jacobianFactor,
176 DofVector &dofs ) const
177 {
178 axpy( x, valueFactor, dofs );
179 axpy( x, jacobianFactor, dofs );
180 }
181
183 template< class QuadratureType, class DofVector, class RangeArray >
184 void evaluateAll ( const QuadratureType &quad, const DofVector &dofs, RangeArray &ranges ) const
185 {
186 // call axpy method for each entry of the given vector, e.g. rangeVector or jacobianVector
187 const unsigned int nop = quad.nop();
188 for( unsigned int qp = 0; qp < nop; ++qp )
189 evaluateAll( quad[ qp ], dofs, ranges[ qp ] );
190 }
191
193 template< class Point, class DofVector >
194 void evaluateAll ( const Point &x, const DofVector &dofs, RangeType &value ) const
195 {
196 value = 0;
197 AxpyFunctor< DofVector, RangeType > f( dofs, value );
198 shapeFunctionSet().evaluateEach( x, f );
199 value = transformation( coordinate( x ) ).apply( value );
200 }
201
203 template< class Point, class RangeArray >
204 void evaluateAll ( const Point &x, RangeArray &values ) const
205 {
206 assert( values.size() >= size() );
207 AssignFunctor< RangeArray > f( values );
208 shapeFunctionSet().evaluateEach( x, f );
209 auto transform = transformation( coordinate( x ) );
210 for( auto &e : values )
211 e = transform.apply( e );
212 }
213
215 template< class QuadratureType, class DofVector, class JacobianArray >
216 void jacobianAll ( const QuadratureType &quad, const DofVector &dofs, JacobianArray &jacobians ) const
217 {
218 // call axpy method for each entry of the given vector, e.g. rangeVector or jacobianVector
219 const unsigned int nop = quad.nop();
220 for( unsigned int qp = 0; qp < nop; ++qp )
221 jacobianAll( quad[ qp ], dofs, jacobians[ qp ] );
222 }
223
225 template< class Point, class DofVector >
226 void jacobianAll ( const Point &x, const DofVector &dofs, JacobianRangeType &jacobian ) const
227 {
228 JacobianRangeType localJacobian( RangeFieldType( 0 ) );
229 AxpyFunctor< DofVector, JacobianRangeType > f( dofs, localJacobian );
230 shapeFunctionSet().jacobianEach( x, f );
231 const Geometry &geo = geometry();
232 JacobianTransformation< GeometryType >( geo, coordinate( x ) )( localJacobian, jacobian );
233 jacobian = transformation( coordinate( x ) ).apply( jacobian );
234 }
235
237 template< class Point, class JacobianRangeArray >
238 void jacobianAll ( const Point &x, JacobianRangeArray &jacobians ) const
239 {
240 assert( jacobians.size() >= size() );
241 const Geometry &geo = geometry();
242 typedef JacobianTransformation< GeometryType > JacobianTransformation;
243 JacobianTransformation jacobianTransformation( geo, coordinate( x ) );
244 AssignFunctor< JacobianRangeArray, JacobianTransformation > f( jacobians, jacobianTransformation );
245 shapeFunctionSet().jacobianEach( x, f );
246 auto transform = transformation( coordinate( x ) );
247 for( auto &jacobian : jacobians )
248 jacobian = transform.apply( jacobian );
249 }
250
252 template< class Point, class DofVector, class HessianRange >
253 void hessianAll ( const Point &x, const DofVector &dofs, HessianRange& hessian ) const
254 {
255 if( ! std::is_same< HessianRange, HessianRangeType >:: value )
256 DUNE_THROW( NotImplemented, "hessianAll: HessianRange mismatch!" );
257
258 DUNE_THROW( NotImplemented, "hessianAll for TransformedBasisFunctionSet not implemented." );
259 }
260
262 template< class Point, class HessianRangeArray >
263 void hessianAll ( const Point &x, HessianRangeArray &hessians ) const
264 {
265 DUNE_THROW( NotImplemented, "hessianAll for TransformedBasisFunctionSet not implemented." );
266 }
267
268 // Non-interface methods
269 // ---------------------
270
271 Transformation transformation ( const DomainType& x ) const
272 {
273 return Transformation( geometry(), x );
274 }
275 };
276
277 } // namespace Fem
278
279} // namespace Dune
280
281#endif // #ifndef DUNE_FEM_SPACE_BASISFUNCTIONSET_TRANSFORMED_HH
Wrapper class for entities.
Definition: entity.hh:66
implementation of entity and geometry storage for basis function set and local functions
Definition: entitygeometry.hh:35
Dune::GeometryType type() const
return geometry type
Definition: entitygeometry.hh:126
const Entity & entity() const
return entity
Definition: entitygeometry.hh:101
const Geometry & geometry() const
return geometry
Definition: entitygeometry.hh:111
const ReferenceElementType & referenceElement() const
return reference element
Definition: entitygeometry.hh:129
Definition: explicitfieldvector.hh:75
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
A vector valued function space.
Definition: functionspace.hh:60
Interface class for shape function sets.
Definition: shapefunctionset.hh:33
implementation of a basis function set for given entity
Definition: transformed.hh:46
void hessianAll(const Point &x, HessianRangeArray &hessians) const
Definition: transformed.hh:263
void jacobianAll(const Point &x, JacobianRangeArray &jacobians) const
Definition: transformed.hh:238
void evaluateAll(const QuadratureType &quad, const DofVector &dofs, RangeArray &ranges) const
Definition: transformed.hh:184
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:107
void evaluateAll(const Point &x, RangeArray &values) const
Definition: transformed.hh:204
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:148
TransformedBasisFunctionSet(const EntityType &entity, const ShapeFunctionSet &shapeFunctionSet=ShapeFunctionSet())
constructor
Definition: transformed.hh:87
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:122
BaseType::Geometry Geometry
geometry
Definition: transformed.hh:55
void hessianAll(const Point &x, const DofVector &dofs, HessianRange &hessian) const
Definition: transformed.hh:253
FunctionSpaceType::RangeType RangeType
range type
Definition: transformed.hh:72
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:137
FunctionSpaceType::JacobianRangeType JacobianRangeType
jacobian range type
Definition: transformed.hh:74
ShapeFunctionSetType::FunctionSpaceType FunctionSpaceType
type of function space
Definition: transformed.hh:67
BaseType::EntityType EntityType
entity type
Definition: transformed.hh:52
FunctionSpaceType::DomainType DomainType
domain type
Definition: transformed.hh:70
void evaluateAll(const Point &x, const DofVector &dofs, RangeType &value) const
Definition: transformed.hh:194
FunctionSpaceType::HessianRangeType HessianRangeType
hessian range type
Definition: transformed.hh:76
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:175
BaseType::ShapeFunctionSetType ShapeFunctionSetType
shape function set type
Definition: transformed.hh:58
TransformedBasisFunctionSet()
constructor
Definition: transformed.hh:84
void jacobianAll(const QuadratureType &quad, const DofVector &dofs, JacobianArray &jacobians) const
Definition: transformed.hh:216
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:166
void jacobianAll(const Point &x, const DofVector &dofs, JacobianRangeType &jacobian) const
Definition: transformed.hh:226
BaseType::ReferenceElementType ReferenceElementType
type of reference element
Definition: transformed.hh:81
Default exception for dummy implementations.
Definition: exceptions.hh:263
#define DUNE_THROW(E, m)
Definition: exceptions.hh:218
Dune namespace.
Definition: alignedallocator.hh:13
A unique label for each type of element that can occur in a grid.
Creative Commons License   |  Legal Statements / Impressum  |  Hosted by TU Dresden  |  generated with Hugo v0.111.3 (Nov 12, 23:30, 2024)