Rheolef  7.2
an efficient C++ finite element environment
 
Loading...
Searching...
No Matches
field_expr_terminal.h
Go to the documentation of this file.
1#ifndef _RHEOLEF_FIELD_EXPR_TERMINAL_H
2#define _RHEOLEF_FIELD_EXPR_TERMINAL_H
23//
24// terminals (leaves) for the non-linear expression tree
25//
26// 0) utilities
27// 0.1) concepts
28// 0.2) terminal wrapper for exprs
29// 0.2) utility for evaluation on sides
30// 1) class-function
31// 1.1) base class for the class-function family
32// 1.2) general function or class-function
33// 1.3) normal to a surface
34// 1.4) h_local
35// 1.5) penalty
36// 2) field and such
37// 2.1) field
38// 2.2) jump of a field
39// 3) convected field, as compose(uh,X) where X is a characteristic
40//
41#include "rheolef/field.h"
42#include "rheolef/field_wdof_sliced.h"
43#include "rheolef/field_expr_utilities.h"
44#include "rheolef/expression.h"
45#include "rheolef/form.h"
46#include "rheolef/basis_on_pointset.h"
47#include "rheolef/test.h"
48#include "rheolef/characteristic.h"
49
50namespace rheolef {
51// -------------------------------------------------------------------
52// 0) utilities
53// -------------------------------------------------------------------
54namespace details {
55
56#ifdef TO_CLEAN
57// -------------------------------------------------------------------
58// 0.1) concepts
59// -------------------------------------------------------------------
60// Define a trait type for detecting field expression valid arguments
61// -> constant and field_convertible are valid terminals:
62template<class Expr, class Sfinae = void> struct is_field_expr_v2_nonlinear_arg : std::false_type {};
63template<class Expr> struct is_field_expr_v2_nonlinear_arg <Expr, typename std::enable_if<
64 is_field_expr_v2_constant<Expr>::value>::type> : std::true_type {};
65template<class Expr> struct is_field_expr_v2_nonlinear_arg <Expr, typename std::enable_if<
66 has_field_rdof_interface<Expr>::value
67 >::type> : std::true_type {};
68
69// Define a trait type for detecting linear & homogeneous expressions
70// these expressions should act homogeneously in the same finite element space
71template<class Expr> struct is_field_expr_affine_homogeneous <Expr, typename std::enable_if<
72 is_field_expr_v2_constant<Expr>::value>::type> : std::true_type {};
73template<class Expr> struct is_field_expr_affine_homogeneous <Expr, typename std::enable_if<
74 has_field_rdof_interface<Expr>::value
75 >::type> : std::true_type {};
76#endif // TO_CLEAN
77// -------------------------------------------------------------------
78// 0.2) terminal wrapper for exprs
79// -------------------------------------------------------------------
80template <class Expr, class Sfinae = void>
82{
83 // catch-all case: non-terminals are not wrapped
84 typedef Expr type;
85};
86// -------------------------------------------------------------------
87// 0.3) utility for evaluation on sides
88// -------------------------------------------------------------------
89template<class T, class M>
90const geo_element& global_get_side (const geo_basic<T,M>& omega, const geo_element& L, const side_information_type& sid);
91
92} // namespace details
93// ---------------------------------------------------------------------------
94// 1. field-function
95// ---------------------------------------------------------------------------
96// 1.1. base class for the field-function family
97// ---------------------------------------------------------------------------
98namespace details {
99
100template<class T>
102public:
103// typedefs:
104
106 typedef rheo_default_memory_model memory_type; // TODO: deduce it
107 typedef T scalar_type;
108 typedef T float_type;
109
110// allocators:
111
114
115// accessors:
116
117 void initialize (
119 const integrate_option& iopt);
120
121 void initialize (
124 const integrate_option& iopt);
125
126 template<class M>
127 const geo_element&
128 get_side (
129 const geo_basic<float_type,M>& omega_K,
130 const geo_element& K,
131 const side_information_type& sid) const;
132
133// data:
134protected:
136};
137// ---------------------------------------------------------------------------
138// 1.2) general function or class-function
139// ---------------------------------------------------------------------------
140template<class Function>
143 <typename float_traits<typename details::function_traits<Function>::result_type>::type>
144{
145public:
146// typedefs:
147
149 typedef rheo_default_memory_model memory_type; // TODO: how to deduce it ?
157
158// alocators:
159
161
162// accessors:
163
165
168
171 const integrate_option& iopt)
172 { base::initialize (pops, iopt); }
173
177 const integrate_option& iopt)
178 { base::initialize (Xh, pops, iopt); }
179
180 void evaluate (
182 const geo_element& K,
183 Eigen::Matrix<result_type,Eigen::Dynamic,1>& value) const
184 {
185 const Eigen::Matrix<piola<float_type>,Eigen::Dynamic,1>& piola = base::_pops.get_piola (omega_K, K);
186 reference_element hat_K = K.variant();
187 size_type loc_nnod = piola.size();
188 value.resize (loc_nnod);
189 for (size_type loc_inod = 0; loc_inod < loc_nnod; ++loc_inod) {
190 const point_basic<float_type>& xi = piola[loc_inod].F;
191 value[loc_inod] = _f (xi);
192 }
193 }
196 const geo_element& L,
197 const side_information_type& sid,
198 Eigen::Matrix<result_type,Eigen::Dynamic,1>& value) const
199 // TODO QUESTION_SIDE: is there some rotation and orientation to apply to the value order ?
200 { evaluate (omega_L, base::get_side(omega_L,L,sid), value); }
201
202 template<class Value>
203 bool valued_check() const {
204 static const bool status = is_equal<Value,result_type>::value;
205 check_macro (status, "unexpected result_type");
206 return status;
207 }
208// data:
209protected:
211// working area:
212public:
213 mutable std::array<
214 Eigen::Matrix<scalar_type,Eigen::Dynamic,1>
216 mutable std::array<
217 Eigen::Matrix<point_basic<scalar_type>,Eigen::Dynamic,1>
219 mutable std::array<
220 Eigen::Matrix<tensor_basic<scalar_type>,Eigen::Dynamic,1>
222 mutable std::array<
223 Eigen::Matrix<tensor3_basic<scalar_type>,Eigen::Dynamic,1>
225 mutable std::array<
226 Eigen::Matrix<tensor4_basic<scalar_type>,Eigen::Dynamic,1>
228};
229
230template<class Function>
232 : base(),
233 _f(f),
234 _scalar_val(),
235 _vector_val(),
236 _tensor_val(),
237 _tensor3_val(),
238 _tensor4_val()
239{
240}
241
242template<class Function>
243class field_expr_v2_nonlinear_terminal_function : public smart_pointer<field_expr_v2_nonlinear_terminal_function_rep<Function> >
244{
245public:
246// typedefs:
247
250 typedef typename rep::size_type size_type;
251 typedef typename rep::memory_type memory_type;
252 typedef typename rep::result_type result_type;
253 typedef typename rep::argument_type argument_type;
254 typedef typename rep::value_type value_type;
255 typedef typename rep::scalar_type scalar_type;
256 typedef typename rep::float_type float_type;
257 static const space_constant::valued_type valued_hint = rep::valued_hint;
258
259// alocators:
260
262 : base(new_macro(rep(f))) {}
263
264 template<class TrueFunction,
265 class Sfinae = typename std::enable_if<std::is_function<TrueFunction>::value, TrueFunction>::type>
267 : base(new_macro(rep(std::ptr_fun(f)))) {}
268
269 template <class Constant,
270 class Sfinae = typename std::enable_if <is_field_expr_v2_constant<Constant>::value, Constant>::type>
271 explicit field_expr_v2_nonlinear_terminal_function (const Constant& c)
272 : base(new_macro(rep(f_constant<point_basic<float_type>,result_type>(c)))) {}
273
274// accessors:
275
276 space_constant::valued_type valued_tag() const { return base::data().valued_tag(); }
277
280 const integrate_option& iopt)
281 { base::data().initialize (pops, iopt); }
282
286 const integrate_option& iopt)
287 { base::data().initialize (Xh, pops, iopt); }
288
289 void evaluate (
291 const geo_element& K,
292 Eigen::Matrix<result_type,Eigen::Dynamic,1>& value) const
293 { return base::data().evaluate (omega_K, K, value); }
294
297 const geo_element& L,
298 const side_information_type& sid,
299 Eigen::Matrix<result_type,Eigen::Dynamic,1>& value) const
300 { return base::data().evaluate_on_side (omega_L, L, sid, value); }
301
302 template<class Value>
303 bool valued_check() const {
304 return base::data().template valued_check<Value>(); }
305};
306// concept::
307template<class F> struct is_field_expr_v2_nonlinear_arg <field_expr_v2_nonlinear_terminal_function<F> > : std::true_type {};
308//template<class F> struct has_field_lazy_interface <field_expr_v2_nonlinear_terminal_function<F> > : std::true_type {};
309template<class F> struct is_field_expr_v2_nonlinear_arg <F,
310 typename std::enable_if<
311 std::conjunction<
312 std::negation<has_field_rdof_interface<F>> // filter: has op()(i_com,j_cmp) inherited from field_rdof
313 ,is_field_function<F>
314 >::value
315 >::type
316> : std::true_type {};
317
318// wrapper:
319template <class Expr>
321 typename std::enable_if<
322 std::conjunction<
323 std::negation<has_field_rdof_interface<Expr>> // filter: has op()(i_com,j_cmp) inherited from field_rdof
324 ,is_field_function<Expr>
325 >::value
326 >::type
327>
328{
330};
331// wrapper for Expr=constant (used by nary compose to wrap constant args)
332template <class Expr>
334 typename std::enable_if<
335 is_field_expr_v2_constant<Expr>::value
336 >::type
337>
338{
339 typedef typename promote<Expr,Float>::type float_type; // promote int to Float, at least
342};
343// ---------------------------------------------------------------------------
344// 1.3) normal to a surface
345// ---------------------------------------------------------------------------
346// the special class-function, used in nonlinear field expressions:
347template<class T>
348struct normal_pseudo_function : std::unary_function <point_basic<T>, point_basic<T> > {
350 fatal_macro ("special normal() class-function should not be directly evaluated");
351 return point_basic<T>();
352 }
353};
354template<class T>
357public:
358// typedefs:
359
362 typedef rheo_default_memory_model memory_type; // TODO: deduce it
366 typedef T scalar_type;
367 typedef T float_type;
369
370// alocators:
371
374
375// accessors:
376
378
380
381 void initialize (
383 const integrate_option& iopt);
384
385 void initialize (
388 const integrate_option& iopt);
389
390 void evaluate (
392 const geo_element& K,
393 Eigen::Matrix<result_type,Eigen::Dynamic,1>& value) const;
394
395 void evaluate_on_side (
397 const geo_element& L,
398 const side_information_type& sid,
399 Eigen::Matrix<result_type,Eigen::Dynamic,1>& value) const;
400
401 void evaluate_internal(
403 const geo_element& K,
404 const T& sign,
405 Eigen::Matrix<result_type,Eigen::Dynamic,1>& value) const;
406
407 template<class Value>
408 bool valued_check() const {
409 static const bool status = details::is_equal<Value,result_type>::value;
410 check_macro (status, "unexpected result_type");
411 return status;
412 }
413// data:
414 mutable bool _is_on_interface, _is_inside_on_local_sides; // normal sign fix
415};
416
417} // namespace details
418
419// the normal() pseudo-function
420template<class T>
421inline
422details::field_expr_v2_nonlinear_terminal_function <details::normal_pseudo_function<T> >
430inline
431details::field_expr_v2_nonlinear_terminal_function <details::normal_pseudo_function<Float> >
433{
434 return normal_basic<Float>();
435}
436// ---------------------------------------------------------------------------
437// 1.4) h_local
438// ---------------------------------------------------------------------------
439namespace details {
440
441// the special class-function, used in nonlinear field expressions:
442template<class T>
443struct h_local_pseudo_function : std::unary_function <point_basic<T>, T> {
445 fatal_macro ("special h_local() class-function should not be directly evaluated");
446 return 0;
447 }
448};
449
450template<class T>
453public:
454// typedefs:
455
458 typedef rheo_default_memory_model memory_type; // TODO: deduce it
460 typedef T result_type;
462 typedef T scalar_type;
463 typedef T float_type;
465
466// alocators:
467
470
471// accessors:
472
474
476
477 void initialize (
479 const integrate_option& iopt);
480
481 void initialize (
484 const integrate_option& iopt);
485
486 void evaluate (
488 const geo_element& K,
489 Eigen::Matrix<result_type,Eigen::Dynamic,1>& value) const;
490
491 void evaluate_on_side (
493 const geo_element& L,
494 const side_information_type& sid,
495 Eigen::Matrix<result_type,Eigen::Dynamic,1>& value) const;
496
497 template<class Value>
498 bool valued_check() const {
499 static const bool status = details::is_equal<Value,result_type>::value;
500 check_macro (status, "unexpected result_type");
501 return status;
502 }
503};
504
505} // namespace details
506
507// the h_local() pseudo-function
508template<class T>
509inline
510details::field_expr_v2_nonlinear_terminal_function <details::h_local_pseudo_function<T> >
518inline
519details::field_expr_v2_nonlinear_terminal_function <details::h_local_pseudo_function<Float> >
521{
522 return h_local_basic<Float>();
523}
524// ---------------------------------------------------------------------------
525// 1.5) penalty
526// ---------------------------------------------------------------------------
527// the special class-function, used in nonlinear field expressions:
528namespace details {
529
530template<class T>
531struct penalty_pseudo_function : std::unary_function <point_basic<T>, T> {
533 fatal_macro ("special penalty() class-function should not be directly evaluated");
534 return 0;
535 }
536};
537
538template<class T>
541public:
542// typedefs:
543
546 typedef rheo_default_memory_model memory_type; // TODO: deduce it
548 typedef T result_type;
550 typedef T scalar_type;
551 typedef T float_type;
553
554// alocators:
555
558
559// accessors:
560
562
564
565 void initialize (
567 const integrate_option& iopt);
568
569 void initialize (
572 const integrate_option& iopt);
573
574 void evaluate (
576 const geo_element& K,
577 Eigen::Matrix<result_type,Eigen::Dynamic,1>& value) const;
578
579 void evaluate_on_side (
581 const geo_element& L,
582 const side_information_type& sid,
583 Eigen::Matrix<result_type,Eigen::Dynamic,1>& value) const;
584
585 template<class Value>
586 bool valued_check() const {
587 static const bool status = details::is_equal<Value,result_type>::value;
588 check_macro (status, "unexpected result_type");
589 return status;
590 }
591protected:
592// internal:
593 T evaluate_measure (
595 const geo_element& K) const;
596
597 void evaluate_internal(
599 const geo_element& K,
600 const geo_element& L,
601 Eigen::Matrix<result_type,Eigen::Dynamic,1>& value) const;
602};
603
604} // namespace details
605
606// the penalty() pseudo-function
607template<class T>
608inline
609details::field_expr_v2_nonlinear_terminal_function <details::penalty_pseudo_function<T> >
617inline
618details::field_expr_v2_nonlinear_terminal_function <details::penalty_pseudo_function<Float> >
620{
621 return penalty_basic<Float>();
622}
623
624// ---------------------------------------------------------------------------
625// 2) field and such
626// ---------------------------------------------------------------------------
627// 2.1) field
628// ---------------------------------------------------------------------------
629namespace details {
630
631template<class T, class M, details::differentiate_option::type Diff>
633public:
634// typedefs:
635
637 using memory_type = M;
638 using scalar_type = T;
640 using result_type = typename std::conditional<
642 ,T // TODO: support also div(tensor) -> vector
644 >::type;
646
647// allocators:
648
649 template <class Expr,
650 class Sfinae = typename std::enable_if<details::has_field_rdof_interface<Expr>::value>::type>
651 explicit field_expr_v2_nonlinear_terminal_field_rep (const Expr& expr, const differentiate_option& gopt);
652
653#ifdef TO_CLEAN
654// --------------------------------------------
655// field_lazy interface:
656// --------------------------------------------
657
658 const geo_type& get_geo() const { return _uh.get_geo(); }
659 const band_type& get_band() const { return _uh.get_band(); }
660 const space_type& get_space() const { return _uh.get_space(); }
661 bool is_on_band () const { return _expr.is_on_band(); }
662 void initialize (const geo_type& omega_K) { _uh.initialize (omega_K); }
663#endif // TO_CLEAN
664
665// --------------------------------------------
666// accessors for the affine & homogeneous case:
667// --------------------------------------------
668
670 Vh = _uh.get_space();
671 return true;
672 }
673 // minimal forward iterator interface:
675 const_iterator begin_dof() const { return _uh.begin_dof(); }
676
677// --------------------------------------------
678// accessors for the general nonlinear case:
679// --------------------------------------------
680
682
684
685 void initialize (
686 const piola_on_pointset<T>& pops,
687 const integrate_option& iopt);
688
689 void initialize (
690 const space_basic<T,M>& Xh,
692 const integrate_option& iopt);
693
694 const geo_element& get_side (const geo_element& K, const side_information_type& sid) const;
695
696 template<class Value>
697 void evaluate (
699 const geo_element& K,
700 Eigen::Matrix<Value,Eigen::Dynamic,1>& value) const;
701
702 template<class Value>
703 void evaluate_on_side (
705 const geo_element& L,
706 const side_information_type& sid,
707 Eigen::Matrix<Value,Eigen::Dynamic,1>& value) const;
708
709 template<class Value>
710 bool valued_check() const;
711
712// data:
716 mutable bool _use_dom2bgd;
717 mutable bool _use_bgd2dom;
718 mutable bool _have_dg_on_sides;
720// working data:
721 mutable basis_on_pointset<T> _piola_on_geo_basis; // for DF, for Hdiv RT ect
722 mutable std::vector<size_type> _dis_inod_geo; // idem
724public:
725 mutable std::array<
726 Eigen::Matrix<T,Eigen::Dynamic,1>
728 mutable std::array<
729 Eigen::Matrix<point_basic<T>,Eigen::Dynamic,1>
731 mutable std::array<
732 Eigen::Matrix<tensor_basic<T>,Eigen::Dynamic,1>
734 mutable std::array<
735 Eigen::Matrix<tensor3_basic<T>,Eigen::Dynamic,1>
737 mutable std::array<
738 Eigen::Matrix<tensor4_basic<T>,Eigen::Dynamic,1>
740};
741// inlined:
742template<class T, class M, details::differentiate_option::type Diff>
743template <class Expr, class Sfinae>
744inline
746 : _gopt(gopt),
747 _uh(expr),
748 _u_test(_uh.get_space()),
749 _use_dom2bgd(false),
750 _use_bgd2dom(false),
751 _have_dg_on_sides(false),
752 _tilde_L(),
753 _piola_on_geo_basis(),
754 _dis_inod_geo(),
755 _need_vector_piola(false),
756 _scalar_val(),
757 _vector_val(),
758 _tensor_val(),
759 _tensor3_val(),
760 _tensor4_val()
761{
762 // e.g. Hdiv RTk: use DF for piola vector/tensor transformation
765 ! _uh.get_space().get_constitution().is_hierarchical();
766}
767template<class T, class M, details::differentiate_option::type Diff>
768template<class Value>
769inline
770bool
772{
773 switch (Diff) {
776 bool status = (_uh.valued_tag() == valued_tag) ||
777 (_uh.valued_tag() == space_constant::unsymmetric_tensor &&
779 check_macro (status,
780 "unexpected "<< _uh.valued()
781 << "-valued field while a " << space_constant::valued_name(valued_tag)
782 << "-valued one is expected in expression");
783 return status;
784 }
786#ifdef TODO
787 base::_fops.template grad_valued_check<Value>();
788#endif // TODO
789 return true;
790 }
792#ifdef TODO
793 base::_fops.template div_valued_check<Value>();
794#endif // TODO
795 return true;
796 }
798#ifdef TODO
799 base::_fops.template curl_valued_check<Value>();
800#endif // TODO
801 return true;
802 }
803 }
804}
805template<class T, class M, details::differentiate_option::type Diff>
806class field_expr_v2_nonlinear_terminal_field : public smart_pointer<field_expr_v2_nonlinear_terminal_field_rep<T,M,Diff> >
807{
808public:
809// typedefs:
810
813 using size_type = typename rep::size_type;
814 using memory_type = typename rep::memory_type;
815 using result_type = typename rep::result_type;
816 using float_type = typename rep::float_type;
817 using scalar_type = typename rep::scalar_type;
818 using value_type = typename rep::value_type;
819
820 static const space_constant::valued_type valued_hint = rep::valued_hint;
821
822// alocators:
823
824 template <class Expr,
825 class Sfinae = typename std::enable_if<details::has_field_rdof_interface<Expr>::value>::type>
827 : base(new_macro(rep(expr,gopt))) {}
828
829#ifdef TO_CLEAN
830// --------------------------------------------
831// field_lazy interface:
832// --------------------------------------------
833
834 const geo_type& get_geo() const { return base::data().get_geo(); }
835 const band_type& get_band() const { return base::data().get_band(); }
836 const space_type& get_space() const { return base::data().get_space(); }
837 bool is_on_band () const { return base::data().is_on_band(); }
838 void initialize (const geo_type& omega_K) { base::data().initialize (omega_K); }
839#endif // TO_CLEAN
840
841// --------------------------------------------
842// accessors for the affine & homogeneous case:
843// --------------------------------------------
844
846 return base::data().have_homogeneous_space (Vh); }
847 // minimal forward iterator interface:
848 using const_iterator = typename rep::const_iterator;
849 const_iterator begin_dof() const { return base::data().begin_dof(); }
850
851// --------------------------------------------
852// accessors for the general nonlinear case:
853// --------------------------------------------
854
855 space_constant::valued_type valued_tag() const { return base::data().valued_tag(); }
856 const geo_basic<T,M>& get_geo() const { return base::data().get_geo(); }
857 const field_basic<T,M>& expr() const { return base::data()._uh; }
858
861 const integrate_option& iopt)
862 { base::data().initialize (pops, iopt); }
863
867 const integrate_option& iopt)
868 { base::data().initialize (Xh, pops, iopt); }
869
870 template<class Value>
871 void evaluate (
873 const geo_element& K,
874 Eigen::Matrix<Value,Eigen::Dynamic,1>& value) const
875 { base::data().evaluate (omega_K, K, value); }
876
877 template<class Value>
880 const geo_element& L,
881 const side_information_type& sid,
882 Eigen::Matrix<Value,Eigen::Dynamic,1>& value) const
883 { base::data().evaluate_on_side (omega_L, L, sid, value); }
884
885 template<class Value>
886 bool valued_check() const {
887 return base::data().template valued_check<Value>();
888 }
889};
890// concepts:
891template<class T, class M, details::differentiate_option::type Diff>
893//template<class T, class M, details::differentiate_option::type Diff>
894//struct has_field_lazy_interface <field_expr_v2_nonlinear_terminal_field<T,M,Diff> > : std::true_type {};
895
896// field terminal is affine-homogeneous whentere is no differentiation
897// e.g. when diff=grad then P1 transforms to (P0)^d and then non-homogeneous space ?
898// TODO: field yh_P0_d = grad(xh_P1); => space should be computed on the fly ?
899
900template<class T, class M>
902
903// wrapper:
904template <class Expr>
906 typename std::enable_if<
907 has_field_rdof_interface<Expr>::value
908 >::type
909>
910{
911 typedef typename Expr::scalar_type T;
912 typedef typename Expr::memory_type M;
914};
915
916} // namespace details
917
919template<class Expr>
920inline
921typename
922std::enable_if<
925 typename Expr::scalar_type
926 ,typename Expr::memory_type
928 >
929>::type
930grad (const Expr& expr)
931{
932 typedef typename Expr::scalar_type T;
933 typedef typename Expr::memory_type M;
935}
937template<class Expr>
938inline
939typename
940std::enable_if<
941 details::has_field_rdof_interface<Expr>::value
942 ,details::field_expr_v2_nonlinear_terminal_field<
943 typename Expr::scalar_type
944 ,typename Expr::memory_type
946 >
947>::type
948grad_s (const Expr& expr)
949{
950 typedef typename Expr::scalar_type T;
951 typedef typename Expr::memory_type M;
952 static details::differentiate_option gopt;
953 gopt.surfacic = true;
955}
957template<class Expr>
958inline
959typename
960std::enable_if<
961 details::has_field_rdof_interface<Expr>::value
962 ,details::field_expr_v2_nonlinear_terminal_field<
963 typename Expr::scalar_type
964 ,typename Expr::memory_type
966 >
967>::type
968grad_h (const Expr& expr)
969{
970 typedef typename Expr::scalar_type T;
971 typedef typename Expr::memory_type M;
972 static details::differentiate_option gopt;
973 gopt.broken = true;
975}
977template<class Expr>
978inline
979typename
980std::enable_if<
981 details::has_field_rdof_interface<Expr>::value
982 ,details::field_expr_v2_nonlinear_terminal_field<
983 typename Expr::scalar_type
984 ,typename Expr::memory_type
986 >
987>::type
988D (const Expr& expr)
989{
990 typedef typename Expr::scalar_type T;
991 typedef typename Expr::memory_type M;
992 details::differentiate_option gopt;
993 gopt.symmetrized = true;
995}
997template<class Expr>
998inline
999typename
1000std::enable_if<
1001 details::has_field_rdof_interface<Expr>::value
1002 ,details::field_expr_v2_nonlinear_terminal_field<
1003 typename Expr::scalar_type
1004 ,typename Expr::memory_type
1006 >
1007>::type
1008Ds (const Expr& expr)
1009{
1010 typedef typename Expr::scalar_type T;
1011 typedef typename Expr::memory_type M;
1012 details::differentiate_option gopt;
1013 gopt.symmetrized = true;
1014 gopt.surfacic = true;
1016}
1018template<class Expr>
1019inline
1020typename
1021std::enable_if<
1022 details::has_field_rdof_interface<Expr>::value
1023 ,details::field_expr_v2_nonlinear_terminal_field<
1024 typename Expr::scalar_type
1025 ,typename Expr::memory_type
1027 >
1028>::type
1029Dh (const Expr& expr)
1030{
1031 typedef typename Expr::scalar_type T;
1032 typedef typename Expr::memory_type M;
1033 details::differentiate_option gopt;
1034 gopt.symmetrized = true;
1035 gopt.broken = true;
1037}
1039template<class Expr>
1040inline
1041typename
1042std::enable_if<
1043 details::has_field_rdof_interface<Expr>::value
1044 ,details::field_expr_v2_nonlinear_terminal_field<
1045 typename Expr::scalar_type
1046 ,typename Expr::memory_type
1048 >
1049>::type
1050div (const Expr& expr)
1051{
1052 typedef typename Expr::scalar_type T;
1053 typedef typename Expr::memory_type M;
1055}
1057template<class Expr>
1058inline
1059typename
1060std::enable_if<
1061 details::has_field_rdof_interface<Expr>::value
1062 ,details::field_expr_v2_nonlinear_terminal_field<
1063 typename Expr::scalar_type
1064 ,typename Expr::memory_type
1066 >
1067>::type
1068div_s (const Expr& expr)
1069{
1070 typedef typename Expr::scalar_type T;
1071 typedef typename Expr::memory_type M;
1072 static details::differentiate_option gopt;
1073 gopt.surfacic = true;
1075}
1077template<class Expr>
1078inline
1079typename
1080std::enable_if<
1081 details::has_field_rdof_interface<Expr>::value
1082 ,details::field_expr_v2_nonlinear_terminal_field<
1083 typename Expr::scalar_type
1084 ,typename Expr::memory_type
1086 >
1087>::type
1088div_h (const Expr& expr)
1089{
1090 typedef typename Expr::scalar_type T;
1091 typedef typename Expr::memory_type M;
1092 static details::differentiate_option gopt;
1093 gopt.broken = true;
1095}
1097template<class Expr>
1098inline
1099typename
1100std::enable_if<
1101 details::has_field_rdof_interface<Expr>::value
1102 ,details::field_expr_v2_nonlinear_terminal_field<
1103 typename Expr::scalar_type
1104 ,typename Expr::memory_type
1106 >
1107>::type
1108curl (const Expr& expr)
1109{
1110 typedef typename Expr::scalar_type T;
1111 typedef typename Expr::memory_type M;
1113}
1114// TODO: bcurl(uh) = Batchelor curl ?
1115// TODO: curl_s(uh) curl_h(uh)...?
1116
1117// ----------------------------------------------------------------------------
1118// 2.2) jump of a field
1119// ----------------------------------------------------------------------------
1120namespace details {
1121
1122template<class T, class M>
1124public:
1125// typedefs:
1126
1133
1134// alocators:
1135
1136 template <class Expr,
1137 class Sfinae = typename std::enable_if<details::has_field_rdof_interface<Expr>::value>::type>
1138 explicit field_expr_v2_nonlinear_terminal_field_dg_rep(const Expr& expr, const float_type& c0, const float_type& c1)
1139 : _expr0(expr), _expr1(expr),
1140 _c0(c0), _c1(c1)
1141 {}
1142
1143// accessors:
1144
1146
1148
1149 void initialize (
1151 const integrate_option& iopt);
1152
1153 void initialize (
1156 const integrate_option& iopt);
1157
1158 template<class Value>
1160 const geo_basic<float_type,memory_type>& omega_K,
1161 const geo_element& K,
1162 Eigen::Matrix<Value,Eigen::Dynamic,1>& value) const;
1163
1164 template<class Value>
1165 bool valued_check() const {
1167 bool status = _expr0.valued_tag() == valued_tag;
1168 check_macro (status,
1170 << "-valued field while a " << space_constant::valued_name(valued_tag)
1171 << "-valued one is expected in expression");
1172 return status;
1173 }
1174// data:
1175protected:
1178};
1179
1180template<class T, class M>
1181class field_expr_v2_nonlinear_terminal_field_dg : public smart_pointer<field_expr_v2_nonlinear_terminal_field_dg_rep<T,M> >
1182{
1183public:
1184// typedefs:
1185
1188 typedef typename rep::size_type size_type;
1189 typedef typename rep::memory_type memory_type;
1190 typedef typename rep::result_type result_type;
1191 typedef typename rep::float_type float_type;
1192 typedef typename rep::scalar_type scalar_type;
1193 typedef typename rep::value_type value_type;
1194 static const space_constant::valued_type valued_hint = rep::valued_hint;
1195
1196// alocators:
1197
1198 template <class Expr,
1199 class Sfinae = typename std::enable_if<details::has_field_rdof_interface<Expr>::value>::type>
1200 explicit field_expr_v2_nonlinear_terminal_field_dg(const Expr& expr, const float_type& c0, const float_type& c1)
1201 : base(new_macro(rep(expr,c0,c1))) {}
1202
1203// accessors:
1204
1205 space_constant::valued_type valued_tag() const { return base::data().valued_tag(); }
1206
1209 const integrate_option& iopt)
1210 { base::data().initialize (pops, iopt); }
1211
1215 const integrate_option& iopt)
1216 { base::data().initialize (Xh, pops, iopt); }
1217
1218 template<class Value>
1220 const geo_basic<float_type,memory_type>& omega_K,
1221 const geo_element& K,
1222 Eigen::Matrix<Value,Eigen::Dynamic,1>& value) const
1223 { base::data().evaluate (omega_K, K, value); }
1224
1225 template<class Value>
1226 bool valued_check() const
1227 { return base::data().template valued_check<Value>(); }
1228};
1229template<class T, class M> struct is_field_expr_v2_nonlinear_arg <field_expr_v2_nonlinear_terminal_field_dg<T,M> > : std::true_type {};
1230//template<class T, class M> struct has_field_lazy_interface <field_expr_v2_nonlinear_terminal_field_dg<T,M> > : std::true_type {};
1231
1232} // namespace details
1233
1234#define _RHEOLEF_make_field_expr_v2_nonlinear_terminal_field_dg(op,c0,c1) \
1235template<class Expr> \
1236inline \
1237typename \
1238std::enable_if< \
1239 details::has_field_rdof_interface<Expr>::value \
1240 ,details::field_expr_v2_nonlinear_terminal_field_dg< \
1241 typename Expr::scalar_type \
1242 ,typename Expr::memory_type \
1243 > \
1244>::type \
1245op (const Expr& expr) \
1246{ \
1247 return details::field_expr_v2_nonlinear_terminal_field_dg \
1248 <typename Expr::scalar_type ,typename Expr::memory_type> \
1249 (expr, c0, c1); \
1250}
1251
1256#undef _RHEOLEF_make_field_expr_v2_nonlinear_terminal_field_dg
1257
1258
1259// ---------------------------------------------------------------------------
1260// 3) convected field, as compose(uh,X) where X is a characteristic
1261// ---------------------------------------------------------------------------
1262namespace details {
1263
1264template<class T, class M>
1266public:
1267// typedefs:
1268
1275
1276// alocators:
1277
1280
1281// accessors:
1282
1284
1286
1287 void initialize (
1289 const integrate_option& iopt);
1290
1291 void initialize (
1294 const integrate_option& iopt);
1295
1296 template<class Value>
1297 void evaluate (
1298 const geo_basic<float_type,memory_type>& omega_K,
1299 const geo_element& K,
1300 Eigen::Matrix<Value,Eigen::Dynamic,1>& value) const;
1301
1302 template<class Value>
1303 bool valued_check() const {
1305 bool status = (_uh.valued_tag() == valued_tag);
1306 check_macro (status,
1307 "unexpected "<<_uh.valued()
1308 << "-valued field while a " << space_constant::valued_name(valued_tag)
1309 << "-valued one is expected in expression");
1310 return status;
1311 }
1312
1313// internal:
1314 template<class Result>
1315 void _check () const;
1316
1317// data:
1327};
1328template<class T, class M>
1330 const field_basic<T,M>& uh,
1332 : _uh (uh),
1333 _X (X),
1334 _fops(),
1335 _scalar_val(),
1336 _vector_val(),
1337 _tensor_val(),
1338 _tensor3_val(),
1339 _tensor4_val(),
1340 _start_q(0)
1341{
1342}
1343template<class T, class M>
1345 const field_o_characteristic<T,M>& uoX)
1346 : _uh (uoX.get_field()),
1347 _X (uoX.get_characteristic()),
1348 _fops(),
1349 _scalar_val(),
1350 _vector_val(),
1351 _tensor_val(),
1352 _tensor3_val(),
1353 _tensor4_val(),
1354 _start_q(0)
1355{
1356}
1357
1358template<class T, class M>
1360 public smart_pointer<field_expr_v2_nonlinear_terminal_field_o_characteristic_rep<T,M> >
1361{
1362public:
1363// typedefs:
1364
1367 typedef typename rep::size_type size_type;
1368 typedef typename rep::memory_type memory_type;
1369 typedef typename rep::result_type result_type;
1370 typedef typename rep::float_type float_type;
1371 typedef typename rep::scalar_type scalar_type;
1372 typedef typename rep::value_type value_type;
1373 static const space_constant::valued_type valued_hint = rep::valued_hint;
1374
1375// alocators:
1376
1379
1382
1383// accessors:
1384
1385 space_constant::valued_type valued_tag() const { return base::data().valued_tag(); }
1386
1389 const integrate_option& iopt)
1390 { base::data().initialize (pops, iopt); }
1391
1395 const integrate_option& iopt)
1396 { base::data().initialize (Xh, pops, iopt); }
1397
1398 template<class Value>
1400 const geo_basic<float_type,memory_type>& omega_K,
1401 const geo_element& K,
1402 Eigen::Matrix<Value,Eigen::Dynamic,1>& value) const
1403 { base::data().evaluate (omega_K, K, value); }
1404
1405 template<class Value>
1406 bool valued_check() const {
1407 return base::data().template valued_check<Value>();
1408 }
1409};
1410template<class T, class M> struct is_field_expr_v2_nonlinear_arg <field_expr_v2_nonlinear_terminal_field_o_characteristic<T,M> > : std::true_type {};
1411//template<class T, class M> struct has_field_lazy_interface <field_expr_v2_nonlinear_terminal_field_o_characteristic<T,M> > : std::true_type {};
1412
1413} // namespace details
1414
1415// compose a field with a characteristic
1416template<class T, class M>
1417inline
1419compose (const field_basic<T,M>& uh, const characteristic_basic<T,M>& X)
1420{
1422}
1423
1424} // namespace rheolef
1425#endif // _RHEOLEF_FIELD_EXPR_TERMINAL_H
field_expr_v2_nonlinear_terminal_field_dg_rep(const Expr &expr, const float_type &c0, const float_type &c1)
void initialize(const piola_on_pointset< float_type > &pops, const integrate_option &iopt)
void evaluate(const geo_basic< float_type, memory_type > &omega_K, const geo_element &K, Eigen::Matrix< Value, Eigen::Dynamic, 1 > &value) const
field_expr_v2_nonlinear_terminal_field< T, M, details::differentiate_option::none > _expr1
field_expr_v2_nonlinear_terminal_field< T, M, details::differentiate_option::none > _expr0
void initialize(const piola_on_pointset< float_type > &pops, const integrate_option &iopt)
field_expr_v2_nonlinear_terminal_field_dg(const Expr &expr, const float_type &c0, const float_type &c1)
void evaluate(const geo_basic< float_type, memory_type > &omega_K, const geo_element &K, Eigen::Matrix< Value, Eigen::Dynamic, 1 > &value) const
field_expr_v2_nonlinear_terminal_field_dg_rep< T, M > rep
void initialize(const space_basic< float_type, memory_type > &Xh, const piola_on_pointset< float_type > &pops, const integrate_option &iopt)
void initialize(const piola_on_pointset< float_type > &pops, const integrate_option &iopt)
void evaluate(const geo_basic< float_type, memory_type > &omega_K, const geo_element &K, Eigen::Matrix< Value, Eigen::Dynamic, 1 > &value) const
field_expr_v2_nonlinear_terminal_field_o_characteristic_rep(const field_o_characteristic< T, M > &uoX)
void initialize(const piola_on_pointset< float_type > &pops, const integrate_option &iopt)
field_expr_v2_nonlinear_terminal_field_o_characteristic(const field_o_characteristic< T, M > &uoX)
void evaluate(const geo_basic< float_type, memory_type > &omega_K, const geo_element &K, Eigen::Matrix< Value, Eigen::Dynamic, 1 > &value) const
field_expr_v2_nonlinear_terminal_field_o_characteristic_rep< T, M > rep
void initialize(const space_basic< float_type, memory_type > &Xh, const piola_on_pointset< float_type > &pops, const integrate_option &iopt)
field_expr_v2_nonlinear_terminal_field_o_characteristic(const field_basic< T, M > &uh, const characteristic_basic< T, M > &X)
std::array< Eigen::Matrix< point_basic< T >, Eigen::Dynamic, 1 >,reference_element::max_variant > _vector_val
const geo_element & get_side(const geo_element &K, const side_information_type &sid) const
typename std::conditional< Diff==details::differentiate_option::divergence,T,undeterminated_basic< T > >::type result_type
std::array< Eigen::Matrix< tensor_basic< T >, Eigen::Dynamic, 1 >,reference_element::max_variant > _tensor_val
void initialize(const piola_on_pointset< T > &pops, const integrate_option &iopt)
field_expr_v2_nonlinear_terminal_field_rep(const Expr &expr, const differentiate_option &gopt)
void evaluate(const geo_basic< float_type, memory_type > &omega_K, const geo_element &K, Eigen::Matrix< Value, Eigen::Dynamic, 1 > &value) const
std::array< Eigen::Matrix< T, Eigen::Dynamic, 1 >,reference_element::max_variant > _scalar_val
bool have_homogeneous_space(space_basic< scalar_type, memory_type > &Vh) const
std::array< Eigen::Matrix< tensor3_basic< T >, Eigen::Dynamic, 1 >,reference_element::max_variant > _tensor3_val
typename field_basic< T, M >::const_iterator const_iterator
std::array< Eigen::Matrix< tensor4_basic< T >, Eigen::Dynamic, 1 >,reference_element::max_variant > _tensor4_val
void evaluate_on_side(const geo_basic< float_type, memory_type > &omega_L, const geo_element &L, const side_information_type &sid, Eigen::Matrix< Value, Eigen::Dynamic, 1 > &value) const
void initialize(const piola_on_pointset< float_type > &pops, const integrate_option &iopt)
field_expr_v2_nonlinear_terminal_field_rep< T, M, Diff > rep
void evaluate(const geo_basic< float_type, memory_type > &omega_K, const geo_element &K, Eigen::Matrix< Value, Eigen::Dynamic, 1 > &value) const
static const space_constant::valued_type valued_hint
bool have_homogeneous_space(space_basic< scalar_type, memory_type > &Vh) const
field_expr_v2_nonlinear_terminal_field(const Expr &expr, const differentiate_option &gopt=differentiate_option())
void evaluate_on_side(const geo_basic< float_type, memory_type > &omega_L, const geo_element &L, const side_information_type &sid, Eigen::Matrix< Value, Eigen::Dynamic, 1 > &value) const
void initialize(const space_basic< float_type, memory_type > &Xh, const piola_on_pointset< float_type > &pops, const integrate_option &iopt)
void initialize(const piola_on_pointset< float_type > &pops, const integrate_option &iopt)
const geo_element & get_side(const geo_basic< float_type, M > &omega_K, const geo_element &K, const side_information_type &sid) const
field_expr_v2_nonlinear_terminal_function_base_rep< float_type > base
std::array< Eigen::Matrix< tensor_basic< scalar_type >, Eigen::Dynamic, 1 >,reference_element::max_variant > _tensor_val
details::function_traits< Function >::copiable_type function_type
details::function_traits< Function >::template arg< 0 >::type argument_type
details::function_traits< Function >::result_type result_type
std::array< Eigen::Matrix< tensor3_basic< scalar_type >, Eigen::Dynamic, 1 >,reference_element::max_variant > _tensor3_val
void initialize(const piola_on_pointset< float_type > &pops, const integrate_option &iopt)
std::array< Eigen::Matrix< tensor4_basic< scalar_type >, Eigen::Dynamic, 1 >,reference_element::max_variant > _tensor4_val
std::array< Eigen::Matrix< scalar_type, Eigen::Dynamic, 1 >,reference_element::max_variant > _scalar_val
std::array< Eigen::Matrix< point_basic< scalar_type >, Eigen::Dynamic, 1 >,reference_element::max_variant > _vector_val
void evaluate(const geo_basic< float_type, memory_type > &omega_K, const geo_element &K, Eigen::Matrix< result_type, Eigen::Dynamic, 1 > &value) const
void initialize(const space_basic< float_type, memory_type > &Xh, const piola_on_pointset< float_type > &pops, const integrate_option &iopt)
void evaluate_on_side(const geo_basic< float_type, memory_type > &omega_L, const geo_element &L, const side_information_type &sid, Eigen::Matrix< result_type, Eigen::Dynamic, 1 > &value) const
void initialize(const piola_on_pointset< float_type > &pops, const integrate_option &iopt)
void evaluate(const geo_basic< float_type, memory_type > &omega_K, const geo_element &K, Eigen::Matrix< result_type, Eigen::Dynamic, 1 > &value) const
field_expr_v2_nonlinear_terminal_function_rep< Function > rep
void initialize(const space_basic< float_type, memory_type > &Xh, const piola_on_pointset< float_type > &pops, const integrate_option &iopt)
void evaluate_on_side(const geo_basic< float_type, memory_type > &omega_L, const geo_element &L, const side_information_type &sid, Eigen::Matrix< result_type, Eigen::Dynamic, 1 > &value) const
see the disarray page for the full documentation
Definition disarray.h:497
valued_type valued_tag() const
Definition field.h:273
iterator begin_dof()
Definition field.h:595
const geo_type & get_geo() const
Definition field.h:271
const space_type & get_space() const
Definition field.h:270
const std::string & valued() const
Definition field.h:274
generic mesh with rerefence counting
Definition geo.h:1089
see the geo_element page for the full documentation
reference_element::size_type size_type
variant_type variant() const
see the integrate_option 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
the finite element space
Definition space.h:382
#define rheo_default_memory_model
#define fatal_macro(message)
Definition dis_macros.h:33
void get_geo(istream &in, my_geo &omega)
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)")
#define _RHEOLEF_make_field_expr_v2_nonlinear_terminal_field_dg(op, c0, c1)
const geo_element & global_get_side(const geo_basic< T, M > &omega_L, const geo_element &L, const side_information_type &sid)
const std::string & valued_name(valued_type valued_tag)
This file is part of Rheolef.
details::field_expr_v2_nonlinear_terminal_function< details::h_local_pseudo_function< T > > h_local_basic()
std::enable_if< details::has_field_rdof_interface< Expr >::value, details::field_expr_v2_nonlinear_terminal_field< typenameExpr::scalar_type, typenameExpr::memory_type, details::differentiate_option::gradient > >::type grad_s(const Expr &expr)
grad_s(uh): see the expression page for the full documentation
std::enable_if< details::has_field_rdof_interface< Expr >::value, details::field_expr_v2_nonlinear_terminal_field< typenameExpr::scalar_type, typenameExpr::memory_type, details::differentiate_option::gradient > >::type grad(const Expr &expr)
grad(uh): see the expression page for the full documentation
std::enable_if< details::has_field_rdof_interface< Expr >::value, details::field_expr_v2_nonlinear_terminal_field< typenameExpr::scalar_type, typenameExpr::memory_type, details::differentiate_option::gradient > >::type D(const Expr &expr)
D(uh): see the expression page for the full documentation.
details::field_expr_v2_nonlinear_terminal_function< details::normal_pseudo_function< T > > normal_basic()
std::enable_if< details::has_field_rdof_interface< Expr >::value, details::field_expr_v2_nonlinear_terminal_field< typenameExpr::scalar_type, typenameExpr::memory_type, details::differentiate_option::divergence > >::type div(const Expr &expr)
div(uh): see the expression page for the full documentation
std::enable_if< details::has_field_rdof_interface< Expr >::value, details::field_expr_v2_nonlinear_terminal_field< typenameExpr::scalar_type, typenameExpr::memory_type, details::differentiate_option::divergence > >::type div_h(const Expr &expr)
div_h(uh): see the expression page for the full documentation
std::enable_if< details::has_field_rdof_interface< Expr >::value, details::field_expr_v2_nonlinear_terminal_field< typenameExpr::scalar_type, typenameExpr::memory_type, details::differentiate_option::curl > >::type curl(const Expr &expr)
curl(uh): see the expression page for the full documentation
details::field_expr_v2_nonlinear_terminal_function< details::penalty_pseudo_function< T > > penalty_basic()
std::enable_if< details::has_field_rdof_interface< Expr >::value, details::field_expr_v2_nonlinear_terminal_field< typenameExpr::scalar_type, typenameExpr::memory_type, details::differentiate_option::gradient > >::type grad_h(const Expr &expr)
grad_h(uh): see the expression page for the full documentation
details::field_expr_v2_nonlinear_terminal_function< details::normal_pseudo_function< Float > > normal()
normal: see the expression page for the full documentation
std::enable_if< details::has_field_rdof_interface< Expr >::value, details::field_expr_v2_nonlinear_terminal_field< typenameExpr::scalar_type, typenameExpr::memory_type, details::differentiate_option::gradient > >::type Dh(const Expr &expr)
Dh(uh): see the expression page for the full documentation.
details::field_expr_v2_nonlinear_terminal_function< details::h_local_pseudo_function< Float > > h_local()
h_local: see the expression page for the full documentation
details::field_expr_v2_nonlinear_terminal_function< details::penalty_pseudo_function< Float > > penalty()
penalty(): see the expression page for the full documentation
std::enable_if< details::has_field_rdof_interface< Expr >::value, details::field_expr_v2_nonlinear_terminal_field< typenameExpr::scalar_type, typenameExpr::memory_type, details::differentiate_option::divergence > >::type div_s(const Expr &expr)
div_s(uh): see the expression page for the full documentation
std::enable_if< details::has_field_rdof_interface< Expr >::value, details::field_expr_v2_nonlinear_terminal_field< typenameExpr::scalar_type, typenameExpr::memory_type, details::differentiate_option::gradient > >::type Ds(const Expr &expr)
Ds(uh): see the expression page for the full documentation.
space_constant::valued_type valued_tag() const
details::field_expr_v2_nonlinear_node_nary< typename details::function_traits< Function >::functor_type, typename details::field_expr_v2_nonlinear_terminal_wrapper_traits< Exprs >::type... > ::type compose(const Function &f, const Exprs &... exprs)
see the compose page for the full documentation
Definition compose.h:247
STL namespace.
Definition cavity_dg.h:29
T operator()(const point_basic< T > &) const
point_basic< T > operator()(const point_basic< T > &) const
T operator()(const point_basic< T > &) const
point_basic< T > F
Definition piola.h:79
promote_not_specialized_for_this_case< T1, T2 > type
Definition promote.h:30
helper for generic field value_type: T, point_basic<T> or tensor_basic<T>
Expr1::memory_type M