Rheolef  7.2
an efficient C++ finite element environment
 
Loading...
Searching...
No Matches
basis_on_pointset.cc
Go to the documentation of this file.
1
21#include "rheolef/basis_on_pointset.h"
22#include "rheolef/memorized_value.h"
23#include "rheolef/pretty_name.h"
25namespace rheolef {
26
27// ----------------------------------------------------------
28// name
29// ----------------------------------------------------------
30template<class T>
31std::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}
45template<class T>
46std::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}
57template<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// ----------------------------------------------------------
80template<class T>
82basis_on_pointset_rep<T>::make_ptr (const std::string& name)
83{
84 return new_macro(basis_on_pointset_rep<T>(name));
85}
86template <class T>
88{
89 if (name() != "") {
91 }
92}
93template<class T>
94void
95basis_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}
103template<class T>
105 : base(new_macro(rep)),
107{
108 reset (name);
109}
110template<class T>
118template<class T>
126template<class T>
127void
129{
131 reset (name);
132}
133template<class T>
134void
136{
138 reset (name);
139}
140// ----------------------------------------------------------
141// cstors
142// ----------------------------------------------------------
143template<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}
165template<class T>
166void
167basis_on_pointset_rep<T>::reset (const std::string& name)
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}
182template<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}
204template<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}
228template<class T>
231{
232 if (_mode == quad_mode) {
233 return _quad.size (hat_K);
234 } else {
235 return _nb.nnod (hat_K);
236 }
237}
238template<class T>
241{
242 return _b.ndof (hat_K);
243}
244template<class T>
245const quadrature<T>&
247{
248 check_macro (_mode == quad_mode, "get_quadrature: pointset mode is not quadrature");
249 return _quad;
250}
251template<class T>
252const 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// -----------------------------------------------------------------------
261template<class T>
262void
264 reference_element hat_K,
265 const Eigen::Matrix<point_basic<T>,Eigen::Dynamic,1>& hat_node) const
266{
267 switch (_b.valued_tag()) {
269 details::basis_on_pointset_evaluate (_b, hat_K, hat_node, _scalar_val[hat_K.variant()]);
270 break;
271 }
273 details::basis_on_pointset_evaluate (_b, hat_K, hat_node, _vector_val[hat_K.variant()]);
274 break;
275 }
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}
285template<class T>
286void
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
306template<class T>
307void
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}
320template<class T>
321void
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
391template<class T>
392void
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()) {
402 if (! _b.option().is_trace_n()) {
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. "trace_n(RTkd)"
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 }
408 break;
409 }
411 if (! _b.option().is_trace_n()) {
412 details::basis_on_pointset_evaluate (_b, tilde_L, hat_node, _sid_vector_val[L_variant][sid.loc_isid][ori_idx][sid.shift]);
413 details::basis_on_pointset_grad_evaluate (_b, tilde_L, hat_node, _sid_tensor_val[L_variant][sid.loc_isid][ori_idx][sid.shift]);
414 } else {
415 fatal_macro ("basis_trace_n(vector): not yet");
416#ifdef TODO
417 details::basis_on_pointset_evaluate_on_sides (_b, tilde_L, hat_node, _sid_vector_val[L_variant][sid.loc_isid][ori_idx][sid.shift]);
418#endif // TODO
419 }
420 break;
421 }
423 if (! _b.option().is_trace_n()) {
424 details::basis_on_pointset_evaluate (_b, tilde_L, hat_node, _sid_tensor_val[L_variant][sid.loc_isid][ori_idx][sid.shift]);
425 details::basis_on_pointset_grad_evaluate (_b, tilde_L, hat_node, _sid_tensor3_val[L_variant][sid.loc_isid][ori_idx][sid.shift]);
426 } else {
427 fatal_macro ("basis_trace_n(tensor): not yet");
428#ifdef TODO
429 details::basis_on_pointset_evaluate_on_sides (_b, tilde_L, hat_node, _sid_tensor_val[L_variant][sid.loc_isid][ori_idx][sid.shift]);
430#endif // TODO
431 }
432 break;
433 }
435 if (! _b.option().is_trace_n()) {
436 details::basis_on_pointset_evaluate (_b, tilde_L, hat_node, _sid_tensor3_val[L_variant][sid.loc_isid][ori_idx][sid.shift]);
437 details::basis_on_pointset_grad_evaluate (_b, tilde_L, hat_node, _sid_tensor4_val[L_variant][sid.loc_isid][ori_idx][sid.shift]);
438 } else {
439 fatal_macro ("basis_trace_n(tensor3): not yet");
440#ifdef TODO
441 details::basis_on_pointset_evaluate_on_sides (_b, tilde_L, hat_node, _sid_tensor_val[L_variant][sid.loc_isid][ori_idx][sid.shift]);
442#endif // TODO
443 }
444 break;
445 }
446 default: error_macro("evaluate_on_sides: unexpected "<<space_constant::valued_name(_b.valued_tag())
447 << "-valued basis \"" << _b.name() << "\"");
448 }
449}
450// TODO: avoid code repetition between initialize & grad_initialize
451template<class T>
452void
454{
455 // transform side hat_node on hat_K into tilde_node on tilde_K, the side of tilde_L:
456 // => orient and shift varies for each side
457 // the transformed nodes from hat_K to tilde_K depends upon (loc_isid,orient,shift)
458 reference_element::variant_type L_variant = tilde_L.variant();
459 Eigen::Matrix<point_basic<T>,Eigen::Dynamic,1> hat_node, tilde_node;
460 if (_mode == quad_mode) {
461 _quad.get_nodes (sid.hat, hat_node);
462 } else {
463 hat_node = _nb.hat_node (sid.hat); // ... on its Lagrange nodes.
464 }
465 tilde_node.resize (hat_node.size());
466 if (! _b.option().is_trace_n()) {
467 // transform hat_xq on side hat_S to tilde_xq on partial tilde_L
468 for (size_type q = 0, nq = hat_node.size(); q < nq; ++q) {
469 tilde_node[q] = reference_element_face_transformation (tilde_L, sid, hat_node[q]);
470 }
471 } else {
472 // else e.g. "trace_n(RTkd)": hat_nodes remains on side hat_S
473 tilde_node = hat_node;
474 }
475 _sid_initialize_continued (tilde_L, sid, tilde_node);
476}
477template<class T>
478void
480{
481 size_type d = tilde_L.dimension();
482 size_type nsid = tilde_L.n_subgeo (d-1);
484 sid.dim = d-1;
485 for (sid.loc_isid = 0; sid.loc_isid < nsid; ++sid.loc_isid) {
486 sid.n_vertex = tilde_L.subgeo_size (sid.dim, sid.loc_isid);
487 sid.hat.set_variant (sid.n_vertex, sid.dim);
488 int last_orient = (sid.dim == 0) ? 1 : -1; // no negative orient for 0d side
489 for (sid.orient = 1; sid.orient >= last_orient; sid.orient -= 2) {
490 size_type n_shift = (d == 3) ? nsid : 1;
491 for (sid.shift = 0; sid.shift < n_shift; ++sid.shift) {
492 _sid_initialize (tilde_L, sid);
493 }
494 }
495 }
496 _sid_initialized [tilde_L.variant()] = true;
497}
498template<class T>
499void
501{
502 fatal_macro ("_sid_grad_initialize: not yet");
503}
504// ----------------------------------------------------------
505// accessors
506// ----------------------------------------------------------
507template<class T>
508template<class Value>
509const Eigen::Matrix<Value,Eigen::Dynamic,Eigen::Dynamic>&
511{
512 if (!_initialized [hat_K.variant()]) _initialize (hat_K);
513 return details::memorized_matrix<T,Value>().get (*this, hat_K);
514}
515template<class T>
516template<class Value>
517const Eigen::Matrix<Value,Eigen::Dynamic,Eigen::Dynamic>&
519{
520 if (!_grad_initialized [hat_K.variant()]) _grad_initialize (hat_K);
521 return details::memorized_matrix<T,Value>().get (*this, hat_K);
522}
523template<class T>
524template<class Value>
525const Eigen::Matrix<Value,Eigen::Dynamic,Eigen::Dynamic>&
527 reference_element hat_K,
528 const side_information_type& sid) const
529{
530 if (!_sid_initialized [hat_K.variant()]) _sid_initialize (hat_K);
531 return details::memorized_side_value<T,Value>().get (*this, hat_K, sid);
532}
533template<class T>
534template<class Value>
535const Eigen::Matrix<Value,Eigen::Dynamic,Eigen::Dynamic>&
537 reference_element hat_K,
538 const side_information_type& sid) const
539{
540 if (!_sid_initialized [hat_K.variant()]) _sid_initialize (hat_K);
541 return details::memorized_side_value<T,Value>().get (*this, hat_K, sid);
542}
543// ----------------------------------------------------------------------------
544// instanciation in library
545// ----------------------------------------------------------------------------
546#define _RHEOLEF_instanciation_value(T,Value) \
547template \
548const Eigen::Matrix<Value,Eigen::Dynamic,Eigen::Dynamic>& \
549basis_on_pointset_rep<T>::evaluate (reference_element hat_K) const; \
550template \
551const Eigen::Matrix<Value,Eigen::Dynamic,Eigen::Dynamic>& \
552basis_on_pointset_rep<T>::evaluate_on_side ( \
553 reference_element hat_K, \
554 const side_information_type& sid) const; \
555template \
556const Eigen::Matrix<Value,Eigen::Dynamic,Eigen::Dynamic>& \
557basis_on_pointset_rep<T>::grad_evaluate (reference_element hat_K) const; \
558template \
559const Eigen::Matrix<Value,Eigen::Dynamic,Eigen::Dynamic>& \
560basis_on_pointset_rep<T>::grad_evaluate_on_side ( \
561 reference_element hat_K, \
562 const side_information_type& sid) const; \
563
564#define _RHEOLEF_instanciation(T) \
565template class basis_on_pointset<T>; \
566template class basis_on_pointset_rep<T>; \
567_RHEOLEF_instanciation_value(T,T) \
568_RHEOLEF_instanciation_value(T,point_basic<T>) \
569_RHEOLEF_instanciation_value(T,tensor_basic<T>) \
570_RHEOLEF_instanciation_value(T,tensor3_basic<T>) \
571_RHEOLEF_instanciation_value(T,tensor4_basic<T>) \
572
574
575} // namespace rheolef
#define _RHEOLEF_instanciation(T, M, A)
Definition asr.cc:223
see the Float page for the full documentation
std::string name() const
Definition basis.h:721
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
basis_on_pointset_rep(const std::string &name="")
basis_on_pointset_rep< T > & operator=(const basis_on_pointset_rep< 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
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 _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
basis_on_pointset(const std::string &name="")
void set(const quadrature< T > &quad, const basis_basic< T > &b)
void reset(const std::string &name)
basis_on_pointset_rep< T > rep
see the persistent_table page for the full documentation
std::string name() const
Definition quadrature.h:213
see the reference_element page for the full documentation
void set_variant(variant_type x)
size_type subgeo_size(size_type subgeo_dim, size_type loc_isid) const
variant_type variant() const
size_type n_subgeo(size_type subgeo_dim) const
#define trace_macro(message)
Definition dis_macros.h:111
#define error_macro(message)
Definition dis_macros.h:49
#define fatal_macro(message)
Definition dis_macros.h:33
Expr1::float_type T
Definition field_expr.h:230
check_macro(expr1.have_homogeneous_space(Xh1), "dual(expr1,expr2); expr1 should have homogeneous space. HINT: use dual(interpolate(Xh, expr1),expr2)")
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)
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)
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)
const std::string & valued_name(valued_type valued_tag)
This file is part of Rheolef.
point_basic< T > reference_element_face_transformation(reference_element tilde_K, const side_information_type &sid, const point_basic< T > &sid_hat_x)
void load(idiststream &in, Float &p, field &uh)
geo_element_indirect::orientation_type orient