Rheolef  7.2
an efficient C++ finite element environment
 
Loading...
Searching...
No Matches
basis_on_pointset.h
Go to the documentation of this file.
1#ifndef _RHEO_BASIS_ON_POINTSET_V2_H
2#define _RHEO_BASIS_ON_POINTSET_V2_H
23
24/*Class:basis_on_pointset
25NAME: @code{basis_on_pointset} - pre-evaluated polynomial basis
26@cindex polynomial basis
27@clindex basis
28@cindex reference element
29@clindex reference_element
30SYNOPSIS:
31 @noindent
32 The @code{basis_on_pointset} class is able to memorize the evaluation
33 of a polynomial basis and its derivatives on a set of point of the reference element
34 (@pxref{reference_element iclass}).
35 The basis is described by the @code{basis} class (@pxref{basis class}).
36 The set of points could be
37 either quadrature nodes on the reference element (@pxref{quadrature iclass})
38 or Lagrange nodes associated to another basis.
39 For application of an integration of on a side, the set of nodes could
40 be defined on a specific side only.
41 In all these cases, the evaluation of polynomials could be performed
42 one time for all on the reference element and its result stored and reused
43 for all elements of the mesh: the speedup is important, espcially for
44 high order polynomials.
45
46AUTHOR: Pierre.Saramito@imag.fr
47DATE: 4 january 2018
48End:
49*/
50
51#include "rheolef/basis.h"
52#include "rheolef/quadrature.h"
53namespace rheolef {
54
55// -----------------------------------------------------------------------
56// basis evaluated on lattice of quadrature formulae
57// -----------------------------------------------------------------------
58template<class T>
60public:
61// typedefs:
62
63 typedef typename std::vector<T>::size_type size_type;
64
65 typedef enum {
68 max_mode = 2
70
71// allocators:
72
74 basis_on_pointset_rep (const std::string& name = "");
77
78// modifiers:
79
80 void reset (const std::string& name);
81
82// accessors:
83
84 std::string name() const;
85 bool is_set() const { return _mode != max_mode; }
86 size_type ndof (reference_element hat_K) const;
87 size_type nnod (reference_element hat_K) const;
88 const basis_basic<T>& get_basis() const { return _b; }
89 bool has_quadrature() const { return _mode == quad_mode; }
90 const quadrature<T>& get_quadrature() const;
91 const basis_basic<T>& get_nodal_basis() const;
92
93 template<class Value>
94 const Eigen::Matrix<Value,Eigen::Dynamic,Eigen::Dynamic>&
95 evaluate (reference_element hat_K) const;
96
97 template<class Value>
98 const Eigen::Matrix<Value,Eigen::Dynamic,Eigen::Dynamic>&
100 reference_element tilde_L,
101 const side_information_type& sid) const;
102
103 template<class Value>
104 const Eigen::Matrix<Value,Eigen::Dynamic,Eigen::Dynamic>&
105 grad_evaluate (reference_element hat_K) const;
106
107 template<class Value>
108 const Eigen::Matrix<Value,Eigen::Dynamic,Eigen::Dynamic>&
110 reference_element tilde_L,
111 const side_information_type& sid) const;
112
113// internal:
114 static basis_on_pointset_rep<T>* make_ptr (const std::string& name);
115 static std::string _make_name(
116 mode_type mode,
117 const std::string& basis_name,
118 const std::string& pointset_name);
119 static mode_type _parse_name (
120 const std::string& name,
121 std::string& basis_name,
122 std::string& node_name);
123protected:
124// data:
126 mutable mode_type _mode;
127 quadrature<T> _quad; // when mode: on quadrature pointset
128 basis_basic<T> _nb; // when mode: on nodal basis pointset
129public:
130
131// _val [tilde_K] (inod,idof)
132#define _RHEOLEF_declare_member(VALUE,MEMBER) \
133 mutable std::array< \
134 Eigen::Matrix<VALUE,Eigen::Dynamic,Eigen::Dynamic> \
135 ,reference_element::max_variant> MEMBER; \
136
142#undef _RHEOLEF_declare_member
143
144// sid_val [tilde_L][loc_isid][orient][shift] (inod,idof)
145#define _RHEOLEF_declare_member(VALUE,MEMBER) \
146 mutable std::array< \
147 std::array< \
148 std::array< \
149 std::array< \
150 Eigen::Matrix<VALUE,Eigen::Dynamic,Eigen::Dynamic>, \
151 8>, \
152 2>, \
153 8>, \
154 reference_element::max_variant> MEMBER;
155
156_RHEOLEF_declare_member(T,_sid_scalar_val)
158_RHEOLEF_declare_member(tensor_basic<T>,_sid_tensor_val)
160_RHEOLEF_declare_member(tensor4_basic<T>,_sid_tensor4_val)
161#undef _RHEOLEF_declare_member
162
163protected:
164 mutable std::array<bool,
166 mutable std::array<bool,
168 mutable std::array<bool,
170// internals:
171 void _initialize (reference_element hat_K) const;
172 void _grad_initialize (reference_element hat_K) const;
174 const Eigen::Matrix<point_basic<T>,Eigen::Dynamic,1>& hat_node) const;
176 const Eigen::Matrix<point_basic<T>,Eigen::Dynamic,1>& hat_node) const;
177 void _sid_initialize (reference_element tilde_L) const;
178 void _sid_grad_initialize (reference_element tilde_L) const;
179 void _sid_initialize (reference_element tilde_L, const side_information_type& sid) const;
182 const Eigen::Matrix<point_basic<T>,Eigen::Dynamic,1>& hat_node) const;
184 const Eigen::Matrix<point_basic<T>,Eigen::Dynamic,1>& hat_node) const;
185};
186// -----------------------------------------------------------------------
187// interface with shallow copy semantic
188// -----------------------------------------------------------------------
189//<verbatim:
190template<class T>
191class basis_on_pointset: public smart_pointer<basis_on_pointset_rep<T>>,
192 public persistent_table<basis_on_pointset<T>> {
193public:
194
197 typedef typename rep::size_type size_type;
198
199// allocators:
200
201 basis_on_pointset (const std::string& name = "");
202 basis_on_pointset (const quadrature<T>& quad, const basis_basic<T>& b);
203 basis_on_pointset (const basis_basic<T>& nb, const basis_basic<T>& b);
204
205// modifiers:
206
207 void reset (const std::string& name);
208 void set (const quadrature<T>& quad, const basis_basic<T>& b);
209 void set (const basis_basic<T>& nb, const basis_basic<T>& b);
210
211// accessors:
212
213 bool is_set() const;
214 const basis_basic<T>& get_basis() const;
215 size_type ndof (reference_element hat_K) const;
216 size_type nnod (reference_element hat_K) const;
217 bool has_quadrature() const;
218 const quadrature<T>& get_quadrature() const;
219 const basis_basic<T>& get_nodal_basis() const;
220
221 template<class Value>
222 const Eigen::Matrix<Value,Eigen::Dynamic,Eigen::Dynamic>&
223 evaluate (reference_element hat_K) const;
224
225 template<class Value>
226 const Eigen::Matrix<Value,Eigen::Dynamic,Eigen::Dynamic>&
228 reference_element tilde_L,
229 const side_information_type& sid) const;
230
231 template<class Value>
232 const Eigen::Matrix<Value,Eigen::Dynamic,Eigen::Dynamic>&
233 grad_evaluate (reference_element hat_K) const;
234
235 template<class Value>
236 const Eigen::Matrix<Value,Eigen::Dynamic,Eigen::Dynamic>&
238 reference_element tilde_L,
239 const side_information_type& sid) const;
240};
241//>verbatim:
242
243// -----------------------------------------------------------------------
244// inlined
245// -----------------------------------------------------------------------
246template<class T>
247inline
248const basis_basic<T>&
250{
251 return base::data().get_basis();
252}
253template<class T>
254inline
257{
258 return base::data().ndof (hat_K);
259}
260template<class T>
261inline
264{
265 return base::data().nnod (hat_K);
266}
267template<class T>
268inline
269bool
271{
272 return base::data().has_quadrature();
273}
274template<class T>
275inline
276const quadrature<T>&
278{
279 return base::data().get_quadrature();
280}
281template<class T>
282inline
283const basis_basic<T>&
285{
286 return base::data().get_nodal_basis();
287}
288template<class T>
289inline
290bool
292{
293 return base::data().is_set();
294}
295template<class T>
296template<class Value>
297inline
298const Eigen::Matrix<Value,Eigen::Dynamic,Eigen::Dynamic>&
300{
301 return base::data().template evaluate<Value> (hat_K);
302}
303template<class T>
304template<class Value>
305inline
306const Eigen::Matrix<Value,Eigen::Dynamic,Eigen::Dynamic>&
308{
309 return base::data().template grad_evaluate<Value> (hat_K);
310}
311template<class T>
312template<class Value>
313inline
314const Eigen::Matrix<Value,Eigen::Dynamic,Eigen::Dynamic>&
316 reference_element tilde_L,
317 const side_information_type& sid) const
318{
319 return base::data().template evaluate_on_side<Value> (tilde_L, sid);
320}
321template<class T>
322template<class Value>
323inline
324const Eigen::Matrix<Value,Eigen::Dynamic,Eigen::Dynamic>&
326 reference_element tilde_L,
327 const side_information_type& sid) const
328{
329 return base::data().template grad_evaluate_on_side<Value> (tilde_L, sid);
330}
331
332}// namespace rheolef
333#endif // _RHEO_BASIS_ON_POINTSET_V2_H
field::size_type size_type
Definition branch.cc:430
size_type ndof(reference_element hat_K) const
Definition basis.h:735
void _sid_grad_initialize(reference_element tilde_L) const
static basis_on_pointset_rep< T > * make_ptr(const std::string &name)
static std::string _make_name(mode_type mode, const std::string &basis_name, const std::string &pointset_name)
const Eigen::Matrix< Value, Eigen::Dynamic, Eigen::Dynamic > & grad_evaluate_on_side(reference_element tilde_L, const side_information_type &sid) const
const basis_basic< T > & get_basis() const
basis_on_pointset_rep< T > & operator=(const basis_on_pointset_rep< T > &)
_RHEOLEF_declare_member(T, _scalar_val) _RHEOLEF_declare_member(point_basic< T >
std::array< bool, reference_element::max_variant > _initialized
void reset(const std::string &name)
static mode_type _parse_name(const std::string &name, std::string &basis_name, std::string &node_name)
void _grad_initialize(reference_element hat_K) const
void _sid_initialize(reference_element tilde_L) const
size_type ndof(reference_element hat_K) const
void _sid_grad_initialize(reference_element tilde_L, const side_information_type &sid) const
size_type nnod(reference_element hat_K) const
void _initialize_continued(reference_element hat_K, const Eigen::Matrix< point_basic< T >, Eigen::Dynamic, 1 > &hat_node) const
void _sid_grad_initialize_continued(reference_element tilde_L, const side_information_type &sid, const Eigen::Matrix< point_basic< T >, Eigen::Dynamic, 1 > &hat_node) const
void _grad_initialize_continued(reference_element hat_K, const Eigen::Matrix< point_basic< T >, Eigen::Dynamic, 1 > &hat_node) const
std::vector< T >::size_type size_type
std::array< bool, reference_element::max_variant > _sid_initialized
const quadrature< T > & get_quadrature() const
std::array< bool, reference_element::max_variant > _grad_initialized
const basis_basic< T > & get_nodal_basis() const
const Eigen::Matrix< Value, Eigen::Dynamic, Eigen::Dynamic > & grad_evaluate(reference_element hat_K) const
void _initialize(reference_element hat_K) const
const Eigen::Matrix< Value, Eigen::Dynamic, Eigen::Dynamic > & evaluate_on_side(reference_element tilde_L, const side_information_type &sid) const
void _sid_initialize_continued(reference_element tilde_L, const side_information_type &sid, const Eigen::Matrix< point_basic< T >, Eigen::Dynamic, 1 > &hat_node) const
const Eigen::Matrix< Value, Eigen::Dynamic, Eigen::Dynamic > & evaluate(reference_element hat_K) const
void set(const quadrature< T > &quad, const basis_basic< T > &b)
const Eigen::Matrix< Value, Eigen::Dynamic, Eigen::Dynamic > & grad_evaluate_on_side(reference_element tilde_L, const side_information_type &sid) const
size_type ndof(reference_element hat_K) const
const basis_basic< T > & get_basis() const
void reset(const std::string &name)
size_type nnod(reference_element hat_K) const
basis_on_pointset_rep< T > rep
const quadrature< T > & get_quadrature() const
const basis_basic< T > & get_nodal_basis() const
const Eigen::Matrix< Value, Eigen::Dynamic, Eigen::Dynamic > & grad_evaluate(reference_element hat_K) const
const Eigen::Matrix< Value, Eigen::Dynamic, Eigen::Dynamic > & evaluate_on_side(reference_element tilde_L, const side_information_type &sid) const
const Eigen::Matrix< Value, Eigen::Dynamic, Eigen::Dynamic > & evaluate(reference_element hat_K) const
see the persistent_table page for the full documentation
see the reference_element page for the full documentation
static const variant_type max_variant
see the smart_pointer page for the full documentation
Expr1::float_type T
Definition field_expr.h:230
This file is part of Rheolef.