1#ifndef _RHEO_BASIS_ON_POINTSET_V2_H
2#define _RHEO_BASIS_ON_POINTSET_V2_H
51#include "rheolef/basis.h"
52#include "rheolef/quadrature.h"
63 typedef typename std::vector<T>::size_type
size_type;
84 std::string
name()
const;
94 const Eigen::Matrix<Value,Eigen::Dynamic,Eigen::Dynamic>&
98 const Eigen::Matrix<Value,Eigen::Dynamic,Eigen::Dynamic>&
103 template<
class Value>
104 const Eigen::Matrix<Value,Eigen::Dynamic,Eigen::Dynamic>&
107 template<
class Value>
108 const Eigen::Matrix<Value,Eigen::Dynamic,Eigen::Dynamic>&
117 const std::string& basis_name,
118 const std::string& pointset_name);
120 const std::string&
name,
121 std::string& basis_name,
122 std::string& node_name);
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; \
142#undef _RHEOLEF_declare_member
145#define _RHEOLEF_declare_member(VALUE,MEMBER) \
146 mutable std::array< \
150 Eigen::Matrix<VALUE,Eigen::Dynamic,Eigen::Dynamic>, \
154 reference_element::max_variant> MEMBER;
161#undef _RHEOLEF_declare_member
164 mutable std::array<bool,
166 mutable std::array<bool,
168 mutable std::array<bool,
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;
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;
207 void reset (
const std::string& name);
221 template<
class Value>
222 const Eigen::Matrix<Value,Eigen::Dynamic,Eigen::Dynamic>&
225 template<
class Value>
226 const Eigen::Matrix<Value,Eigen::Dynamic,Eigen::Dynamic>&
231 template<
class Value>
232 const Eigen::Matrix<Value,Eigen::Dynamic,Eigen::Dynamic>&
235 template<
class Value>
236 const Eigen::Matrix<Value,Eigen::Dynamic,Eigen::Dynamic>&
251 return base::data().get_basis();
258 return base::data().
ndof (hat_K);
265 return base::data().nnod (hat_K);
272 return base::data().has_quadrature();
279 return base::data().get_quadrature();
286 return base::data().get_nodal_basis();
293 return base::data().is_set();
298const Eigen::Matrix<Value,Eigen::Dynamic,Eigen::Dynamic>&
301 return base::data().template evaluate<Value> (hat_K);
306const Eigen::Matrix<Value,Eigen::Dynamic,Eigen::Dynamic>&
309 return base::data().template grad_evaluate<Value> (hat_K);
314const Eigen::Matrix<Value,Eigen::Dynamic,Eigen::Dynamic>&
319 return base::data().template evaluate_on_side<Value> (tilde_L, sid);
324const Eigen::Matrix<Value,Eigen::Dynamic,Eigen::Dynamic>&
329 return base::data().template grad_evaluate_on_side<Value> (tilde_L, sid);
field::size_type size_type
size_type ndof(reference_element hat_K) const
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 >
bool has_quadrature() const
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
bool has_quadrature() const
void reset(const std::string &name)
size_type nnod(reference_element hat_K) const
smart_pointer< rep > base
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
This file is part of Rheolef.