Rheolef  7.1
an efficient C++ finite element environment
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
3 
24 /*Class:basis_on_pointset
25 NAME: @code{basis_on_pointset} - pre-evaluated polynomial basis
26 @cindex polynomial basis
27 @clindex basis
28 @cindex reference element
29 @clindex reference_element
30 SYNOPSIS:
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 
46 AUTHOR: Pierre.Saramito@imag.fr
47 DATE: 4 january 2018
48 End:
49 */
50 
51 #include "rheolef/basis.h"
52 #include "rheolef/quadrature.h"
53 namespace rheolef {
54 
55 // -----------------------------------------------------------------------
56 // basis evaluated on lattice of quadrature formulae
57 // -----------------------------------------------------------------------
58 template<class T>
60 public:
61 // typedefs:
62 
64 
65  typedef enum {
66  quad_mode = 0,
69  } mode_type;
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);
123 protected:
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
129 public:
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 
137 _RHEOLEF_declare_member(T,_scalar_val)
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)
161 #undef _RHEOLEF_declare_member
162 
163 protected:
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;
180  void _sid_grad_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:
190 template<class T>
191 class basis_on_pointset: public smart_pointer<basis_on_pointset_rep<T>>,
192  public persistent_table<basis_on_pointset<T>> {
193 public:
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 // -----------------------------------------------------------------------
246 template<class T>
247 inline
248 const basis_basic<T>&
250 {
251  return base::data().get_basis();
252 }
253 template<class T>
254 inline
257 {
258  return base::data().ndof (hat_K);
259 }
260 template<class T>
261 inline
264 {
265  return base::data().nnod (hat_K);
266 }
267 template<class T>
268 inline
269 bool
271 {
272  return base::data().has_quadrature();
273 }
274 template<class T>
275 inline
276 const quadrature<T>&
278 {
279  return base::data().get_quadrature();
280 }
281 template<class T>
282 inline
283 const basis_basic<T>&
285 {
286  return base::data().get_nodal_basis();
287 }
288 template<class T>
289 inline
290 bool
292 {
293  return base::data().is_set();
294 }
295 template<class T>
296 template<class Value>
297 inline
298 const Eigen::Matrix<Value,Eigen::Dynamic,Eigen::Dynamic>&
300 {
301  return base::data().template evaluate<Value> (hat_K);
302 }
303 template<class T>
304 template<class Value>
305 inline
306 const Eigen::Matrix<Value,Eigen::Dynamic,Eigen::Dynamic>&
308 {
309  return base::data().template grad_evaluate<Value> (hat_K);
310 }
311 template<class T>
312 template<class Value>
313 inline
314 const 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 }
321 template<class T>
322 template<class Value>
323 inline
324 const 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
rheolef::basis_on_pointset_rep::get_nodal_basis
const basis_basic< T > & get_nodal_basis() const
Definition: basis_on_pointset.cc:253
rheolef::basis_on_pointset::has_quadrature
bool has_quadrature() const
Definition: basis_on_pointset.h:270
rheolef::basis_on_pointset_rep::_initialize_continued
void _initialize_continued(reference_element hat_K, const Eigen::Matrix< point_basic< T >, Eigen::Dynamic, 1 > &hat_node) const
Definition: basis_on_pointset.cc:263
rheolef::basis_on_pointset_rep::grad_evaluate
const Eigen::Matrix< Value, Eigen::Dynamic, Eigen::Dynamic > & grad_evaluate(reference_element hat_K) const
Definition: basis_on_pointset.cc:520
rheolef::basis_on_pointset_rep::_make_name
static std::string _make_name(mode_type mode, const std::string &basis_name, const std::string &pointset_name)
Definition: basis_on_pointset.cc:32
rheolef::basis_on_pointset_rep::_mode
mode_type _mode
Definition: basis_on_pointset.h:126
rheolef::point_basic
Definition: point.h:87
rheolef::basis_on_pointset::basis_on_pointset
basis_on_pointset(const std::string &name="")
Definition: basis_on_pointset.cc:104
rheolef::tensor4_basic
Definition: tensor4.h:80
rheolef::basis_on_pointset::evaluate
const Eigen::Matrix< Value, Eigen::Dynamic, Eigen::Dynamic > & evaluate(reference_element hat_K) const
Definition: basis_on_pointset.h:299
rheolef::basis_on_pointset::get_nodal_basis
const basis_basic< T > & get_nodal_basis() const
Definition: basis_on_pointset.h:284
rheolef::basis_on_pointset_rep::_initialize
void _initialize(reference_element hat_K) const
Definition: basis_on_pointset.cc:308
rheolef::basis_on_pointset_rep::_initialized
std::array< bool, reference_element::max_variant > _initialized
Definition: basis_on_pointset.h:165
rheolef::basis_on_pointset_rep::quad_mode
@ quad_mode
Definition: basis_on_pointset.h:66
rheolef::basis_on_pointset_rep::max_mode
@ max_mode
Definition: basis_on_pointset.h:68
rheolef::basis_on_pointset_rep::ndof
size_type ndof(reference_element hat_K) const
Definition: basis_on_pointset.cc:240
rheolef::tensor_basic
Definition: tensor.h:90
rheolef::basis_on_pointset_rep::reset
void reset(const std::string &name)
Definition: basis_on_pointset.cc:167
rheolef::basis_on_pointset_rep::_b
basis_basic< T > _b
Definition: basis_on_pointset.h:125
rheolef::basis_on_pointset_rep::has_quadrature
bool has_quadrature() const
Definition: basis_on_pointset.h:89
rheolef::basis_on_pointset::get_basis
const basis_basic< T > & get_basis() const
Definition: basis_on_pointset.h:249
rheolef::side_information_type
Definition: reference_element_face_transformation.h:37
rheolef::basis_basic
Definition: basis.h:206
rheolef::persistent_table
see the persistent_table page for the full documentation
Definition: persistent_table.h:84
rheolef::basis_on_pointset::is_set
bool is_set() const
Definition: basis_on_pointset.h:291
rheolef::reference_element
see the reference_element page for the full documentation
Definition: reference_element.h:66
rheolef::smart_pointer
see the smart_pointer page for the full documentation
Definition: smart_pointer.h:351
rheolef::basis_on_pointset_rep::~basis_on_pointset_rep
~basis_on_pointset_rep()
Definition: basis_on_pointset.cc:87
rheolef::basis_on_pointset::size_type
rep::size_type size_type
Definition: basis_on_pointset.h:197
rheolef::basis_on_pointset_rep::_sid_initialized
std::array< bool, reference_element::max_variant > _sid_initialized
Definition: basis_on_pointset.h:169
rheolef::basis_on_pointset::grad_evaluate_on_side
const Eigen::Matrix< Value, Eigen::Dynamic, Eigen::Dynamic > & grad_evaluate_on_side(reference_element tilde_L, const side_information_type &sid) const
Definition: basis_on_pointset.h:325
rheolef::basis_on_pointset_rep::_grad_initialized
std::array< bool, reference_element::max_variant > _grad_initialized
Definition: basis_on_pointset.h:167
rheolef::basis_on_pointset_rep::_sid_grad_initialize
void _sid_grad_initialize(reference_element tilde_L) const
Definition: basis_on_pointset.cc:502
rheolef::basis_on_pointset::reset
void reset(const std::string &name)
Definition: basis_on_pointset.cc:95
rheolef::basis_on_pointset_rep::evaluate
const Eigen::Matrix< Value, Eigen::Dynamic, Eigen::Dynamic > & evaluate(reference_element hat_K) const
Definition: basis_on_pointset.cc:512
rheolef::basis_on_pointset_rep::evaluate_on_side
const Eigen::Matrix< Value, Eigen::Dynamic, Eigen::Dynamic > & evaluate_on_side(reference_element tilde_L, const side_information_type &sid) const
Definition: basis_on_pointset.cc:528
rheolef::basis_on_pointset::rep
basis_on_pointset_rep< T > rep
Definition: basis_on_pointset.h:195
rheolef::basis_on_pointset_rep::_parse_name
static mode_type _parse_name(const std::string &name, std::string &basis_name, std::string &node_name)
Definition: basis_on_pointset.cc:59
rheolef::basis_on_pointset_rep::_sid_initialize
void _sid_initialize(reference_element tilde_L) const
Definition: basis_on_pointset.cc:481
rheolef::tensor3_basic
Definition: tensor3.h:73
rheolef::basis_on_pointset_rep::_grad_initialize_continued
void _grad_initialize_continued(reference_element hat_K, const Eigen::Matrix< point_basic< T >, Eigen::Dynamic, 1 > &hat_node) const
Definition: basis_on_pointset.cc:287
rheolef::basis_on_pointset_rep::_grad_initialize
void _grad_initialize(reference_element hat_K) const
Definition: basis_on_pointset.cc:322
rheolef::basis_on_pointset_rep::get_basis
const basis_basic< T > & get_basis() const
Definition: basis_on_pointset.h:88
rheolef::basis_on_pointset_rep::_RHEOLEF_declare_member
_RHEOLEF_declare_member(T, _scalar_val) _RHEOLEF_declare_member(point_basic< T >
rheolef
This file is part of Rheolef.
Definition: compiler_eigen.h:37
rheolef::basis_on_pointset::ndof
size_type ndof(reference_element hat_K) const
Definition: basis_on_pointset.h:256
rheolef::basis_on_pointset::base
smart_pointer< rep > base
Definition: basis_on_pointset.h:196
rheolef::reference_element::max_variant
static const variant_type max_variant
Definition: reference_element.h:82
rheolef::basis_on_pointset_rep::_sid_initialize_continued
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
Definition: basis_on_pointset.cc:393
rheolef::basis_on_pointset_rep::size_type
std::vector< T >::size_type size_type
Definition: basis_on_pointset.h:63
rheolef::basis_on_pointset_rep::name
std::string name() const
Definition: basis_on_pointset.cc:47
rheolef::basis_on_pointset_rep::mode_type
mode_type
Definition: basis_on_pointset.h:65
rheolef::basis_on_pointset
Definition: basis_on_pointset.h:191
rheolef::basis_on_pointset::grad_evaluate
const Eigen::Matrix< Value, Eigen::Dynamic, Eigen::Dynamic > & grad_evaluate(reference_element hat_K) const
Definition: basis_on_pointset.h:307
mkgeo_ball.b
b
Definition: mkgeo_ball.sh:152
rheolef::basis_on_pointset::nnod
size_type nnod(reference_element hat_K) const
Definition: basis_on_pointset.h:263
rheolef::basis_on_pointset::set
void set(const quadrature< T > &quad, const basis_basic< T > &b)
Definition: basis_on_pointset.cc:128
size_type
field::size_type size_type
Definition: branch.cc:425
rheolef::basis_on_pointset_rep::nnod
size_type nnod(reference_element hat_K) const
Definition: basis_on_pointset.cc:230
rheolef::basis_on_pointset_rep::operator=
basis_on_pointset_rep< T > & operator=(const basis_on_pointset_rep< T > &)
Definition: basis_on_pointset.cc:206
rheolef::basis_on_pointset_rep::is_set
bool is_set() const
Definition: basis_on_pointset.h:85
rheolef::basis_on_pointset_rep::_quad
quadrature< T > _quad
Definition: basis_on_pointset.h:127
rheolef::basis_on_pointset_rep
Definition: basis_on_pointset.h:59
rheolef::quadrature
Definition: quadrature.h:185
rheolef::basis_on_pointset_rep::_nb
basis_basic< T > _nb
Definition: basis_on_pointset.h:128
mkgeo_contraction.name
name
Definition: mkgeo_contraction.sh:133
rheolef::basis_on_pointset::get_quadrature
const quadrature< T > & get_quadrature() const
Definition: basis_on_pointset.h:277
rheolef::basis_on_pointset_rep::basis_on_pointset_rep
basis_on_pointset_rep(const std::string &name="")
Definition: basis_on_pointset.cc:144
rheolef::basis_on_pointset_rep::grad_evaluate_on_side
const Eigen::Matrix< Value, Eigen::Dynamic, Eigen::Dynamic > & grad_evaluate_on_side(reference_element tilde_L, const side_information_type &sid) const
Definition: basis_on_pointset.cc:538
rheolef::basis_on_pointset_rep::_sid_grad_initialize_continued
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
T
Expr1::float_type T
Definition: field_expr.h:218
rheolef::basis_on_pointset_rep::make_ptr
static basis_on_pointset_rep< T > * make_ptr(const std::string &name)
Definition: basis_on_pointset.cc:82
rheolef::basis_on_pointset_rep::get_quadrature
const quadrature< T > & get_quadrature() const
Definition: basis_on_pointset.cc:246
rheolef::basis_on_pointset_rep::nodal_mode
@ nodal_mode
Definition: basis_on_pointset.h:67
rheolef::basis_on_pointset::evaluate_on_side
const Eigen::Matrix< Value, Eigen::Dynamic, Eigen::Dynamic > & evaluate_on_side(reference_element tilde_L, const side_information_type &sid) const
Definition: basis_on_pointset.h:315