20 #include "libmesh/fem_context.h" 22 #include "libmesh/boundary_info.h" 23 #include "libmesh/diff_system.h" 24 #include "libmesh/dof_map.h" 25 #include "libmesh/elem.h" 26 #include "libmesh/fe_base.h" 27 #include "libmesh/fe_interface.h" 28 #include "libmesh/libmesh_logging.h" 29 #include "libmesh/mesh_base.h" 30 #include "libmesh/numeric_vector.h" 31 #include "libmesh/quadrature.h" 32 #include "libmesh/system.h" 33 #include "libmesh/time_solver.h" 34 #include "libmesh/unsteady_solver.h" 40 const std::vector<unsigned int> * active_vars,
41 bool allocate_local_matrices)
42 :
FEMContext(sys, sys.extra_quadrature_order, active_vars,
43 allocate_local_matrices)
49 int extra_quadrature_order,
50 const std::vector<unsigned int> * active_vars,
51 bool allocate_local_matrices)
59 _custom_solution(nullptr),
60 _boundary_info(sys.get_mesh().get_boundary_info()),
62 _dim(
cast_int<unsigned char>(sys.get_mesh().mesh_dimension())),
64 _elem_dims(sys.get_mesh().elem_dimensions()),
67 _extra_quadrature_order(extra_quadrature_order)
73 std::make_unique<std::vector<unsigned int>>(*active_vars);
76 std::sort(vars_copy->begin(), vars_copy->end());
93 auto check_var = [&hardest_fe_type, &sys](
unsigned int v)
113 hardest_fe_type = fe_type;
123 return hardest_fe_type;
131 auto attach_rules = [
this, &sys](
unsigned int v)
207 _element_fe = std::vector<std::map<FEType, std::unique_ptr<FEAbstract>>>(4);
208 _side_fe = std::vector<std::map<FEType, std::unique_ptr<FEAbstract>>>(4);
215 unsigned int nv = sys.
n_vars();
218 bool have_scalar =
false;
243 auto build_var_fe = [
this, &sys](
unsigned int dim,
248 const bool add_p_level = dof_map.
should_p_refine(dof_map.var_group_from_var_number(i));
255 element_fe->add_p_level_in_reinit(add_p_level);
257 side_fe->add_p_level_in_reinit(add_p_level);
263 edge_fe->add_p_level_in_reinit(add_p_level);
283 build_var_fe(
dim, v);
286 build_var_fe(
dim, v);
312 template<
typename OutputType,
318 const unsigned int n_dofs = cast_int<unsigned int>
329 const std::vector<std::vector
335 for (
unsigned int l=0; l != n_dofs; l++)
336 u += phi[l][qp] * coef(l);
341 template<
typename OutputType,
347 const unsigned int n_dofs = cast_int<unsigned int>
358 const std::vector<std::vector
365 for (
unsigned int l=0; l != n_dofs; l++)
366 du.add_scaled(dphi[l][qp], coef(l));
373 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 374 template<
typename OutputType,
380 const unsigned int n_dofs = cast_int<unsigned int>
391 const std::vector<std::vector
398 for (
unsigned int l=0; l != n_dofs; l++)
399 d2u.add_scaled(d2phi[l][qp], coef(l));
416 template<
typename OutputType>
418 OutputType & u)
const 421 &FEMContext::get_element_fe<typename TensorTools::MakeReal<OutputType>::type>,
426 template<
typename OutputType>
429 std::vector<OutputType> & u_vals)
const 434 const unsigned int n_dofs = cast_int<unsigned int>
442 this->get_element_fe<OutputShape>( var, fe, this->
get_elem_dim() );
445 const std::vector<std::vector<OutputShape>> & phi = fe->
get_phi();
450 OutputType & u = u_vals[qp];
455 for (
unsigned int l=0; l != n_dofs; l++)
456 u += phi[l][qp] * coef(l);
463 unsigned int qp)
const 474 template<
typename OutputType>
477 OutputType & du)
const 482 <OutputType>::type>::type>,
488 template<
typename OutputType>
491 std::vector<OutputType> & du_vals)
const 498 const unsigned int n_dofs = cast_int<unsigned int>
506 this->get_element_fe<OutputShape>( var, fe, this->
get_elem_dim() );
509 const std::vector<std::vector<typename FEGenericBase<OutputShape>::OutputGradient>> & dphi = fe->
get_dphi();
514 OutputType & du = du_vals[qp];
519 for (
unsigned int l=0; l != n_dofs; l++)
520 du.add_scaled(dphi[l][qp], coef(l));
526 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 536 template<
typename OutputType>
538 OutputType & d2u)
const 545 <OutputType>::type>::type>::type>,
550 template<
typename OutputType>
553 std::vector<OutputType> & d2u_vals)
const 560 const unsigned int n_dofs = cast_int<unsigned int>
568 this->get_element_fe<OutputShape>( var, fe, this->
get_elem_dim() );
571 const std::vector<std::vector<typename FEGenericBase<OutputShape>::OutputTensor>> & d2phi = fe->
get_d2phi();
576 OutputType & d2u = d2u_vals[qp];
581 for (
unsigned int l=0; l != n_dofs; l++)
582 d2u.add_scaled(d2phi[l][qp], coef(l));
589 #endif // ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 592 template<
typename OutputType>
594 OutputType & curl_u)
const 599 const unsigned int n_dofs = cast_int<unsigned int>
608 this->get_element_fe<OutputShape>( var, fe, this->
get_elem_dim() );
611 const std::vector<std::vector<typename FEGenericBase<OutputShape>::OutputShape>> & curl_phi = fe->
get_curl_phi();
616 for (
unsigned int l=0; l != n_dofs; l++)
617 curl_u.add_scaled(curl_phi[l][qp], coef(l));
623 template<
typename OutputType>
625 OutputType & div_u)
const 632 const unsigned int n_dofs = cast_int<unsigned int>
641 this->get_element_fe<OutputShape>( var, fe, this->
get_elem_dim() );
644 const std::vector<std::vector<typename FEGenericBase<OutputShape>::OutputDivergence>> & div_phi = fe->
get_div_phi();
649 for (
unsigned int l=0; l != n_dofs; l++)
650 div_u += div_phi[l][qp] * coef(l);
657 unsigned int qp)
const 667 template<
typename OutputType>
670 OutputType & u)
const 673 &FEMContext::get_side_fe<typename TensorTools::MakeReal<OutputType>::type>,
678 template<
typename OutputType>
681 std::vector<OutputType> & u_vals)
const 686 const unsigned int n_dofs = cast_int<unsigned int>
694 this->get_side_fe<OutputShape>( var, the_side_fe, this->
get_elem_dim() );
697 const std::vector<std::vector<OutputShape>> & phi = the_side_fe->
get_phi();
702 OutputType & u = u_vals[qp];
707 for (
unsigned int l=0; l != n_dofs; l++)
708 u += phi[l][qp] * coef(l);
724 template<
typename OutputType>
726 OutputType & du)
const 733 const unsigned int n_dofs = cast_int<unsigned int>
742 this->get_side_fe<OutputShape>( var, the_side_fe, this->
get_elem_dim() );
745 const std::vector<std::vector<typename FEGenericBase<OutputShape>::OutputGradient>> & dphi = the_side_fe->
get_dphi();
750 for (
unsigned int l=0; l != n_dofs; l++)
751 du.add_scaled(dphi[l][qp], coef(l));
758 template<
typename OutputType>
761 std::vector<OutputType> & du_vals)
const 768 const unsigned int n_dofs = cast_int<unsigned int>
776 this->get_side_fe<OutputShape>( var, the_side_fe, this->
get_elem_dim() );
779 const std::vector<std::vector<typename FEGenericBase<OutputShape>::OutputGradient>> & dphi = the_side_fe->
get_dphi();
784 OutputType & du = du_vals[qp];
789 for (
unsigned int l=0; l != n_dofs; l++)
790 du.add_scaled(dphi[l][qp], coef(l));
796 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 798 unsigned int qp)
const 809 template<
typename OutputType>
812 OutputType & d2u)
const 819 <OutputType>::type>::type>::type>,
825 template<
typename OutputType>
828 std::vector<OutputType> & d2u_vals)
const 835 const unsigned int n_dofs = cast_int<unsigned int>
843 this->get_side_fe<OutputShape>( var, the_side_fe, this->
get_elem_dim() );
846 const std::vector<std::vector<typename FEGenericBase<OutputShape>::OutputTensor>> & d2phi = the_side_fe->
get_d2phi();
851 OutputType & d2u = d2u_vals[qp];
856 for (
unsigned int l=0; l != n_dofs; l++)
857 d2u.add_scaled(d2phi[l][qp], coef(l));
865 #endif // ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 878 template<
typename OutputType>
882 const Real tolerance)
const 887 const unsigned int n_dofs = cast_int<unsigned int>
896 this->get_element_fe<OutputShape>( var, fe, this->
get_elem_dim() );
903 const std::vector<std::vector<OutputShape>> & phi = fe_new->get_phi();
907 for (
unsigned int l=0; l != n_dofs; l++)
908 u += phi[l][0] * coef(l);
926 template<
typename OutputType>
930 const Real tolerance)
const 937 const unsigned int n_dofs = cast_int<unsigned int>
946 this->get_element_fe<OutputShape>( var, fe, this->
get_elem_dim() );
953 const std::vector<std::vector<typename FEGenericBase<OutputShape>::OutputGradient>> & dphi = fe_new->get_dphi();
957 for (
unsigned int l=0; l != n_dofs; l++)
958 grad_u.add_scaled(dphi[l][0], coef(l));
965 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 977 template<
typename OutputType>
981 const Real tolerance)
const 988 const unsigned int n_dofs = cast_int<unsigned int>
997 this->get_element_fe<OutputShape>( var, fe, this->
get_elem_dim() );
1004 const std::vector<std::vector<typename FEGenericBase<OutputShape>::OutputTensor>> & d2phi = fe_new->get_d2phi();
1008 for (
unsigned int l=0; l != n_dofs; l++)
1009 hess_u.add_scaled(d2phi[l][0], coef(l));
1014 #endif // LIBMESH_ENABLE_SECOND_DERIVATIVES 1017 template<
typename OutputType>
1020 OutputType & curl_u,
1021 const Real tolerance)
const 1026 const unsigned int n_dofs = cast_int<unsigned int>
1035 this->get_element_fe<OutputShape>( var, fe, this->
get_elem_dim() );
1042 const std::vector<std::vector<typename FEGenericBase<OutputShape>::OutputShape>> & curl_phi = fe_new->get_curl_phi();
1046 for (
unsigned int l=0; l != n_dofs; l++)
1047 curl_u.add_scaled(curl_phi[l][0], coef(l));
1065 template<
typename OutputType>
1067 OutputType & u)
const 1087 template<
typename OutputType>
1089 OutputType & du)
const 1096 <OutputType>::type>::type>,
1103 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 1114 template<
typename OutputType>
1116 OutputType & d2u)
const 1123 <OutputType>::type>::type>::type>,
1126 #endif // ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 1140 template<
typename OutputType>
1142 OutputType & u)
const 1164 template<
typename OutputType>
1166 OutputType & du)
const 1172 <OutputType>::type>::type>,
1178 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 1188 template<
typename OutputType>
1190 OutputType & d2u)
const 1197 <OutputType>::type>::type>::type>,
1200 #endif // ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 1213 template<
typename OutputType>
1217 const Real tolerance)
const 1222 const unsigned int n_dofs = cast_int<unsigned int>
1231 this->get_element_fe<OutputShape>( var, fe, this->
get_elem_dim() );
1238 const std::vector<std::vector<OutputShape>> & phi = fe_new->get_phi();
1242 for (
unsigned int l=0; l != n_dofs; l++)
1243 u += phi[l][0] * coef(l);
1261 template<
typename OutputType>
1264 OutputType & grad_u,
1265 const Real tolerance)
const 1272 const unsigned int n_dofs = cast_int<unsigned int>
1281 this->get_element_fe<OutputShape>( var, fe, this->
get_elem_dim() );
1288 const std::vector<std::vector<typename FEGenericBase<OutputShape>::OutputGradient>> & dphi = fe_new->get_dphi();
1292 for (
unsigned int l=0; l != n_dofs; l++)
1293 grad_u.add_scaled(dphi[l][0], coef(l));
1299 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 1312 template<
typename OutputType>
1315 OutputType & hess_u,
1316 const Real tolerance)
const 1323 const unsigned int n_dofs = cast_int<unsigned int>
1332 this->get_element_fe<OutputShape>( var, fe, this->
get_elem_dim() );
1339 const std::vector<std::vector<typename FEGenericBase<OutputShape>::OutputTensor>> & d2phi = fe_new->get_d2phi();
1343 for (
unsigned int l=0; l != n_dofs; l++)
1344 hess_u.add_scaled(d2phi[l][0], coef(l));
1349 #endif // LIBMESH_ENABLE_SECOND_DERIVATIVES 1353 template<
typename OutputType>
1355 OutputType & u)
const 1363 template<
typename OutputType>
1365 OutputType & dudot)
const 1370 <OutputType>::type>::type>,
1374 template<
typename OutputType>
1376 OutputType & u)
const 1384 template<
typename OutputType>
1386 OutputType & u)
const 1396 template<
typename OutputType>
1398 OutputType & u)
const 1485 pr.second->reinit(&(this->
get_elem()), pts);
1488 pr.second->reinit(
nullptr);
1560 for (
unsigned int i=0; i !=
n_nodes; ++i)
1564 for (
unsigned int i=0; i !=
n_nodes; ++i)
1568 for (
unsigned int i=0; i !=
n_nodes; ++i)
1585 pr.second->get_fe_map().set_jacobian_tolerance(tol);
1589 pr.second->get_fe_map().set_jacobian_tolerance(tol);
1592 pr.second->get_fe_map().set_jacobian_tolerance(tol);
1638 for (
unsigned int i=0; i !=
n_nodes; ++i)
1639 const_cast<Elem &>(this->
get_elem()).point(i)(0) =
1643 for (
unsigned int i=0; i !=
n_nodes; ++i)
1644 const_cast<Elem &>(this->
get_elem()).point(i)(1) =
1648 for (
unsigned int i=0; i !=
n_nodes; ++i)
1649 const_cast<Elem &>(this->
get_elem()).point(i)(2) =
1693 #ifdef LIBMESH_ENABLE_AMR 1705 #endif // LIBMESH_ENABLE_AMR 1707 const unsigned int n_dofs = cast_int<unsigned int>
1709 const unsigned int n_qoi = sys.
n_qois();
1751 for (std::size_t q=0; q != n_qoi; ++q)
1758 unsigned int sub_dofs = 0;
1771 #ifdef LIBMESH_ENABLE_AMR 1782 #endif // LIBMESH_ENABLE_AMR 1788 const unsigned int n_dofs_var = cast_int<unsigned int>
1797 (sub_dofs, n_dofs_var);
1808 (sub_dofs, n_dofs_var);
1815 (sub_dofs, n_dofs_var);
1821 (sub_dofs, n_dofs_var);
1826 (sub_dofs, n_dofs_var);
1828 for (std::size_t q=0; q != n_qoi; ++q)
1830 (sub_dofs, n_dofs_var);
1834 for (
unsigned int j=0; j != i; ++j)
1836 const unsigned int n_dofs_var_j =
1837 cast_int<unsigned int>
1842 n_dofs_var, n_dofs_var_j);
1845 n_dofs_var_j, n_dofs_var);
1848 (sub_dofs, sub_dofs,
1855 sub_dofs += n_dofs_var;
1863 libmesh_assert_equal_to (sub_dofs, n_dofs);
1874 for (; localized_vec_it != localized_vec_end; ++localized_vec_it)
1882 unsigned int sub_dofs = 0;
1883 auto init_localized_var_data = [
this, localized_vec_it, &sub_dofs](
unsigned int i)
1885 const unsigned int n_dofs_var = cast_int<unsigned int>
1891 localized_vec_it->second.second[i].reposition
1892 (sub_dofs, n_dofs_var);
1894 sub_dofs += n_dofs_var;
1899 init_localized_var_data(v);
1902 init_localized_var_data(v);
1904 libmesh_assert_equal_to (sub_dofs, n_dofs);
1915 cast_int<unsigned char>(this->
_elem ? this->
_elem->
dim() : 0);
1935 const int get_derivative_level )
const 1937 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 1938 const bool fe_needs_inf =
1944 fe_type !=
_real_fe->get_fe_type() ||
1950 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 1957 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 1977 const int get_derivative_level )
const 1979 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 1980 const bool fe_needs_inf =
1992 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 1999 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 2016 template<
typename OutputShape>
2020 const Real tolerance,
2021 const int get_derivative_level)
const 2030 #ifdef LIBMESH_ENABLE_AMR 2038 fe_type.
order = static_cast<Order>(fe_type.
order + add_p_level);
2040 #endif // LIBMESH_ENABLE_AMR 2045 cached_fe<OutputShape>(elem_dim, fe_type, get_derivative_level);
2046 #ifdef LIBMESH_ENABLE_AMR 2048 #endif // LIBMESH_ENABLE_AMR 2056 std::vector<Point> coor(1, master_point);
2058 switch (get_derivative_level)
2063 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 2075 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 2079 libmesh_not_implemented();
2094 fe_new->
reinit (
nullptr, &coor);
2104 template LIBMESH_EXPORT
void FEMContext::interior_value<Number>(
unsigned int,
unsigned int,
Number &)
const;
2106 std::vector<Number> &)
const;
2107 template LIBMESH_EXPORT
void FEMContext::interior_value<Gradient>(
unsigned int,
unsigned int,
Gradient &)
const;
2109 std::vector<Gradient> &)
const;
2111 template LIBMESH_EXPORT
void FEMContext::interior_gradient<Gradient>(
unsigned int,
unsigned int,
Gradient &)
const;
2113 std::vector<Gradient> &)
const;
2114 template LIBMESH_EXPORT
void FEMContext::interior_gradient<Tensor>(
unsigned int,
unsigned int,
Tensor &)
const;
2116 std::vector<Tensor> &)
const;
2118 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 2119 template LIBMESH_EXPORT
void FEMContext::interior_hessian<Tensor>(
unsigned int,
unsigned int,
Tensor &)
const;
2121 std::vector<Tensor> &)
const;
2128 template LIBMESH_EXPORT
void FEMContext::interior_curl<Gradient>(
unsigned int,
unsigned int,
Gradient &)
const;
2130 template LIBMESH_EXPORT
void FEMContext::interior_div<Number>(
unsigned int,
unsigned int,
Number &)
const;
2132 template LIBMESH_EXPORT
void FEMContext::side_value<Number>(
unsigned int,
unsigned int,
Number &)
const;
2133 template LIBMESH_EXPORT
void FEMContext::side_value<Gradient>(
unsigned int,
unsigned int,
Gradient &)
const;
2135 std::vector<Number> &)
const;
2137 std::vector<Gradient> &)
const;
2139 template LIBMESH_EXPORT
void FEMContext::side_gradient<Gradient>(
unsigned int,
unsigned int,
Gradient &)
const;
2141 std::vector<Gradient> &)
const;
2142 template LIBMESH_EXPORT
void FEMContext::side_gradient<Tensor>(
unsigned int,
unsigned int,
Tensor &)
const;
2144 std::vector<Tensor> &)
const;
2147 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 2148 template LIBMESH_EXPORT
void FEMContext::side_hessian<Tensor>(
unsigned int,
unsigned int,
Tensor &)
const;
2150 std::vector<Tensor> &)
const;
2158 template LIBMESH_EXPORT
void FEMContext::point_value<Number>(
unsigned int,
const Point &,
Number &,
const Real)
const;
2159 template LIBMESH_EXPORT
void FEMContext::point_value<Gradient>(
unsigned int,
const Point &,
Gradient &,
const Real)
const;
2161 template LIBMESH_EXPORT
void FEMContext::point_gradient<Gradient>(
unsigned int,
const Point &,
Gradient &,
const Real)
const;
2162 template LIBMESH_EXPORT
void FEMContext::point_gradient<Tensor>(
unsigned int,
const Point &,
Tensor &,
const Real)
const;
2164 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 2165 template LIBMESH_EXPORT
void FEMContext::point_hessian<Tensor>(
unsigned int,
const Point &,
Tensor &,
const Real)
const;
2170 template LIBMESH_EXPORT
void FEMContext::point_curl<Gradient>(
unsigned int,
const Point &,
Gradient &,
const Real)
const;
2172 template LIBMESH_EXPORT
void FEMContext::fixed_interior_value<Number>(
unsigned int,
unsigned int,
Number &)
const;
2173 template LIBMESH_EXPORT
void FEMContext::fixed_interior_value<Gradient>(
unsigned int,
unsigned int,
Gradient &)
const;
2175 template LIBMESH_EXPORT
void FEMContext::fixed_interior_gradient<Gradient>(
unsigned int,
unsigned int,
Gradient &)
const;
2176 template LIBMESH_EXPORT
void FEMContext::fixed_interior_gradient<Tensor>(
unsigned int,
unsigned int,
Tensor &)
const;
2178 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 2179 template LIBMESH_EXPORT
void FEMContext::fixed_interior_hessian<Tensor>(
unsigned int,
unsigned int,
Tensor &)
const;
2184 template LIBMESH_EXPORT
void FEMContext::fixed_side_value<Number>(
unsigned int,
unsigned int,
Number &)
const;
2185 template LIBMESH_EXPORT
void FEMContext::fixed_side_value<Gradient>(
unsigned int,
unsigned int,
Gradient &)
const;
2187 template LIBMESH_EXPORT
void FEMContext::fixed_side_gradient<Gradient>(
unsigned int,
unsigned int,
Gradient &)
const;
2188 template LIBMESH_EXPORT
void FEMContext::fixed_side_gradient<Tensor>(
unsigned int,
unsigned int,
Tensor &)
const;
2190 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 2191 template LIBMESH_EXPORT
void FEMContext::fixed_side_hessian<Tensor>(
unsigned int,
unsigned int,
Tensor &)
const;
2196 template LIBMESH_EXPORT
void FEMContext::fixed_point_value<Number>(
unsigned int,
const Point &,
Number &,
const Real)
const;
2197 template LIBMESH_EXPORT
void FEMContext::fixed_point_value<Gradient>(
unsigned int,
const Point &,
Gradient &,
const Real)
const;
2199 template LIBMESH_EXPORT
void FEMContext::fixed_point_gradient<Gradient>(
unsigned int,
const Point &,
Gradient &,
const Real)
const;
2200 template LIBMESH_EXPORT
void FEMContext::fixed_point_gradient<Tensor>(
unsigned int,
const Point &,
Tensor &,
const Real)
const;
2202 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 2203 template LIBMESH_EXPORT
void FEMContext::fixed_point_hessian<Tensor>(
unsigned int,
const Point &,
Tensor &,
const Real)
const;
2208 template LIBMESH_EXPORT
void FEMContext::interior_rate<Number>(
unsigned int,
unsigned int,
Number &)
const;
2209 template LIBMESH_EXPORT
void FEMContext::interior_rate<Gradient>(
unsigned int,
unsigned int,
Gradient &)
const;
2211 template LIBMESH_EXPORT
void FEMContext::interior_rate_gradient<Gradient>(
unsigned int,
unsigned int,
Gradient &)
const;
2212 template LIBMESH_EXPORT
void FEMContext::interior_rate_gradient<Tensor>(
unsigned int,
unsigned int,
Tensor &)
const;
2214 template LIBMESH_EXPORT
void FEMContext::side_rate<Number>(
unsigned int,
unsigned int,
Number &)
const;
2215 template LIBMESH_EXPORT
void FEMContext::side_rate<Gradient>(
unsigned int,
unsigned int,
Gradient &)
const;
2217 template LIBMESH_EXPORT
void FEMContext::interior_accel<Number>(
unsigned int,
unsigned int,
Number &)
const;
2218 template LIBMESH_EXPORT
void FEMContext::interior_accel<Gradient>(
unsigned int,
unsigned int,
Gradient &)
const;
2220 template LIBMESH_EXPORT
void FEMContext::side_accel<Number>(
unsigned int,
unsigned int,
Number &)
const;
2221 template LIBMESH_EXPORT
void FEMContext::side_accel<Gradient>(
unsigned int,
unsigned int,
Gradient &)
const;
class FEType hides (possibly multiple) FEFamily and approximation orders, thereby enabling specialize...
void interior_curl(unsigned int var, unsigned int qp, OutputType &curl_u) const
FEFamily family
The type of finite element.
unsigned char get_edge() const
Accessor for current edge of Elem object.
virtual void nonlocal_reinit(Real theta) override
Gives derived classes the opportunity to reinitialize data needed for nonlocal calculations at a new ...
FEGenericBase< OutputShape > * cached_fe(const unsigned int elem_dim, const FEType fe_type, const int get_derivative_level) const
virtual_for_inffe const std::vector< std::vector< OutputDivergence > > & get_div_phi() const
const DenseMatrix< Number > & get_elem_jacobian() const
Const accessor for element Jacobian.
Order
defines an enum for polynomial orders.
Number fixed_side_value(unsigned int var, unsigned int qp) const
Number side_value(unsigned int var, unsigned int qp) const
Tensor fixed_point_hessian(unsigned int var, const Point &p) const
Number interior_value(unsigned int var, unsigned int qp) const
void elem_position_set(Real theta)
Uses the coordinate data specified by mesh_*_position configuration to set the geometry of elem to th...
std::vector< DenseSubVector< Number > > _elem_fixed_subsolutions
This class provides all data required for a physics package (e.g.
void get_side_fe(unsigned int var, FEGenericBase< OutputShape > *&fe) const
Accessor for edge/face (2D/3D) finite element object for variable var for the largest dimension in th...
virtual void pre_fe_reinit(const System &, const Elem *e)
Reinitializes local data vectors/matrices on the current geometric element.
bool has_boundary_id(const Node *const node, const boundary_id_type id) const
const unsigned int invalid_uint
A number which is used quite often to represent an invalid or uninitialized value for an unsigned int...
unsigned int get_mesh_x_var() const
Accessor for x-variable of moving mesh System.
std::unique_ptr< const std::vector< unsigned int > > _active_vars
Variables on which to enable calculations, or nullptr if all variables in the System are to be enable...
const NumericVector< Number > * _custom_solution
Data with which to do algebra reinitialization.
void set_elem(const Elem *e)
Helper function to promote accessor usage.
std::unique_ptr< FEGenericBase< RealGradient > > _real_grad_fe
virtual void get(const std::vector< numeric_index_type > &index, T *values) const
Access multiple components at once.
void dof_indices(const Elem *const elem, std::vector< dof_id_type > &di) const
Fills the vector di with the global degree of freedom indices for the element.
const DenseVector< Number > & get_elem_fixed_solution() const
Accessor for element fixed solution.
const DenseVector< Number > & get_elem_solution_rate() const
Accessor for element solution rate of change w.r.t.
unsigned char get_elem_dim() const
const Elem & get_elem() const
Accessor for current Elem object.
static Point inverse_map(const unsigned int dim, const Elem *elem, const Point &p, const Real tolerance=TOLERANCE, const bool secure=true, const bool extra_checks=true)
FEMContext(const System &sys, const std::vector< unsigned int > *active_vars=nullptr, bool allocate_local_matrices=true)
Constructor.
bool _real_grad_fe_is_inf
unsigned int n_qois() const
Number of currently active quantities of interest.
void side_hessians(unsigned int var, const NumericVector< Number > &_system_vector, std::vector< OutputType > &d2u_vals) const
Fills a vector of hessians of the _system_vector at the all the quadrature points on the current elem...
void interior_accel(unsigned int var, unsigned int qp, OutputType &u) const
TensorTools::MakeReal< OutputType >::type value_shape
const DenseSubVector< Number > &(DiffContext::* diff_subsolution_getter)(unsigned int) const
Helper typedef to simplify refactoring.
void resize(const unsigned int n)
Resize the vector.
RefinementState p_refinement_flag() const
void interior_gradients(unsigned int var, const NumericVector< Number > &_system_vector, std::vector< OutputType > &interior_gradients_vector) const
Fills a vector with the gradient of the solution variable var at all the quadrature points in the cur...
std::vector< T > & get_values()
virtual void elem_edge_reinit(Real theta) override
Resets the current time in the context.
This is the base class from which all geometric element types are derived.
OrderWrapper order
The approximation order of the element.
int _extra_quadrature_order
The extra quadrature order for this context.
void boundary_ids(const Node *node, std::vector< boundary_id_type > &vec_to_fill) const
Fills a user-provided std::vector with the boundary ids associated with Node node.
void should_p_refine(unsigned int g, bool p_refine)
Describe whether the given variable group should be p-refined.
This class defines a vector in LIBMESH_DIM dimensional Real or Complex space.
std::vector< std::vector< DenseSubVector< Number > > > _elem_qoi_subderivatives
virtual bool is_steady() const =0
Is this effectively a steady-state solver?
void elem_position_get()
Uses the geometry of elem to set the coordinate data specified by mesh_*_position configuration...
The libMesh namespace provides an interface to certain functionality in the library.
void interior_values(unsigned int var, const NumericVector< Number > &_system_vector, std::vector< OutputType > &interior_values_vector) const
Fills a vector of values of the _system_vector at the all the quadrature points in the current elemen...
bool has_side_boundary_id(boundary_id_type id) const
Reports if the boundary id is found on the current side.
void use_default_quadrature_rules(int extra_quadrature_order=0)
Use quadrature rules designed to over-integrate a mass matrix, plus extra_quadrature_order.
virtual ~FEMContext()
Destructor.
void side_boundary_ids(std::vector< boundary_id_type > &vec_to_fill) const
As above, but fills in the std::set provided by the user.
unsigned char get_side() const
Accessor for current side of Elem object.
Number point_value(unsigned int var, const Point &p) const
std::vector< std::map< FEType, std::unique_ptr< FEAbstract > > > _element_fe
Finite element objects for each variable's interior, sides and edges.
virtual void elem_fe_reinit(const std::vector< Point > *const pts=nullptr)
Reinitializes interior FE objects on the current geometric element.
const bool _have_local_matrices
Whether we have local matrices allocated/initialized.
std::unique_ptr< FEGenericBase< Real > > _real_fe
Number fixed_interior_value(unsigned int var, unsigned int qp) const
Tnew cast_int(Told oldvar)
void interior_rate_gradient(unsigned int var, unsigned int qp, OutputType &u) const
Defines a dense subvector for use in finite element computations.
FEType get_fe_type() const
This class provides a specific system class.
Tensor fixed_side_hessian(unsigned int var, unsigned int qp) const
std::unique_ptr< QBase > default_quadrature_rule(const unsigned int dim, const int extraorder=0) const
const std::vector< std::vector< OutputGradient > > & get_dphi() const
std::set< unsigned char > _elem_dims
Cached dimensions of elements in the mesh, plus dimension 0 if SCALAR variables are in use...
Gradient interior_gradient(unsigned int var, unsigned int qp) const
void point_curl(unsigned int var, const Point &p, OutputType &curl_u, const Real tolerance=TOLERANCE) const
const dof_id_type n_nodes
virtual unsigned int time_order() const =0
const std::set< unsigned int > & get_second_order_vars() const
const DenseVector< Number > & get_elem_solution() const
Accessor for element solution.
unsigned char _elem_dim
Cached dimension of this->_elem.
void side_gradients(unsigned int var, const NumericVector< Number > &_system_vector, std::vector< OutputType > &side_gradients_vector) const
Fills a vector with the gradient of the solution variable var at all the quadrature points on the cur...
Tensor side_hessian(unsigned int var, unsigned int qp) const
virtual void side_fe_reinit()
Reinitializes side FE objects on the current geometric element.
std::map< FEType, std::unique_ptr< FEAbstract > > _edge_fe
const System & get_system() const
Accessor for associated system.
void _update_time_from_system(Real theta)
Update the time in the context object for the given value of theta, based on the values of "time" and...
virtual unsigned int n_nodes() const =0
void init_internal_data(const System &sys)
Helper function used in constructors to set up internal data.
Gradient fixed_point_gradient(unsigned int var, const Point &p) const
bool use_fixed_solution
A boolean to be set to true by systems using elem_fixed_solution, for optional use by e...
void side_accel(unsigned int var, unsigned int qp, OutputType &u) const
Manages consistently variables, degrees of freedom, and coefficient vectors.
void set_jacobian_tolerance(Real tol)
Calls set_jacobian_tolerance() on all the FE objects controlled by this class.
System * _mesh_sys
System from which to acquire moving mesh information.
void interior_hessians(unsigned int var, const NumericVector< Number > &_system_vector, std::vector< OutputType > &d2u_vals) const
Fills a vector of hessians of the _system_vector at the all the quadrature points in the current elem...
void some_value(unsigned int var, unsigned int qp, OutputType &u) const
Helper function to reduce some code duplication in the *interior_value methods.
void some_hessian(unsigned int var, unsigned int qp, OutputType &u) const
Helper function to reduce some code duplication in the *interior_hessian methods. ...
FEType find_hardest_fe_type()
Helper function for creating quadrature rules.
void add_p_level_in_reinit(bool value)
Indicate whether to add p-refinement levels in init/reinit methods.
Tensor interior_hessian(unsigned int var, unsigned int qp) const
const std::vector< unsigned int > * active_vars() const
Return a pointer to the vector of active variables being computed for, or a null pointer if all varia...
Tensor point_hessian(unsigned int var, const Point &p) const
Gradient fixed_interior_gradient(unsigned int var, unsigned int qp) const
static std::unique_ptr< FEGenericBase > build(const unsigned int dim, const FEType &type)
Builds a specific finite element type.
void use_unweighted_quadrature_rules(int extra_quadrature_order=0)
Use quadrature rules designed to exactly integrate unweighted undistorted basis functions, plus extra_quadrature_order.
This class provides all data required for a physics package (e.g.
virtual void reinit(const Elem *elem, const std::vector< Point > *const pts=nullptr, const std::vector< Real > *const weights=nullptr)=0
This is at the core of this class.
void side_values(unsigned int var, const NumericVector< Number > &_system_vector, std::vector< OutputType > &side_values_vector) const
Fills a vector of values of the _system_vector at the all the quadrature points on the current elemen...
const std::vector< dof_id_type > & get_dof_indices() const
Accessor for element dof indices.
const BoundaryInfo & _boundary_info
Saved reference to BoundaryInfo on the mesh for this System.
void set_time(Real time_in)
Set the time for which the current nonlinear_solution is defined.
std::vector< std::vector< FEAbstract * > > _element_fe_var
Pointers to the same finite element objects, but indexed by variable number.
const std::vector< std::vector< OutputTensor > > & get_d2phi() const
std::map< const NumericVector< Number > *, std::pair< DenseVector< Number >, std::vector< DenseSubVector< Number > > > >::iterator localized_vectors_iterator
Typedef for the localized_vectors iterator.
std::map< const NumericVector< Number > *, std::pair< DenseVector< Number >, std::vector< DenseSubVector< Number > > > > _localized_vectors
Contains pointers to vectors the user has asked to be localized, keyed with pairs of element localize...
int _real_grad_fe_derivative_level
virtual void elem_side_reinit(Real theta) override
Resets the current time in the context.
DenseSubVector< Number > & get_localized_subvector(const NumericVector< Number > &localized_vector, unsigned int var)
Return a reference to DenseSubVector localization of localized_vector at variable var contained in th...
Real get_system_time() const
Accessor for the time variable stored in the system class.
Number fixed_point_value(unsigned int var, const Point &p) const
virtual void edge_fe_reinit()
Reinitializes edge FE objects on the current geometric element.
std::vector< DenseSubVector< Number > > _elem_subsolutions
std::vector< FEAbstract * > _edge_fe_var
const DenseVector< Number > & get_elem_residual() const
Const accessor for element residual.
This is a generic class that defines a solver to handle time integration of DifferentiableSystems.
std::vector< std::unique_ptr< QBase > > _side_qrule
Quadrature rules for element sides The FEM context will try to find a quadrature rule that correctly ...
std::vector< std::map< FEType, std::unique_ptr< FEAbstract > > > _side_fe
const FEType & variable_type(const unsigned int i) const
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
std::vector< std::unique_ptr< QBase > > _element_qrule
Quadrature rule for element interior.
std::unique_ptr< QBase > unweighted_quadrature_rule(const unsigned int dim, const int extraorder=0) const
unsigned char side
Current side for side_* to examine.
virtual unsigned short dim() const =0
int _real_fe_derivative_level
const DenseVector< Number > & get_elem_solution_accel() const
Accessor for element solution accel of change w.r.t.
static std::unique_ptr< FEAbstract > build(const unsigned int dim, const FEType &type)
Builds a specific finite element type.
void resize(const unsigned int new_m, const unsigned int new_n)
Resizes the matrix to the specified size and calls zero().
IntRange< T > make_range(T beg, T end)
The 2-parameter make_range() helper function returns an IntRange<T> when both input parameters are of...
std::unique_ptr< QBase > _edge_qrule
Quadrature rules for element edges.
std::unique_ptr< NumericVector< Number > > current_local_solution
All the values I need to compute my contribution to the simulation at hand.
Helper nested class for C++03-compatible "template typedef".
void get_element_fe(unsigned int var, FEGenericBase< OutputShape > *&fe) const
Accessor for interior finite element object for variable var for the largest dimension in the mesh...
static std::unique_ptr< FEGenericBase > build_InfFE(const unsigned int dim, const FEType &type)
Builds a specific infinite element type.
Gradient fixed_side_gradient(unsigned int var, unsigned int qp) const
unsigned int get_mesh_y_var() const
Accessor for y-variable of moving mesh System.
AlgebraicType algebraic_type() const
FEGenericBase< OutputShape > * build_new_fe(const FEGenericBase< OutputShape > *fe, const Point &p, const Real tolerance=TOLERANCE, const int get_derivative_level=-1) const
Helper function to reduce some code duplication in the *_point_* methods.
void interior_div(unsigned int var, unsigned int qp, OutputType &div_u) const
virtual bool infinite() const =0
void some_gradient(unsigned int var, unsigned int qp, OutputType &u) const
Helper function to reduce some code duplication in the *interior_gradient methods.
std::vector< std::vector< FEAbstract * > > _side_fe_var
unsigned int n_vars() const
Gradient point_gradient(unsigned int var, const Point &p) const
const Elem * _elem
Current element for element_* to examine.
void attach_quadrature_rules()
Helper function for attaching quadrature rules.
void _do_elem_position_set(Real theta)
Uses the coordinate data specified by mesh_*_position configuration to set the geometry of elem to th...
virtual Order default_order() const =0
const std::vector< DenseVector< Number > > & get_qoi_derivatives() const
Const accessor for QoI derivatives.
bool has_elem() const
Test for current Elem object.
Gradient side_gradient(unsigned int var, unsigned int qp) const
const DofMap & get_dof_map() const
void interior_rate(unsigned int var, unsigned int qp, OutputType &u) const
A Point defines a location in LIBMESH_DIM dimensional Real space.
const Point & point(const unsigned int i) const
void ErrorVector unsigned int
auto index_range(const T &sizable)
Helper function that returns an IntRange<std::size_t> representing all the indices of the passed-in v...
void side_rate(unsigned int var, unsigned int qp, OutputType &u) const
This class forms the foundation from which generic finite elements may be derived.
virtual void elem_reinit(Real theta) override
Resets the current time in the context.
unsigned int get_mesh_z_var() const
Accessor for z-variable of moving mesh System.
This class defines a tensor in LIBMESH_DIM dimensional Real or Complex space.
void old_dof_indices(const Elem &elem, unsigned int n, std::vector< dof_id_type > &di, const unsigned int vn) const
Appends to the vector di the old global degree of freedom indices for elem.node_ref(n), for one variable vn.
virtual_for_inffe const std::vector< std::vector< OutputShape > > & get_curl_phi() const
const std::vector< std::vector< OutputShape > > & get_phi() const
TimeSolver & get_time_solver()
static FEFamily map_fe_type(const Elem &elem)
Tensor fixed_interior_hessian(unsigned int var, unsigned int qp) const