Go to the documentation of this file.
25 #ifndef NIHU_SINGLE_INTEGRAL_HPP_INCLUDED
26 #define NIHU_SINGLE_INTEGRAL_HPP_INCLUDED
29 #include "../util/store_pattern.hpp"
30 #include "../util/plain_type.hpp"
31 #include "../util/product_type.hpp"
32 #include "../util/block_product.hpp"
33 #include "../util/materialize_sequence.hpp"
37 #include "jacobian_computer.hpp"
46 template <
class TestField,
class TrialField>
51 typename TestField::nset_t::shape_t,
52 Eigen::Matrix<double, TestField::quantity_dimension, TrialField::quantity_dimension>,
53 typename TrialField::nset_t::shape_t
63 template <class TestField, class TrialField, class Formalism = typename get_formalism<TestField, TrialField>::type>
78 typedef typename test_field_t::elem_t
elem_t;
80 typedef typename elem_t::lset_t
lset_t;
91 static unsigned const degree
92 = test_nset_t::polynomial_order
93 + trial_nset_t::polynomial_order
94 + lset_t::jacobian_order;
116 auto acc = create_dual_field_type_accelerator(
117 test_store_t::get_data()[
degree],
118 trial_store_t::get_data()[
degree],
122 enum {Rows = TestField::quantity_dimension, Cols = TrialField::quantity_dimension};
123 auto mat = Eigen::Matrix<double, Rows, Cols>::Identity();
125 static bool printed =
false;
128 std::cout <<
"single_integral_impl::eval() called" << std::endl;
132 for (
auto it = acc.begin(); it != acc.end(); ++it)
143 test_field.eval(it.get_first()->get_xi()),
144 (mat * it.get_first()->get_w()*jac).eval(),
145 trial_field.eval(it.get_second()->get_xi())
159 template <
class TestField,
class TrialField>
168 typedef typename TrialField::elem_t trial_elem_t;
186 enum {Rows = TestField::quantity_dimension, Cols = TrialField::quantity_dimension};
187 auto mat = Eigen::Matrix<double, Rows, Cols>::Identity();
189 for (
unsigned row = 0; row < test_nset_t::num_nodes; ++row)
191 auto N = trial_field.eval(test_nset_t::corner_at(row));
192 double weight = test_field.get_corner_angle(row);
194 for (
unsigned col = 0; col < trial_nset_t::num_nodes; ++col)
195 result.template block<Rows, Cols>(row*Rows, col*Cols) = weight * mat * N(col);
208 template <
class TestField,
class TrialField,
class =
void>
222 return result_t::Zero();
232 template <
class TestField,
class TrialField>
234 typename std::enable_if<
236 typename TestField::elem_t,
237 typename TrialField::elem_t
traits_t::result_t result_t
the result matrix type of the single integral
implementation of Gaussian quadratures
single_integral_traits< TestField, TrialField > traits_t
the traits class
block_product_result_type< typename TestField::nset_t::shape_t, Eigen::Matrix< double, TestField::quantity_dimension, TrialField::quantity_dimension >, typename TrialField::nset_t::shape_t >::type result_t
result type of the single integral residual
Storage class with a static member.
single_integral_traits< TestField, TrialField >::result_t result_t
the result matrix type
trial_field_t::nset_t trial_nset_t
N-set of the trial field.
return weak form formalism from the test and trial field types
traits class of NiHu::single_integral
TestField::nset_t test_nset_t
N-set of the test field.
global constants governing some accuracy parameters
static const unsigned degree
the polynomial degree of the product of the n-sets and the jacobian
gauss_family_tag quadrature_family_t
the quadrature family
elem_t::lset_t lset_t
L-set of the elem.
static result_t eval(field_base< test_field_t > const &test_field, field_base< trial_field_t > const &trial_field)
evaluate single integral on a given field
The geometrical shape set of the element.
single integral for different element types
test_field_t::nset_t test_nset_t
N-set of the test field.
TestField test_field_t
template parameter as nested type
single_integral_traits< test_field_t, trial_field_t > traits_t
the traits class
single integral over an element for the general case
const elem_t & get_elem() const
return underlying element
CRTP base class of all fields.
TrialField::nset_t trial_nset_t
N-set of the trial field.
block_product_result_type< leftDerived, mat, rightDerived >::type block_product(Eigen::MatrixBase< leftDerived > const &l, mat const &m, Eigen::MatrixBase< rightDerived > const &r)
compute a block product l * m * r^T
const unsigned GLOBAL_MAX_ORDER
the maximal order of accelerated quadratures and field type accelerators
metafunction returning the value type of a block product
static result_t eval(field_base< TestField > const &test_field, field_base< TrialField > const &trial_field)
evaluate collocational integral on a given field
TrialField trial_field_t
template parameter as nested type
static constexpr result_t eval(field_base< TestField > const &, field_base< TrialField > const &)
specialisation of single_integral::eval for the empty case
traits_t::result_t result_t
the result matrix type
test_field_t::elem_t elem_t
Elem type.
tag for the family of Gaussian quadratures
Acceleration of field types by storing shape functions for each quadrature point.