Rheolef  7.1
an efficient C++ finite element environment
basis_on_pointset.cc
Go to the documentation of this file.
1 #include "rheolef/basis_on_pointset.h"
22 #include "rheolef/memorized_value.h"
23 #include "rheolef/pretty_name.h"
25 namespace rheolef {
26 
27 // ----------------------------------------------------------
28 // name
29 // ----------------------------------------------------------
30 template<class T>
31 std::string
33  mode_type mode,
34  const std::string& basis_name,
35  const std::string& pointset_name)
36 {
37  if (mode == max_mode) {
38  return "";
39  } else if (mode == quad_mode) {
40  return basis_name + "@q(" + pointset_name + ")";
41  } else {
42  return basis_name + "@b(" + pointset_name + ")";
43  }
44 }
45 template<class T>
46 std::string
48 {
49  if (_mode == max_mode) {
50  return "";
51  } else if (_mode == quad_mode) {
52  return _make_name (_mode, _b.name(), _quad.name());
53  } else {
54  return _make_name (_mode, _b.name(), _nb.name());
55  }
56 }
57 template<class T>
60  const std::string& name,
61  std::string& basis_name,
62  std::string& pointset_name)
63 {
64  long i_pos_at = name.find ('@');
65  check_macro (i_pos_at > 0 && i_pos_at+1 < long(name.size()),
66  "invalid basis_on_pointset name, expect e.g. \"P1@b(gauss(3))\"");
67  basis_name = name.substr (0, i_pos_at);
68  check_macro (name[i_pos_at+1] == 'q' || name[i_pos_at+1] == 'b',
69  "invalid basis_on_pointset name, expect e.g. \"P1@b(gauss(3))\"");
70  mode_type mode = name[i_pos_at+1] == 'q' ? quad_mode : nodal_mode;
71  long l = name.size() - (i_pos_at+4);
72  check_macro (l > 0 && i_pos_at+3 < long(name.size()),
73  "invalid basis_on_pointset name, expect e.g. \"P1@b(gauss(3))\"");
74  pointset_name = name.substr (i_pos_at+3, l);
75  return mode;
76 }
77 // ----------------------------------------------------------
78 // persistent table
79 // ----------------------------------------------------------
80 template<class T>
83 {
84  return new_macro(basis_on_pointset_rep<T>(name));
85 }
86 template <class T>
88 {
89  if (name() != "") {
91  }
92 }
93 template<class T>
94 void
95 basis_on_pointset<T>::reset (const std::string& name)
96 {
97  if (name == "") {
98  base::operator= (new_macro(rep()));
99  } else {
100  base::operator= (persistent_table<basis_on_pointset<T>>::load (name));
101  }
102 }
103 template<class T>
105  : base(new_macro(rep)),
107 {
108  reset (name);
109 }
110 template<class T>
112  : base(),
114 {
116  reset (name);
117 }
118 template<class T>
120  : base(),
122 {
124  reset (name);
125 }
126 template<class T>
127 void
129 {
131  reset (name);
132 }
133 template<class T>
134 void
136 {
138  reset (name);
139 }
140 // ----------------------------------------------------------
141 // cstors
142 // ----------------------------------------------------------
143 template<class T>
145  : _b(),
146  _mode(max_mode),
147  _quad(),
148  _nb(),
149  _scalar_val(),
150  _vector_val(),
151  _tensor_val(),
152  _tensor3_val(),
153  _tensor4_val(),
154  _sid_scalar_val(),
155  _sid_vector_val(),
156  _sid_tensor_val(),
157  _sid_tensor3_val(),
158  _sid_tensor4_val(),
159  _initialized(),
160  _grad_initialized(),
161  _sid_initialized()
162 {
163  reset (name);
164 }
165 template<class T>
166 void
168 {
169  _initialized.fill(false);
170  _grad_initialized.fill(false);
171  _sid_initialized.fill(false);
172  if (name == "") return;
173  std::string basis_name, pointset_name;
174  _mode = _parse_name (name, basis_name, pointset_name);
175  _b = basis_basic<T> (basis_name);
176  if (_mode == quad_mode) {
177  _quad = quadrature<T> (pointset_name);
178  } else {
179  _nb = basis_basic<T> (pointset_name);
180  }
181 }
182 template<class T>
184  : _b(x._b),
185  _mode(x._mode),
186  _quad(x._quad),
187  _nb(x._nb),
188  _scalar_val(x._scalar_val),
189  _vector_val(x._vector_val),
190  _tensor_val(x._tensor_val),
191  _tensor3_val(x._tensor3_val),
192  _tensor4_val(x._tensor4_val),
193  _sid_scalar_val(x._sid_scalar_val),
194  _sid_vector_val(x._sid_vector_val),
195  _sid_tensor_val(x._sid_tensor_val),
196  _sid_tensor3_val(x._sid_tensor3_val),
197  _sid_tensor4_val(x._sid_tensor4_val),
198  _initialized (x._initialized),
199  _grad_initialized(x._grad_initialized),
200  _sid_initialized (x._sid_initialized)
201 {
202  trace_macro("PHYSICAL COPY***************");
203 }
204 template<class T>
207 {
208  trace_macro("ASSIGN**************");
209  _b = x._b;
210  _mode = x._mode;
211  _quad = x._quad;
212  _nb = x._nb;
213  _scalar_val = x._scalar_val;
214  _vector_val = x._vector_val;
215  _tensor_val = x._tensor_val;
216  _tensor3_val = x._tensor3_val;
217  _tensor4_val = x._tensor4_val;
218  _sid_scalar_val = x._sid_scalar_val;
219  _sid_vector_val = x._sid_vector_val;
220  _sid_tensor_val = x._sid_tensor_val;
221  _sid_tensor3_val = x._sid_tensor3_val;
222  _sid_tensor4_val = x._sid_tensor4_val;
223  _initialized = x._initialized;
224  _grad_initialized = x._grad_initialized;
225  _sid_initialized = x._sid_initialized;
226  return *this;
227 }
228 template<class T>
231 {
232  if (_mode == quad_mode) {
233  return _quad.size (hat_K);
234  } else {
235  return _nb.nnod (hat_K);
236  }
237 }
238 template<class T>
241 {
242  return _b.ndof (hat_K);
243 }
244 template<class T>
245 const quadrature<T>&
247 {
248  check_macro (_mode == quad_mode, "get_quadrature: pointset mode is not quadrature");
249  return _quad;
250 }
251 template<class T>
252 const basis_basic<T>&
254 {
255  check_macro (_mode == nodal_mode, "get_nodal_basis: pointset mode is not nodal");
256  return _nb;
257 }
258 // -----------------------------------------------------------------------
259 // basis evaluated on lattice of quadrature formulae
260 // -----------------------------------------------------------------------
261 template<class T>
262 void
264  reference_element hat_K,
265  const Eigen::Matrix<point_basic<T>,Eigen::Dynamic,1>& hat_node) const
266 {
267  switch (_b.valued_tag()) {
268  case space_constant::scalar: {
269  details::basis_on_pointset_evaluate (_b, hat_K, hat_node, _scalar_val[hat_K.variant()]);
270  break;
271  }
272  case space_constant::vector: {
273  details::basis_on_pointset_evaluate (_b, hat_K, hat_node, _vector_val[hat_K.variant()]);
274  break;
275  }
276  case space_constant::tensor: {
277  // note: not yet available in basis ; for specific tensor basis, e.g. winter
278  details::basis_on_pointset_evaluate (_b, hat_K, hat_node, _tensor_val[hat_K.variant()]);
279  break;
280  }
281  default: error_macro("evaluate: unexpected "<<space_constant::valued_name(_b.valued_tag())
282  << "-valued basis \"" << _b.name() << "\"");
283  }
284 }
285 template<class T>
286 void
288  reference_element hat_K,
289  const Eigen::Matrix<point_basic<T>,Eigen::Dynamic,1>& hat_node) const
290 {
291  switch (_b.valued_tag()) {
293  details::basis_on_pointset_grad_evaluate (_b, hat_K, hat_node, _vector_val[hat_K.variant()]);
294  break;
296  details::basis_on_pointset_grad_evaluate (_b, hat_K, hat_node, _tensor_val[hat_K.variant()]);
297  break;
299  details::basis_on_pointset_grad_evaluate (_b, hat_K, hat_node, _tensor3_val[hat_K.variant()]);
300  break;
301  default: error_macro("grad_evaluate: unexpected "<<space_constant::valued_name(_b.valued_tag())
302  << "-valued basis \"" << _b.name() << "\"");
303  }
304 }
305 // TODO: avoid code repetition between initialize & grad_initialize
306 template<class T>
307 void
309 {
310  if (_mode == quad_mode) {
311  Eigen::Matrix<point_basic<T>,Eigen::Dynamic,1> hat_node;
312  _quad.get_nodes (hat_K, hat_node);
313  _initialize_continued (hat_K, hat_node);
314  } else {
315  const Eigen::Matrix<point_basic<T>,Eigen::Dynamic,1>& hat_node = _nb.hat_node (hat_K);
316  _initialize_continued (hat_K, hat_node);
317  }
318  _initialized [hat_K.variant()] = true;
319 }
320 template<class T>
321 void
323 {
324  reference_element::variant_type K_variant = hat_K.variant();
325  if (_mode == quad_mode) {
326  Eigen::Matrix<point_basic<T>,Eigen::Dynamic,1> hat_node;
327  _quad.get_nodes (hat_K, hat_node);
328  _grad_initialize_continued (hat_K, hat_node);
329  } else {
330  const Eigen::Matrix<point_basic<T>,Eigen::Dynamic,1>& hat_node = _nb.hat_node (hat_K);
331  _grad_initialize_continued (hat_K, hat_node);
332  }
333  _grad_initialized [K_variant] = true;
334 }
335 // ----------------------------------------------------------
336 // side restriction (for DG)
337 // ----------------------------------------------------------
338 // We have a side tilde_K on an element tilde_L
339 // and a Piola transform from hat_K to tilde_L
340 // We have a quadrature formula or a point-set on hat_K:
341 // (hat_xq) q=0..nq
342 // and we are able to transform it on tilde_K subset partial tilde_L:
343 // (tilde_xq) q=0..nq
344 // We want to evaluate on this pointset a basis "b" of tilde_L.
345 //
346 // basis_on_pointset<T>::init_all_sides (reference_element tilde_L)
347 // => loop on all possible tilde_K side of tilde_L
348 // build the (tilde_xq)q set and evaluate on it the basis "b"
349 // (possibly, have to reverse or reorder (ie rotate in 3d) the tilde_xq set)
350 // store all results contiguously, side by side, in the same array
351 //
352 // basis_on_pointset<T>::evaluate_on_side (reference_element tilde_L, const side_information_type& sid)
353 // => if not initialized, call init_all_sides
354 // then, from the current side index, compute begin and end pointers
355 // that will be used by accessors
356 //
357 // AVANTAGES:
358 // - transparent : no changes in downstream codes (field_evaluate.cc etc)
359 // - fast : the basis is evaluate only once for each reference element
360 // TECHNIQUE:
361 // possibly, have to reverse or reorder (ie rotate in 3d) the tilde_xq set
362 // - first, do it in 2d, with only an optional reverse ordering
363 // when the reference side has the opposite orientation.
364 // This is not known at "init_all_sides()" call, as it is given later
365 // by the "side_information_type" argument of "restrict_on_side()" call.
366 // So, we cannot use begin() and end() iterators.
367 // Thus, the random accessor is preferable: "value()"
368 // and any call to iterators should check that we are not in the "side" mode.
369 // - second, extend it to 3d: the reordering is more complex when there is shift
370 // * nodes on the vertices are rotated
371 // * nodes on the boundary edges of a face are shifted: iedge=(iedge0+shift)%nedge
372 // * nodes inside the face are completely reordered
373 // => a table of reordering is required, it should be build by restrict_on_side()
374 // Note that there is a finite number of possible shift+reverse combinations
375 // and we can enumerate it and compute one time for all the all reordering
376 // at init_all_sides(). Then, at restrict_on_side(), we have just to compute
377 // the index associated to the reordering.
378 // IMPLEMENTATION:
379 // - may have a different data structure for eval on volume or side
380 // _scalar_val(), _vector_val(),
381 // _sid_scalar_val(), _sid_vector_val(),
382 // or with a vol/side flag:
383 // mutable std::array<2, array<Eigen::Matrix<T,Eigen::Dynamic,1>,
384 // reference_element::max_variant> > _scalar_val;
385 // =>
386 // _scalar_val [vol_side_flag] [hat_K] [idof]
387 // motivation: when the test::_bops is used for volume, then for side,
388 // and then volume again, we do not have to reinitialize.
389 // perhaps, it should solve the BUG_TEST_CONST ?
390 
391 template<class T>
392 void
394  reference_element tilde_L,
395  const side_information_type& sid,
396  const Eigen::Matrix<point_basic<T>,Eigen::Dynamic,1>& hat_node) const
397 {
398  reference_element::variant_type L_variant = tilde_L.variant();
399  size_type ori_idx = (sid.orient == 1) ? 0 : 1;
400  switch (_b.valued_tag()) {
401  case space_constant::scalar: {
402  if (!_b.option().is_restricted_to_sides()) {
403  details::basis_on_pointset_evaluate (_b, tilde_L, hat_node, _sid_scalar_val[L_variant][sid.loc_isid][ori_idx][sid.shift]);
404  details::basis_on_pointset_grad_evaluate (_b, tilde_L, hat_node, _sid_vector_val[L_variant][sid.loc_isid][ori_idx][sid.shift]);
405  } else { // e.g. "Pkd[sides]"
406  details::basis_on_pointset_evaluate_on_side (_b, tilde_L, sid, hat_node, _sid_scalar_val[L_variant][sid.loc_isid][ori_idx][sid.shift]);
407 #ifdef TODO
408 #endif // TODO
409  }
410  break;
411  }
412  case space_constant::vector: {
413  if (!_b.option().is_restricted_to_sides()) {
414  details::basis_on_pointset_evaluate (_b, tilde_L, hat_node, _sid_vector_val[L_variant][sid.loc_isid][ori_idx][sid.shift]);
415  details::basis_on_pointset_grad_evaluate (_b, tilde_L, hat_node, _sid_tensor_val[L_variant][sid.loc_isid][ori_idx][sid.shift]);
416  } else {
417  fatal_macro ("basis_sides: not yet");
418 #ifdef TODO
419  details::basis_on_pointset_evaluate_on_sides (_b, tilde_L, hat_node, _sid_vector_val[L_variant][sid.loc_isid][ori_idx][sid.shift]);
420 #endif // TODO
421  }
422  break;
423  }
424  case space_constant::tensor: {
425  if (!_b.option().is_restricted_to_sides()) {
426  details::basis_on_pointset_evaluate (_b, tilde_L, hat_node, _sid_tensor_val[L_variant][sid.loc_isid][ori_idx][sid.shift]);
427  details::basis_on_pointset_grad_evaluate (_b, tilde_L, hat_node, _sid_tensor3_val[L_variant][sid.loc_isid][ori_idx][sid.shift]);
428  } else {
429  fatal_macro ("basis_sides: not yet");
430 #ifdef TODO
431  details::basis_on_pointset_evaluate_on_sides (_b, tilde_L, hat_node, _sid_tensor_val[L_variant][sid.loc_isid][ori_idx][sid.shift]);
432 #endif // TODO
433  }
434  break;
435  }
437  if (!_b.option().is_restricted_to_sides()) {
438  details::basis_on_pointset_evaluate (_b, tilde_L, hat_node, _sid_tensor3_val[L_variant][sid.loc_isid][ori_idx][sid.shift]);
439  details::basis_on_pointset_grad_evaluate (_b, tilde_L, hat_node, _sid_tensor4_val[L_variant][sid.loc_isid][ori_idx][sid.shift]);
440  } else {
441  fatal_macro ("basis_sides: not yet");
442 #ifdef TODO
443  details::basis_on_pointset_evaluate_on_sides (_b, tilde_L, hat_node, _sid_tensor_val[L_variant][sid.loc_isid][ori_idx][sid.shift]);
444 #endif // TODO
445  }
446  break;
447  }
448  default: error_macro("evaluate_on_sides: unexpected "<<space_constant::valued_name(_b.valued_tag())
449  << "-valued basis \"" << _b.name() << "\"");
450  }
451 }
452 // TODO: avoid code repetition between initialize & grad_initialize
453 template<class T>
454 void
456 {
457  // transform side hat_node on hat_K into tilde_node on tilde_K, the side of tilde_L:
458  // => orient and shift varies for each side
459  // the transformed nodes from hat_K to tilde_K depends upon (loc_isid,orient,shift)
460  reference_element::variant_type L_variant = tilde_L.variant();
461  Eigen::Matrix<point_basic<T>,Eigen::Dynamic,1> hat_node, tilde_node;
462  if (_mode == quad_mode) {
463  _quad.get_nodes (sid.hat, hat_node);
464  } else {
465  hat_node = _nb.hat_node (sid.hat); // ... on its Lagrange nodes.
466  }
467  tilde_node.resize (hat_node.size());
468  if (!_b.option().is_restricted_to_sides()) {
469  // transform hat_xq on side hat_S to tilde_xq on partial tilde_L
470  for (size_type q = 0, nq = hat_node.size(); q < nq; ++q) {
471  tilde_node[q] = reference_element_face_transformation (tilde_L, sid, hat_node[q]);
472  }
473  } else {
474  // else e.g. "Pkd[sides]": hat_nodes remains on side hat_S
475  tilde_node = hat_node;
476  }
477  _sid_initialize_continued (tilde_L, sid, tilde_node);
478 }
479 template<class T>
480 void
482 {
483  size_type d = tilde_L.dimension();
484  size_type nsid = tilde_L.n_subgeo (d-1);
486  sid.dim = d-1;
487  for (sid.loc_isid = 0; sid.loc_isid < nsid; ++sid.loc_isid) {
488  sid.n_vertex = tilde_L.subgeo_size (sid.dim, sid.loc_isid);
489  sid.hat.set_variant (sid.n_vertex, sid.dim);
490  int last_orient = (sid.dim == 0) ? 1 : -1; // no negative orient for 0d side
491  for (sid.orient = 1; sid.orient >= last_orient; sid.orient -= 2) {
492  size_type n_shift = (d == 3) ? nsid : 1;
493  for (sid.shift = 0; sid.shift < n_shift; ++sid.shift) {
494  _sid_initialize (tilde_L, sid);
495  }
496  }
497  }
498  _sid_initialized [tilde_L.variant()] = true;
499 }
500 template<class T>
501 void
503 {
504  fatal_macro ("_sid_grad_initialize: not yet");
505 }
506 // ----------------------------------------------------------
507 // accessors
508 // ----------------------------------------------------------
509 template<class T>
510 template<class Value>
511 const Eigen::Matrix<Value,Eigen::Dynamic,Eigen::Dynamic>&
513 {
514  if (!_initialized [hat_K.variant()]) _initialize (hat_K);
515  return details::memorized_matrix<T,Value>().get (*this, hat_K);
516 }
517 template<class T>
518 template<class Value>
519 const Eigen::Matrix<Value,Eigen::Dynamic,Eigen::Dynamic>&
521 {
522  if (!_grad_initialized [hat_K.variant()]) _grad_initialize (hat_K);
523  return details::memorized_matrix<T,Value>().get (*this, hat_K);
524 }
525 template<class T>
526 template<class Value>
527 const Eigen::Matrix<Value,Eigen::Dynamic,Eigen::Dynamic>&
529  reference_element hat_K,
530  const side_information_type& sid) const
531 {
532  if (!_sid_initialized [hat_K.variant()]) _sid_initialize (hat_K);
533  return details::memorized_side_value<T,Value>().get (*this, hat_K, sid);
534 }
535 template<class T>
536 template<class Value>
537 const Eigen::Matrix<Value,Eigen::Dynamic,Eigen::Dynamic>&
539  reference_element hat_K,
540  const side_information_type& sid) const
541 {
542  if (!_sid_initialized [hat_K.variant()]) _sid_initialize (hat_K);
543  return details::memorized_side_value<T,Value>().get (*this, hat_K, sid);
544 }
545 // ----------------------------------------------------------------------------
546 // instanciation in library
547 // ----------------------------------------------------------------------------
548 #define _RHEOLEF_instanciation_value(T,Value) \
549 template \
550 const Eigen::Matrix<Value,Eigen::Dynamic,Eigen::Dynamic>& \
551 basis_on_pointset_rep<T>::evaluate (reference_element hat_K) const; \
552 template \
553 const Eigen::Matrix<Value,Eigen::Dynamic,Eigen::Dynamic>& \
554 basis_on_pointset_rep<T>::evaluate_on_side ( \
555  reference_element hat_K, \
556  const side_information_type& sid) const; \
557 template \
558 const Eigen::Matrix<Value,Eigen::Dynamic,Eigen::Dynamic>& \
559 basis_on_pointset_rep<T>::grad_evaluate (reference_element hat_K) const; \
560 template \
561 const Eigen::Matrix<Value,Eigen::Dynamic,Eigen::Dynamic>& \
562 basis_on_pointset_rep<T>::grad_evaluate_on_side ( \
563  reference_element hat_K, \
564  const side_information_type& sid) const; \
565 
566 #define _RHEOLEF_instanciation(T) \
567 template class basis_on_pointset<T>; \
568 template class basis_on_pointset_rep<T>; \
569 _RHEOLEF_instanciation_value(T,T) \
570 _RHEOLEF_instanciation_value(T,point_basic<T>) \
571 _RHEOLEF_instanciation_value(T,tensor_basic<T>) \
572 _RHEOLEF_instanciation_value(T,tensor3_basic<T>) \
573 _RHEOLEF_instanciation_value(T,tensor4_basic<T>) \
574 
576 
577 } // namespace rheolef
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_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
check_macro
check_macro(expr1.have_homogeneous_space(Xh1), "dual(expr1,expr2); expr1 should have homogeneous space. HINT: use dual(interpolate(Xh, expr1),expr2)")
rheolef::details::basis_on_pointset_evaluate_on_side
void basis_on_pointset_evaluate_on_side(const Basis &b, const reference_element &tilde_K, const side_information_type &sid, const Eigen::Matrix< point_basic< T >, Eigen::Dynamic, 1 > &hat_x, Eigen::Matrix< Value, Eigen::Dynamic, Eigen::Dynamic > &vdm)
Definition: basis_on_pointset_evaluate.icc:65
rheolef::details::basis_on_pointset_evaluate
void basis_on_pointset_evaluate(const Basis &b, const reference_element &hat_K, const Eigen::Matrix< point_basic< T >, Eigen::Dynamic, 1 > &hat_x, Eigen::Matrix< Value, Eigen::Dynamic, Eigen::Dynamic > &vdm)
Definition: basis_on_pointset_evaluate.icc:31
rheolef::details::reset
void reset(T &x)
Definition: keller_details.h:31
rheolef::reference_element_face_transformation
point_basic< T > reference_element_face_transformation(reference_element tilde_K, const side_information_type &sid, const point_basic< T > &sid_hat_x)
Definition: reference_element_face_transformation.cc:33
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::_RHEOLEF_instanciation
_RHEOLEF_instanciation(Float, sequential, std::allocator< Float >) _RHEOLEF_instanciation(Float
rheolef::reference_element::subgeo_size
size_type subgeo_size(size_type subgeo_dim, size_type loc_isid) const
Definition: reference_element.h:120
rheolef::reference_element::n_subgeo
size_type n_subgeo(size_type subgeo_dim) const
Definition: reference_element.h:110
rheolef::basis_on_pointset_rep::ndof
size_type ndof(reference_element hat_K) const
Definition: basis_on_pointset.cc:240
rheolef::space_constant::tensor3
Definition: space_constant.h:140
rheolef::side_information_type::n_vertex
size_t n_vertex
Definition: reference_element_face_transformation.h:42
rheolef::basis_on_pointset_rep::reset
void reset(const std::string &name)
Definition: basis_on_pointset.cc:167
rheolef::side_information_type::dim
size_t dim
Definition: reference_element_face_transformation.h:41
rheolef::side_information_type::orient
geo_element_indirect::orientation_type orient
Definition: reference_element_face_transformation.h:40
rheolef::basis_on_pointset_rep::_b
basis_basic< T > _b
Definition: basis_on_pointset.h:125
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::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::space_constant::tensor
Definition: space_constant.h:138
rheolef::space_constant::scalar
Definition: space_constant.h:136
rheolef::basis_on_pointset_rep::~basis_on_pointset_rep
~basis_on_pointset_rep()
Definition: basis_on_pointset.cc:87
rheolef::quadrature::name
std::string name() const
Definition: quadrature.h:213
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_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::reference_element::variant_type
size_type variant_type
Definition: reference_element.h:72
rheolef::basis_on_pointset::reset
void reset(const std::string &name)
Definition: basis_on_pointset.cc:95
rheolef::basis_basic::name
std::string name() const
Definition: basis.h:682
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::details::basis_on_pointset_grad_evaluate
void basis_on_pointset_grad_evaluate(const Basis &b, const reference_element &hat_K, const Eigen::Matrix< point_basic< T >, Eigen::Dynamic, 1 > &hat_x, Eigen::Matrix< Value, Eigen::Dynamic, Eigen::Dynamic > &vdm_grad)
Definition: basis_on_pointset_evaluate.icc:48
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::reference_element::variant
variant_type variant() const
Definition: reference_element.h:99
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::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
fatal_macro
#define fatal_macro(message)
Definition: dis_macros.h:33
rheolef::basis_on_pointset_rep::_grad_initialize
void _grad_initialize(reference_element hat_K) const
Definition: basis_on_pointset.cc:322
rheolef::side_information_type::loc_isid
size_t loc_isid
Definition: reference_element_face_transformation.h:38
rheolef::details::memorized_matrix
Definition: memorized_value.h:70
rheolef::space_constant::valued_name
const std::string & valued_name(valued_type valued_tag)
Definition: space_constant.cc:43
rheolef
This file is part of Rheolef.
Definition: compiler_eigen.h:37
error_macro
#define error_macro(message)
Definition: dis_macros.h:49
Float
see the Float page for the full documentation
basis_on_pointset_evaluate.icc
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
mkgeo_ball.d
d
Definition: mkgeo_ball.sh:154
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::reference_element::dimension
size_type dimension() const
Definition: reference_element.h:101
mkgeo_ball.b
b
Definition: mkgeo_ball.sh:152
rheolef::side_information_type::shift
size_t shift
Definition: reference_element_face_transformation.h:39
load
void load(idiststream &in, Float &p, field &uh)
Definition: p_laplacian_post.cc:64
rheolef::space_constant::vector
Definition: space_constant.h:137
rheolef::reference_element::set_variant
void set_variant(variant_type x)
Definition: reference_element.h:93
rheolef::basis_on_pointset::set
void set(const quadrature< T > &quad, const basis_basic< T > &b)
Definition: basis_on_pointset.cc:128
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::_quad
quadrature< T > _quad
Definition: basis_on_pointset.h:127
rheolef::basis_on_pointset_rep
Definition: basis_on_pointset.h:59
rheolef::side_information_type::hat
reference_element hat
Definition: reference_element_face_transformation.h:43
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_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
T
Expr1::float_type T
Definition: field_expr.h:218
rheolef::details::memorized_side_value
Definition: memorized_value.h:92
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
trace_macro
#define trace_macro(message)
Definition: dis_macros.h:111