Rheolef  7.1
an efficient C++ finite element environment
basis_fem_vector.cc
Go to the documentation of this file.
1 //
22 // vectors
23 //
24 #include "basis_fem_vector.h"
25 #include "rheolef/rheostream.h"
26 
27 namespace rheolef {
28 
29 using namespace std;
30 
31 // =========================================================================
32 // basis members
33 // =========================================================================
34 template<class T>
36 {
37 }
38 template<class T>
40  : basis_rep<T>(sopt),
41  _n_comp(0),
42  _scalar_basis(scalar_basis),
43  _scalar_value(),
44  _vector_value()
45 {
48  base::_piola_fem = _scalar_basis.get_piola_fem();
49  check_macro (base::option().dimension() != std::numeric_limits<basis_option::size_type>::max(),
50  "vector(basis): basis.option.map_dimension should be initialized for component number");
53 }
54 #ifdef TO_CLEAN
55 template<class T>
58 {
59  const size_type unset = std::numeric_limits<basis_option::size_type>::max();
60  return (base::option().dimension() == unset) ? map_d : base::option().dimension();
61 }
62 #endif // TO_CLEAN
63 template<class T>
64 void
66 {
67  for (size_type map_d = 0; map_d < 4; ++map_d) {
68  for (size_type subgeo_variant = 0; subgeo_variant < reference_element::max_variant; ++subgeo_variant) {
69  base::_ndof_on_subgeo [map_d][subgeo_variant] = _n_comp*_scalar_basis.ndof_on_subgeo (map_d, subgeo_variant);
70  base::_nnod_on_subgeo [map_d][subgeo_variant] = _scalar_basis.nnod_on_subgeo (map_d, subgeo_variant);
71  }
72  }
74  reference_element hat_K (variant);
75  for (size_type subgeo_d = 0; subgeo_d < 5; ++subgeo_d) {
76  base::_first_idof_by_dimension [variant][subgeo_d] = _n_comp*_scalar_basis.first_idof_by_dimension (hat_K, subgeo_d);
77  base::_first_inod_by_dimension [variant][subgeo_d] = _scalar_basis.first_inod_by_dimension (hat_K, subgeo_d);
78  }
79  }
80 }
81 template<class T>
82 void
84 {
85 }
86 template<class T>
87 const Eigen::Matrix<point_basic<T>,Eigen::Dynamic,1>&
89 {
90  return _scalar_basis.hat_node (hat_K);
91 }
92 template<class T>
93 void
95  reference_element hat_K,
96  const point_basic<T>& hat_x,
97  Eigen::Matrix<point_basic<T>,Eigen::Dynamic,1>& value) const
98 {
99  base::_initialize_data_guard (hat_K);
100  _scalar_basis.evaluate (hat_K, hat_x, _scalar_value);
101  size_type loc_comp_ndof = _scalar_value.size();
102  size_type loc_ndof = _n_comp*loc_comp_ndof;
103  value.resize (loc_ndof);
104  value.fill (point_basic<T>()); // do not remove !
105  for (size_type loc_comp_idof = 0; loc_comp_idof < loc_comp_ndof; ++loc_comp_idof) {
106  for (size_type i_comp = 0; i_comp < _n_comp; ++i_comp) {
107  size_type loc_idof = _n_comp*loc_comp_idof + i_comp;
108  value[loc_idof][i_comp] = _scalar_value[loc_comp_idof];
109  }
110  }
111 }
112 template<class T>
113 void
115  reference_element hat_K,
116  const point_basic<T>& hat_x,
117  Eigen::Matrix<tensor_basic<T>,Eigen::Dynamic,1>& value) const
118 {
119  base::_initialize_data_guard (hat_K);
120  _scalar_basis.grad_evaluate (hat_K, hat_x, _vector_value);
121  size_type loc_comp_ndof = _vector_value.size();
122  size_type loc_ndof = _n_comp*loc_comp_ndof;
123  value.resize (loc_ndof);
124  value.fill (tensor_basic<T>()); // do not remove !
125  for (size_type loc_comp_idof = 0; loc_comp_idof < loc_comp_ndof; ++loc_comp_idof) {
126  for (size_type i_comp = 0; i_comp < _n_comp; ++i_comp) {
127  size_type loc_idof = _n_comp*loc_comp_idof + i_comp;
128  for (size_type j_comp = 0; j_comp < _n_comp; ++j_comp) {
129  value[loc_idof] (i_comp,j_comp) = _vector_value[loc_comp_idof][j_comp];
130  }
131  }
132  }
133 }
134 template<class T>
135 void
137  reference_element hat_K,
138  const Eigen::Matrix<point_basic<T>,Eigen::Dynamic,1>& f_xnod,
139  Eigen::Matrix<T,Eigen::Dynamic,1>& dof) const
140 {
141  base::_initialize_data_guard (hat_K);
142  size_type loc_comp_ndof = _scalar_basis.ndof (hat_K);
143  size_type loc_comp_nnod = _scalar_basis.nnod (hat_K);
144  size_type loc_ndof = _n_comp*loc_comp_ndof;
145  Eigen::Matrix<T,Eigen::Dynamic,1> f_comp_xnod (loc_comp_nnod); // TODO: class working array
146  Eigen::Matrix<T,Eigen::Dynamic,1> comp_dof (loc_comp_ndof); // TODO: class working array
147  dof.resize (loc_ndof);
148  for (size_type i_comp = 0; i_comp < _n_comp; ++i_comp) {
149  for (size_type loc_comp_inod = 0; loc_comp_inod < loc_comp_nnod; ++loc_comp_inod) {
150  f_comp_xnod [loc_comp_inod] = f_xnod [loc_comp_inod] [i_comp];
151  }
152  _scalar_basis.compute_dofs (hat_K, f_comp_xnod, comp_dof);
153  for (size_type loc_comp_idof = 0; loc_comp_idof < loc_comp_ndof; ++loc_comp_idof) {
154  size_type loc_idof = _n_comp*loc_comp_idof + i_comp;
155  dof [loc_idof] = comp_dof [loc_comp_idof];
156  }
157  }
158 }
159 // ----------------------------------------------------------------------------
160 // instanciation in library
161 // ----------------------------------------------------------------------------
162 #define _RHEOLEF_instanciation(T) \
163 template class basis_fem_vector<T>;
164 
166 
167 }// namespace rheolef
field::size_type size_type
Definition: branch.cc:425
see the Float page for the full documentation
basis_fem_vector(const basis_basic< T > &scalar_basis, const basis_option &sopt)
void grad_evaluate(reference_element hat_K, const point_basic< T > &hat_x, Eigen::Matrix< tensor_basic< T >, Eigen::Dynamic, 1 > &value) const
void _initialize_cstor_sizes() const
basis_basic< T > _scalar_basis
void evaluate(reference_element hat_K, const point_basic< T > &hat_x, Eigen::Matrix< point_basic< T >, Eigen::Dynamic, 1 > &value) const
void _compute_dofs(reference_element hat_K, const Eigen::Matrix< point_basic< T >, Eigen::Dynamic, 1 > &f_xnod, Eigen::Matrix< T, Eigen::Dynamic, 1 > &dof) const
size_type family_index() const
const Eigen::Matrix< point_basic< T >, Eigen::Dynamic, 1 > & hat_node(reference_element hat_K) const
std::string family_name() const
void _initialize_data(reference_element hat_K) const
see the basis_option page for the full documentation
Definition: basis_option.h:93
size_type dimension() const
Definition: basis_option.h:138
void set_valued_tag(valued_type v)
Definition: basis_option.h:149
reference_element::size_type size_type
Definition: basis.h:214
piola_fem< T > _piola_fem
Definition: basis.h:394
basis_option _sopt
Definition: basis.h:393
const basis_option & option() const
Definition: basis.h:238
static std::string standard_naming(std::string family_name, size_t degree, const basis_option &sopt)
Definition: basis_rep.cc:44
std::string _name
Definition: basis.h:392
see the reference_element page for the full documentation
static const variant_type max_variant
rheolef::std value
const size_t dimension
Definition: edge.icc:64
Expr1::float_type T
Definition: field_expr.h:261
check_macro(expr1.have_homogeneous_space(Xh1), "dual(expr1,expr2); expr1 should have homogeneous space. HINT: use dual(interpolate(Xh, expr1),expr2)")
size_type n_component(valued_type valued_tag, size_type d, coordinate_type sys_coord)
This file is part of Rheolef.
_RHEOLEF_instanciation(Float) _RHEOLEF_instanciation_evaluate(Float