1 #ifndef _RHEO_FIELD_ELEMENT_H
2 #define _RHEO_FIELD_ELEMENT_H
31 #include "rheolef/quadrature.h"
32 #include "rheolef/basis.h"
33 #include "rheolef/field_expr_utilities.h"
58 template <
class Function>
68 Eigen::Matrix<Value,Eigen::Dynamic,1>
phi;
71 Value result = Value(0);
72 for (
size_t loc_idof = 0, loc_ndof =
_dof.size(); loc_idof < loc_ndof; ++loc_idof) {
73 result = result +
_dof[loc_idof]*
phi[loc_idof];
80 Eigen::Matrix<Value,Eigen::Dynamic,1> grad_phi;
83 Value result = Value();
84 for (
size_t loc_idof = 0, loc_ndof =
_dof.size(); loc_idof < loc_ndof; ++loc_idof) {
85 result = result +
_dof[loc_idof]*grad_phi[loc_idof];
91 return (*this).template evaluate<Float> (hat_x);
97 Eigen::Matrix<Float,Eigen::Dynamic,1>
_dof;
100 template<
class Value>
103 {
return _uk.template evaluate<Value> (hat_x); }
115 (*this).interpolate (
f);
120 (*this).interpolate (
f);
123 default: error_macro (
"unsupported valued field_element::interpolate(field_element)");
133 template <
class Function>
137 if (qorder ==
size_t(0)) qorder = max(
size_t(5), 2*uk.
_V.
_b.degree() + 2);
139 qopt.set_order (qorder);;
143 last_q = quad.
end(uk.
_V.
_hat_K); iter_q != last_q; iter_q++) {
144 const point& hat_xq = (*iter_q).x;
145 Float hat_wq = (*iter_q).w;
146 err2 +=
generic_norm2(
u(hat_xq) - uk.template evaluate<result_type>(hat_xq))*hat_wq;
153 template <
class Function>
157 if (nsub == 0) nsub = max(
size_t(5), 2*uk.
_V.
_b.degree() + 2);
158 Eigen::Matrix<point_basic<Float>,Eigen::Dynamic,1> hat_xnod;
161 for (
size_t loc_inod = 0, loc_nnod = hat_xnod.size(); loc_inod < loc_nnod; ++loc_inod) {
162 err2 = max (err2,
generic_norm2(uk.template evaluate<result_type>(hat_xnod[loc_inod]) -
u(hat_xnod[loc_inod])));
169 template <
class Function>
175 if (qorder ==
size_t(0)) qorder = max(
size_t(5), 2*uk.
_V.
_b.degree() + 2);
177 qopt.set_order (qorder);
181 last_q = quad.
end(uk.
_V.
_hat_K); iter_q != last_q; iter_q++) {
182 const point& hat_xq = (*iter_q).x;
183 Float hat_wq = (*iter_q).w;
184 grad_type grad_u_xq = uk.template grad_evaluate<grad_type>(hat_xq);
185 div_type div_u_xq =
tr(grad_u_xq);
193 template <
class Function>
199 if (qorder ==
size_t(0)) qorder = max(
size_t(5), 2*uk.
_V.
_b.degree() + 2);
201 qopt.set_order (qorder);
205 last_q = quad.
end(uk.
_V.
_hat_K); iter_q != last_q; iter_q++) {
206 const point& hat_xq = (*iter_q).x;
207 Float hat_wq = (*iter_q).w;
208 grad_type grad_u_xq = uk.template grad_evaluate<grad_type>(hat_xq);
209 div_type div_u_xq =
tr(grad_u_xq);
Float phi(const point &nu, Float a, Float b)
see the Float page for the full documentation
see the basis page for the full documentation
see the point page for the full documentation
rep::const_iterator const_iterator
const_iterator end(reference_element hat_K) const
const_iterator begin(reference_element hat_K) const
see the reference_element page for the full documentation
integrate_option quadrature_option
typename details::generic_binary_traits< BinaryFunction >::template result_hint< typename Expr1::result_type, typename Expr2::result_type >::type result_type
see the tensor page for the full documentation
check_macro(expr1.have_homogeneous_space(Xh1), "dual(expr1,expr2); expr1 should have homogeneous space. HINT: use dual(interpolate(Xh, expr1),expr2)")
This file is part of Rheolef.
Float error_div_linf(const field_element &uk, const Function &u, size_t qorder=size_t(0))
void pointset_lagrange_equispaced(reference_element hat_K, size_t order_in, Eigen::Matrix< point_basic< T >, Eigen::Dynamic, 1 > &hat_xnod, size_t internal=0)
Float error_l2(const field_element &uk, const Function &u, size_t qorder=size_t(0))
Float error_linf(const field_element &uk, const Function &u, size_t nsub=0)
Float measure(reference_element hat_K)
T norm2(const vec< T, M > &x)
norm2(x): see the expression page for the full documentation
U tr(const tensor_basic< U > &a, size_t d=3)
Float generic_norm2(const Float &x)
Float error_div_l2(const field_element &uk, const Function &u, size_t qorder=size_t(0))
Value operator()(const point &hat_x) const
field_element_function_wrapper(const field_element &uk)
Eigen::Matrix< Float, Eigen::Dynamic, 1 > _dof
Value grad_evaluate(const point &hat_x) const
Float operator()(const point &hat_x) const
Float dof(size_t i) const
void interpolate(const Function &f)
Value evaluate(const point &hat_x) const
field_element(const space_element &V)
space_element(std::string name, reference_element hat_K)
valued_type valued_tag() const
space_constant::valued_type valued_type
space_element(const basis &b, reference_element hat_K)
helper for generic field value_type: T, point_basic<T> or tensor_basic<T>