1#ifndef _RHEOLEF_FIELD_EXPR_QUADRATURE_H
2#define _RHEOLEF_FIELD_EXPR_QUADRATURE_H
43#include "rheolef/field_expr_variational.h"
44#include "rheolef/field_expr_variational_terminal.h"
45#include "rheolef/init_expr_quadrature.h"
78 typedef geo_basic <float_type,memory_type>
geo_type;
89 template<class Sfinae = typename std::enable_if<is_field_expr_v2_variational_arg<Expr>::value, Expr>::type>
102 template<
class Value>
106 Eigen::Matrix<Value,Eigen::Dynamic,1>& lk)
const;
108 template<
class Value>
123template<
class Sfinae>
138 _pops.initialize (omega_K.get_piola_basis(), quad, new_iopt);
139 _expr.initialize (_pops, new_iopt);
149 _pops.initialize (
gh.band().get_piola_basis(), quad, new_iopt);
150 _expr.initialize (
gh, _pops, new_iopt);
158 Eigen::Matrix<Value,Eigen::Dynamic,1>& lk)
const
160 Eigen::Matrix<float_type,Eigen::Dynamic,Eigen::Dynamic> phij_xi;
161 const Eigen::Matrix<float_type,Eigen::Dynamic,1>& w = _pops.get_weight (omega_K, K);
162 _expr.evaluate (omega_K, K, phij_xi);
164 size_t ni = phij_xi.rows();
165 size_t nj = phij_xi.cols();
167 for (
size_t j = 0; j < nj; ++j) {
169 for (
size_t i = 0; i < ni; ++i) {
170 sum += w[i] * phij_xi(i,j);
195 typedef geo_basic <float_type,memory_type>
geo_type;
206 template<class Sfinae = typename std::enable_if<is_field_expr_v2_variational_arg<Expr>::value, Expr>::type>
219 template<
class Value>
223 Eigen::Matrix<Value,Eigen::Dynamic,1>& value)
const;
225 template<
class Value>
236 mutable Eigen::Matrix<float_type,Eigen::Dynamic,Eigen::Dynamic>
_value_i;
246template<
class Sfinae>
251 _ignore_sys_coord(false),
265 _pops.initialize (omega_K.get_piola_basis(), quad, new_iopt);
266 _expr.initialize (_pops, new_iopt);
276 _pops.initialize (
gh.band().get_piola_basis(), quad, new_iopt);
277 _expr.initialize (
gh, _pops, new_iopt);
278 fatal_macro(
"on_local_sides: banded level set not yet supported, sorry");
286 Eigen::Matrix<Value,Eigen::Dynamic,1>& value)
const
288 bool do_local_component_assembly =
true;
292 const geo_element& S = omega_K.dis_get_geo_element (sid_dim, dis_isid);
295 _expr.evaluate_on_side (omega_K, K, sid, _value_i, do_local_component_assembly);
296 const Eigen::Matrix<float_type,Eigen::Dynamic,1>& w = _pops.get_weight (omega_K, S);
298 size_t ni = _value_i.rows();
299 size_t nj = _value_i.cols();
301 value = Eigen::Matrix<Value,Eigen::Dynamic,1>::Zero(nj, 1);
303 for (
size_t j = 0; j < nj; ++j) {
305 for (
size_t i = 0; i < ni; ++i) {
306 sum += w[i] * _value_i(i,j);
336template<
class UnaryFunction,
class Expr>
337class field_expr_quadrature_unary {
341 typedef geo_element::size_type
size_type;
342 typedef typename Expr::memory_type memory_type;
343 typedef typename details::generic_unary_traits<UnaryFunction>::template result_hint<
344 typename Expr::value_type>::type result_hint;
345 typedef typename details::generic_unary_traits<UnaryFunction>::template hint<
346 typename Expr::value_type
347 ,result_hint>::result_type value_type;
348 typedef typename scalar_traits<value_type>::type scalar_type;
349 typedef typename float_traits<value_type>::type float_type;
350 typedef space_basic<float_type,memory_type> space_type;
351 typedef geo_basic <float_type,memory_type> geo_type;
352 typedef typename Expr::vf_tag_type vf_tag_type;
353 typedef typename details::dual_vf_tag<vf_tag_type>::type
355 typedef field_expr_quadrature_unary<UnaryFunction,Expr> self_type;
356 typedef field_expr_quadrature_unary<UnaryFunction, typename Expr::dual_self_type>
359 static const space_constant::valued_type valued_hint = space_constant::valued_tag_traits<value_type>::value;
363 template<class Sfinae = typename std::enable_if<is_field_expr_quadrature_arg<Expr>::value, Expr>::type>
364 field_expr_quadrature_unary (
const UnaryFunction&
f,
const Expr& expr)
365 : _f(
f), _expr(expr) {}
369 const space_type& get_vf_space()
const {
return _expr.get_vf_space(); }
370 size_type n_derivative()
const {
return _expr.n_derivative(); }
374 void initialize (
const geo_basic<float_type,memory_type>& omega_K,
const integrate_option& iopt) {
375 return _expr.initialize (omega_K, iopt); }
376 void initialize (
const band_basic<float_type,memory_type>&
gh,
const integrate_option& iopt) {
377 return _expr.initialize (
gh, iopt); }
379 template<
class Value>
380 void evaluate (
const geo_element& K, Eigen::Matrix<Value,Eigen::Dynamic,1>& value)
const {
382 _expr.evaluate (K, value);
383 for (
size_type i = 0, ni = value.rows(); i < ni; ++i) {
384 for (
size_type j = 0, nj = value.cols(); j < nj; ++j) {
385 value(i,j) = _f (value(i,j));
388 template<
class Value>
391 if (! is_undeterminated<A1>::value) { _expr.template valued_check<A1>(); }
398template<
class F,
class Expr>
struct is_field_expr_quadrature_arg <field_expr_quadrature_unary<F,Expr> > : std::true_type {};
405#define _RHEOLEF_make_field_expr_quadrature_unary(FUNCTION,FUNCTOR) \
406template<class Expr> \
410 details::is_field_expr_quadrature_arg<Expr>::value \
411 ,details::field_expr_quadrature_unary< \
416FUNCTION (const Expr& expr) \
418 return details::field_expr_quadrature_unary <FUNCTOR,Expr> (FUNCTOR(), expr); \
421_RHEOLEF_make_field_expr_quadrature_unary (
operator+, details::unary_plus)
422_RHEOLEF_make_field_expr_quadrature_unary (
operator-, details::negate)
423#undef _RHEOLEF_make_field_expr_quadrature_unary
436template<
class BinaryFunction,
class Expr1,
class Expr2>
445 typename Expr1::value_type
448 typename Expr1::value_type
449 ,
typename Expr2::value_type
454 typedef geo_basic <float_type,memory_type>
geo_type;
456 typename Expr1::vf_tag_type,
462 typename Expr2::dual_self_type>
468 template<
class Sfinae
469 =
typename std::enable_if<
487 _expr1.initialize (omega_K, iopt);
488 _expr2.initialize (omega_K, iopt);
494 template<
class Value>
498 Eigen::Matrix<Value,Eigen::Dynamic,1>& value)
const
505 Eigen::Matrix<A2,Eigen::Dynamic,1> value2 (value.size());
506 _expr1.evaluate (omega_K, K, value);
507 _expr2.evaluate (omega_K, K, value2);
508 for (
size_type i = 0, ni = value.size(); i < ni; ++i) {
509 value[i] =
_f (value[i], value2[i]);
512 template<
class Value>
536#define _RHEOLEF_field_expr_quadrature_binary(FUNCTION,FUNCTOR) \
537template <class Expr1, class Expr2> \
541 details::is_field_expr_quadrature_arg <Expr1>::value \
542 && details::is_field_expr_quadrature_arg <Expr2>::value \
543 ,details::field_expr_quadrature_binary< \
549FUNCTION (const Expr1& expr1, const Expr2& expr2) \
551 return details::field_expr_quadrature_binary \
552 <FUNCTOR, Expr1, Expr2> \
553 (FUNCTOR(), expr1, expr2); \
555template <class Expr1, class ExprVar2> \
559 details::is_field_expr_quadrature_arg <Expr1>::value \
560 && details::is_field_expr_v2_variational_arg <ExprVar2>::value \
561 ,details::field_expr_quadrature_binary< \
564 ,details::field_expr_quadrature_on_element<ExprVar2> \
567FUNCTION (const Expr1& expr1, const ExprVar2& expr_var2) \
569 using Expr2 = details::field_expr_quadrature_on_element<ExprVar2>; \
570 return details::field_expr_quadrature_binary \
571 <FUNCTOR, Expr1, Expr2> \
572 (FUNCTOR(), expr1, Expr2(expr_var2)); \
574template <class ExprVar1, class Expr2> \
578 details::is_field_expr_v2_variational_arg <ExprVar1>::value \
579 && details::is_field_expr_quadrature_arg <Expr2>::value \
580 ,details::field_expr_quadrature_binary< \
582 ,details::field_expr_quadrature_on_element<ExprVar1> \
586FUNCTION (const ExprVar1& expr_var1, const Expr2& expr2) \
588 using Expr1 = details::field_expr_quadrature_on_element<ExprVar1>; \
589 return details::field_expr_quadrature_binary \
590 <FUNCTOR, Expr1, Expr2> \
591 (FUNCTOR(), Expr1(expr_var1), expr2); \
597#undef _RHEOLEF_field_expr_quadrature_binary
609template<
class Expr1,
class Expr2,
class Sfinae =
void>
610struct is_field_expr_quadrature_binary_multiplies_divides_constant_left : std::false_type {};
612template<
class Expr1,
class Expr2>
613struct is_field_expr_quadrature_binary_multiplies_divides_constant_left <
618 is_rheolef_arithmetic <Expr1>::value
619 && is_field_expr_quadrature_arg<Expr2>::value
625template<
class Expr1,
class Expr2>
626struct is_field_expr_quadrature_binary_multiplies_divides_constant_right
627: is_field_expr_quadrature_binary_multiplies_divides_constant_left <Expr2,Expr1> {};
631#define _RHEOLEF_make_field_expr_quadrature_binary_operator_multiplies_divides_constant_left(FUNCTION,FUNCTOR) \
632template<class Expr1, class Expr2> \
636 details::is_field_expr_quadrature_binary_multiplies_divides_constant_left <Expr1,Expr2>::value \
637 ,details::field_expr_quadrature_unary< \
638 details::binder_first <FUNCTOR, Expr1> \
642FUNCTION (const Expr1& expr1, const Expr2& expr2) \
644 return details::field_expr_quadrature_unary \
645 <details::binder_first <FUNCTOR,Expr1>, Expr2> \
646 (details::binder_first <FUNCTOR,Expr1> (FUNCTOR(), expr1), expr2); \
649#define _RHEOLEF_make_field_expr_quadrature_binary_operator_multiplies_divides_constant_right(FUNCTION,FUNCTOR) \
650template<class Expr1, class Expr2> \
654 details::is_field_expr_quadrature_binary_multiplies_divides_constant_right <Expr1,Expr2>::value \
655 ,details::field_expr_quadrature_unary< \
656 details::binder_second <FUNCTOR, Expr2> \
660FUNCTION (const Expr1& expr1, const Expr2& expr2) \
662 return details::field_expr_quadrature_unary \
663 <details::binder_second <FUNCTOR,Expr2>, Expr1> \
664 (details::binder_second <FUNCTOR,Expr2> (FUNCTOR(), expr2), expr1); \
667#define _RHEOLEF_make_field_expr_quadrature_binary_operator_multiplies_divides_constant(FUNCTION,FUNCTOR) \
668 _RHEOLEF_make_field_expr_quadrature_binary_operator_multiplies_divides_constant_left (FUNCTION,FUNCTOR) \
669 _RHEOLEF_make_field_expr_quadrature_binary_operator_multiplies_divides_constant_right (FUNCTION,FUNCTOR)
672_RHEOLEF_make_field_expr_quadrature_binary_operator_multiplies_divides_constant (
operator*, details::multiplies)
673_RHEOLEF_make_field_expr_quadrature_binary_operator_multiplies_divides_constant_right (
operator/, details::divides)
675#undef _RHEOLEF_make_field_expr_quadrature_binary_operator_multiplies_divides_constant_right
676#undef _RHEOLEF_make_field_expr_quadrature_binary_operator_multiplies_divides_constant_left
677#undef _RHEOLEF_make_field_expr_quadrature_binary_operator_multiplies_divides_constant
field::size_type size_type
field gh(Float epsilon, Float t, const field &uh, const test &v)
void initialize(const band_basic< float_type, memory_type > &gh, const integrate_option &iopt)
geo_basic< float_type, memory_type > geo_type
field_expr_quadrature_binary< BinaryFunction, typename Expr1::dual_self_type, typename Expr2::dual_self_type > dual_self_type
space_basic< float_type, memory_type > space_type
details::generic_binary_traits< BinaryFunction >::template result_hint< typenameExpr1::value_type, typenameExpr2::value_type >::type result_hint
void initialize(const geo_basic< float_type, memory_type > &omega_K, 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
void valued_check() const
field_expr_quadrature_binary< BinaryFunction, Expr1, Expr2 > self_type
details::bf_vf_tag< BinaryFunction, typenameExpr1::vf_tag_type, typenameExpr2::vf_tag_type >::type vf_tag_type
static const space_constant::valued_type valued_hint
size_type n_derivative() const
details::dual_vf_tag< vf_tag_type >::type vf_dual_tag_type
field_expr_quadrature_binary(const BinaryFunction &f, const Expr1 &expr1, const Expr2 &expr2)
scalar_traits< value_type >::type scalar_type
details::generic_binary_traits< BinaryFunction >::template hint< typenameExpr1::value_type, typenameExpr2::value_type, result_hint >::result_type value_type
promote_memory< typenameExpr1::memory_type, typenameExpr2::memory_type >::type memory_type
float_traits< value_type >::type float_type
geo_element::size_type size_type
const space_type & get_vf_space() const
geo_basic< float_type, memory_type > geo_type
space_basic< float_type, memory_type > space_type
void initialize(const geo_basic< float_type, memory_type > &omega_K, const integrate_option &iopt)
void valued_check() const
field_expr_quadrature_on_element< typename Expr::dual_self_type > dual_self_type
static const space_constant::valued_type valued_hint
size_type n_derivative() const
details::dual_vf_tag< vf_tag_type >::type vf_dual_tag_type
scalar_traits< value_type >::type scalar_type
field_expr_quadrature_on_element< Expr > self_type
Expr::value_type value_type
Expr::memory_type memory_type
void evaluate(const geo_basic< float_type, memory_type > &omega_K, const geo_element &K, Eigen::Matrix< Value, Eigen::Dynamic, 1 > &lk) const
Expr::value_type result_hint
piola_on_pointset< float_type > _pops
field_expr_quadrature_on_element(const Expr &expr)
Expr::vf_tag_type vf_tag_type
float_traits< value_type >::type float_type
geo_element::size_type size_type
const space_type & get_vf_space() const
geo_basic< float_type, memory_type > geo_type
field_expr_quadrature_on_sides(const Expr &expr)
space_basic< float_type, memory_type > space_type
void initialize(const geo_basic< float_type, memory_type > &omega_K, 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
std::vector< size_type > _dis_inod_S
void valued_check() const
static const space_constant::valued_type valued_hint
size_type n_derivative() const
details::dual_vf_tag< vf_tag_type >::type vf_dual_tag_type
field_expr_quadrature_on_sides< Expr > self_type
tensor_basic< float_type > _DF
scalar_traits< value_type >::type scalar_type
Expr::value_type value_type
Expr::memory_type memory_type
Expr::value_type result_hint
piola_on_pointset< float_type > _pops
Expr::vf_tag_type vf_tag_type
Eigen::Matrix< float_type, Eigen::Dynamic, Eigen::Dynamic > _value_i
field_expr_quadrature_on_sides< typename Expr::dual_self_type > dual_self_type
float_traits< value_type >::type float_type
geo_element::size_type size_type
const space_type & get_vf_space() const
generic mesh with rerefence counting
see the geo_element page for the full documentation
size_type edge(size_type i) const
size_type face(size_type i) const
reference_element::size_type size_type
size_type dimension() const
orientation_type get_side_informations(const geo_element &S, size_type &loc_isid, size_type &shift) const
size_type n_subgeo(size_type subgeo_dim) const
see the integrate_option page for the full documentation
bool _is_inside_on_local_sides
#define fatal_macro(message)
#define _RHEOLEF_field_expr_quadrature_binary(FUNCTION, FUNCTOR)
integrate_option expr_quadrature_init_iopt(const geo_basic< T, M > &omega_K, const space_basic< T, M > &X, size_t n_derivative, const integrate_option &iopt)
This file is part of Rheolef.
bool valued_check() const
void initialize(const piola_on_pointset< float_type > &pops, const integrate_option &iopt)
void evaluate(const geo_basic< float_type, M > &omega_K, const geo_element &K, Eigen::Matrix< Result, Eigen::Dynamic, 1 > &value) const
std::enable_if< details::is_field_expr_v2_variational_arg< Expr >::value, details::field_expr_quadrature_on_sides< Expr > >::type on_local_sides(const Expr &expr)
on_local_sides(expr): see the expression page for the full documentation