16 #ifndef dealii_fe_values_h
17 #define dealii_fe_values_h
46 #include <type_traits>
52 #ifdef DEAL_II_WITH_PETSC
60 template <
int dim,
int spacedim = dim>
70 template <
int dim,
class NumberType =
double>
79 template <
class NumberType>
91 template <
class NumberType>
103 template <
class NumberType>
146 template <
int dim,
int spacedim = dim>
184 template <
typename Number>
193 template <
typename Number>
203 template <
typename Number>
213 template <
typename Number>
223 template <
typename Number>
233 template <
typename Number>
361 value(const
unsigned int shape_function, const
unsigned int q_point) const;
375 const
unsigned int q_point) const;
389 const
unsigned int q_point) const;
403 const
unsigned int q_point) const;
422 template <class InputVector>
425 const InputVector &fe_function,
463 template <class InputVector>
466 const InputVector &dof_values,
487 template <class InputVector>
490 const InputVector &fe_function,
500 template <class InputVector>
503 const InputVector &dof_values,
524 template <class InputVector>
527 const InputVector &fe_function,
537 template <class InputVector>
540 const InputVector &dof_values,
563 template <class InputVector>
566 const InputVector &fe_function,
576 template <class InputVector>
579 const InputVector &dof_values,
602 template <class InputVector>
605 const InputVector &fe_function,
608 &third_derivatives) const;
616 template <class InputVector>
619 const InputVector &dof_values,
622 &third_derivatives) const;
674 template <
int dim,
int spacedim = dim>
720 using curl_type = typename ::internal::CurlType<spacedim>::type;
742 template <
typename Number>
751 template <
typename Number>
761 template <
typename Number>
771 template <
typename Number>
781 template <
typename Number>
791 template <
typename Number>
800 template <
typename Number>
810 template <
typename Number>
820 template <
typename Number>
941 const unsigned int first_vector_component);
991 value(
const unsigned int shape_function,
const unsigned int q_point)
const;
1008 const unsigned int q_point)
const;
1027 const unsigned int q_point)
const;
1041 const unsigned int q_point)
const;
1064 curl(
const unsigned int shape_function,
const unsigned int q_point)
const;
1078 const unsigned int q_point)
const;
1092 const unsigned int q_point)
const;
1111 template <
class InputVector>
1114 const InputVector &fe_function,
1152 template <
class InputVector>
1155 const InputVector &dof_values,
1176 template <
class InputVector>
1179 const InputVector &fe_function,
1189 template <
class InputVector>
1192 const InputVector &dof_values,
1219 template <
class InputVector>
1221 get_function_symmetric_gradients(
1222 const InputVector &fe_function,
1225 &symmetric_gradients)
const;
1233 template <
class InputVector>
1235 get_function_symmetric_gradients_from_local_dof_values(
1236 const InputVector &dof_values,
1239 &symmetric_gradients)
const;
1259 template <
class InputVector>
1261 get_function_divergences(
1262 const InputVector &fe_function,
1264 &divergences)
const;
1272 template <
class InputVector>
1274 get_function_divergences_from_local_dof_values(
1275 const InputVector &dof_values,
1277 &divergences)
const;
1297 template <
class InputVector>
1300 const InputVector &fe_function,
1310 template <
class InputVector>
1312 get_function_curls_from_local_dof_values(
1313 const InputVector &dof_values,
1334 template <
class InputVector>
1337 const InputVector &fe_function,
1347 template <
class InputVector>
1350 const InputVector &dof_values,
1372 template <
class InputVector>
1375 const InputVector &fe_function,
1385 template <
class InputVector>
1388 const InputVector &dof_values,
1410 template <
class InputVector>
1413 const InputVector &fe_function,
1416 &third_derivatives)
const;
1424 template <
class InputVector>
1427 const InputVector &dof_values,
1430 &third_derivatives)
const;
1451 template <
int rank,
int dim,
int spacedim = dim>
1476 template <
int dim,
int spacedim>
1505 template <
typename Number>
1514 template <
typename Number>
1525 template <
typename Number>
1549 struct ShapeFunctionData
1559 bool is_nonzero_shape_function_component
1560 [value_type::n_independent_components];
1571 unsigned int row_index[value_type::n_independent_components];
1604 const unsigned int first_tensor_component);
1649 value(const
unsigned int shape_function, const
unsigned int q_point) const;
1665 divergence(const
unsigned int shape_function,
1666 const
unsigned int q_point) const;
1685 template <class InputVector>
1688 const InputVector &fe_function,
1726 template <class InputVector>
1729 const InputVector &dof_values,
1754 template <class InputVector>
1756 get_function_divergences(
1757 const InputVector &fe_function,
1759 &divergences) const;
1767 template <class InputVector>
1769 get_function_divergences_from_local_dof_values(
1770 const InputVector &dof_values,
1772 &divergences) const;
1784 const
unsigned int first_tensor_component;
1793 template <
int rank,
int dim,
int spacedim = dim>
1814 template <
int dim,
int spacedim>
1841 template <
typename Number>
1850 template <
typename Number>
1860 template <
typename Number>
1871 template <
typename Number>
1903 struct ShapeFunctionData
1913 bool is_nonzero_shape_function_component
1914 [value_type::n_independent_components];
1925 unsigned int row_index[value_type::n_independent_components];
1975 const unsigned int first_tensor_component);
2008 value(
const unsigned int shape_function,
const unsigned int q_point)
const;
2025 const unsigned int q_point)
const;
2042 const unsigned int q_point)
const;
2061 template <
class InputVector>
2064 const InputVector &fe_function,
2102 template <
class InputVector>
2105 const InputVector &dof_values,
2130 template <
class InputVector>
2132 get_function_divergences(
2133 const InputVector &fe_function,
2135 &divergences)
const;
2143 template <
class InputVector>
2145 get_function_divergences_from_local_dof_values(
2146 const InputVector &dof_values,
2148 &divergences)
const;
2166 template <
class InputVector>
2169 const InputVector &fe_function,
2179 template <
class InputVector>
2182 const InputVector &dof_values,
2215 template <
int dim,
int spacedim,
typename Extractor>
2226 template <
int dim,
int spacedim>
2229 using type = typename ::FEValuesViews::Scalar<dim, spacedim>;
2239 template <
int dim,
int spacedim>
2242 using type = typename ::FEValuesViews::Vector<dim, spacedim>;
2252 template <
int dim,
int spacedim,
int rank>
2255 using type = typename ::FEValuesViews::Tensor<rank, dim, spacedim>;
2265 template <
int dim,
int spacedim,
int rank>
2269 typename ::FEValuesViews::SymmetricTensor<rank, dim, spacedim>;
2279 template <
int dim,
int spacedim>
2286 std::vector<::FEValuesViews::Scalar<dim, spacedim>>
scalars;
2287 std::vector<::FEValuesViews::Vector<dim, spacedim>>
vectors;
2288 std::vector<::FEValuesViews::SymmetricTensor<2, dim, spacedim>>
2290 std::vector<::FEValuesViews::Tensor<2, dim, spacedim>>
2307 template <
int dim,
int spacedim,
typename Extractor>
2308 using View = typename ::internal::FEValuesViews::
2309 ViewType<dim, spacedim, Extractor>::type;
2412 template <
int dim,
int spacedim>
2419 static constexpr
unsigned int dimension = dim;
2424 static constexpr
unsigned int space_dimension = spacedim;
2462 const unsigned int dofs_per_cell,
2536 const unsigned int q_point,
2565 shape_grad(
const unsigned int i,
const unsigned int q_point)
const;
2585 const unsigned int q_point,
2628 const unsigned int q_point,
2671 const unsigned int q_point,
2714 template <
class InputVector>
2717 const InputVector & fe_function,
2718 std::vector<typename InputVector::value_type> &
values)
const;
2733 template <
class InputVector>
2736 const InputVector & fe_function,
2795 template <
class InputVector>
2798 const InputVector & fe_function,
2800 std::vector<typename InputVector::value_type> &
values)
const;
2810 template <
class InputVector>
2813 const InputVector & fe_function,
2839 template <
class InputVector>
2842 const InputVector & fe_function,
2845 const bool quadrature_points_fastest)
const;
2887 template <
class InputVector>
2890 const InputVector &fe_function,
2910 template <
class InputVector>
2913 const InputVector &fe_function,
2926 template <
class InputVector>
2929 const InputVector & fe_function,
2942 template <
class InputVector>
2945 const InputVector & fe_function,
2950 const bool quadrature_points_fastest =
false)
const;
2995 template <
class InputVector>
2998 const InputVector &fe_function,
3019 template <
class InputVector>
3022 const InputVector &fe_function,
3026 const bool quadrature_points_fastest =
false)
const;
3036 template <
class InputVector>
3039 const InputVector & fe_function,
3052 template <
class InputVector>
3055 const InputVector & fe_function,
3060 const bool quadrature_points_fastest =
false)
const;
3102 template <
class InputVector>
3105 const InputVector & fe_function,
3106 std::vector<typename InputVector::value_type> &laplacians)
const;
3127 template <
class InputVector>
3130 const InputVector & fe_function,
3141 template <
class InputVector>
3144 const InputVector & fe_function,
3146 std::vector<typename InputVector::value_type> & laplacians)
const;
3156 template <
class InputVector>
3159 const InputVector & fe_function,
3171 template <
class InputVector>
3174 const InputVector & fe_function,
3176 std::vector<std::vector<typename InputVector::value_type>> &laplacians,
3177 const bool quadrature_points_fastest =
false)
const;
3221 template <
class InputVector>
3224 const InputVector &fe_function,
3226 &third_derivatives)
const;
3246 template <
class InputVector>
3249 const InputVector &fe_function,
3252 & third_derivatives,
3253 const bool quadrature_points_fastest =
false)
const;
3263 template <
class InputVector>
3266 const InputVector & fe_function,
3269 &third_derivatives)
const;
3279 template <
class InputVector>
3282 const InputVector & fe_function,
3287 const bool quadrature_points_fastest =
false)
const;
3429 const std::vector<Point<spacedim>> &
3449 JxW(
const unsigned int q_point)
const;
3454 const std::vector<double> &
3472 const std::vector<DerivativeForm<1, dim, spacedim>> &
3491 const std::vector<DerivativeForm<2, dim, spacedim>> &
3511 const std::vector<Tensor<3, spacedim>> &
3530 const std::vector<DerivativeForm<3, dim, spacedim>> &
3551 const std::vector<Tensor<4, spacedim>> &
3571 const std::vector<DerivativeForm<4, dim, spacedim>> &
3592 const std::vector<Tensor<5, spacedim>> &
3610 const std::vector<DerivativeForm<1, spacedim, dim>> &
3642 const std::vector<Tensor<1, spacedim>> &
3643 get_normal_vectors()
const;
3731 get_cell_similarity()
const;
3751 <<
"You are requesting information from an FEValues/FEFaceValues/FESubfaceValues "
3752 <<
"object for which this kind of information has not been computed. What "
3753 <<
"information these objects compute is determined by the update_* flags you "
3754 <<
"pass to the constructor. Here, the operation you are attempting requires "
3756 <<
"> flag to be set, but it was apparently not specified "
3757 <<
"upon construction.");
3765 "FEValues object is not reinit'ed to any cell");
3775 "The FiniteElement you provided to FEValues and the FiniteElement that belongs "
3776 "to the DoFHandler that provided the cell iterator do not match.");
3784 <<
"The shape function with index " << arg1
3785 <<
" is not primitive, i.e. it is vector-valued and "
3786 <<
"has more than one non-zero vector component. This "
3787 <<
"function cannot be called for these shape functions. "
3788 <<
"Maybe you want to use the same function with the "
3789 <<
"_component suffix?");
3799 "The given FiniteElement is not a primitive element but the requested operation "
3800 "only works for those. See FiniteElement::is_primitive() for more information.");
3814 "You have previously called the FEValues::reinit() function with a "
3815 "cell iterator of type Triangulation<dim,spacedim>::cell_iterator. However, "
3816 "when you do this, you cannot call some functions in the FEValues "
3817 "class, such as the get_function_values/gradients/hessians/third_derivatives "
3818 "functions. If you need these functions, then you need to call "
3819 "FEValues::reinit() with an iterator type that allows to extract "
3820 "degrees of freedom, such as DoFHandler<dim,spacedim>::cell_iterator.");
3844 is_initialized()
const;
3859 n_dofs_for_dof_handler()
const;
3865 template <
typename VectorType>
3867 get_interpolated_dof_values(
3868 const VectorType & in,
3876 get_interpolated_dof_values(
const IndexSet & in,
3917 invalidate_present_cell();
3929 maybe_invalidate_previous_present_cell(
3943 std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
3966 std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
3992 compute_update_flags(
const UpdateFlags update_flags)
const;
4007 check_cell_similarity(
4022 template <
int,
int,
int>
4024 template <
int,
int,
int>
4039 template <
int dim,
int spacedim = dim>
4047 static constexpr
unsigned int integral_dimension = dim;
4094 template <
bool level_dof_access>
4178 template <
int dim,
int spacedim = dim>
4186 static constexpr
unsigned int integral_dimension = dim - 1;
4233 const std::vector<Tensor<1, spacedim>> &
4234 get_boundary_forms()
const;
4298 template <
int dim,
int spacedim = dim>
4306 static constexpr
unsigned int dimension = dim;
4308 static constexpr
unsigned int space_dimension = spacedim;
4314 static constexpr
unsigned int integral_dimension = dim - 1;
4359 template <
bool level_dof_access>
4363 const unsigned int face_no);
4371 template <
bool level_dof_access>
4392 const unsigned int face_no);
4444 do_reinit(
const unsigned int face_no);
4464 template <
int dim,
int spacedim = dim>
4471 static constexpr
unsigned int dimension = dim;
4476 static constexpr
unsigned int space_dimension = spacedim;
4482 static constexpr
unsigned int integral_dimension = dim - 1;
4529 template <
bool level_dof_access>
4533 const unsigned int face_no,
4534 const unsigned int subface_no);
4540 template <
bool level_dof_access>
4562 const unsigned int face_no,
4563 const unsigned int subface_no);
4634 do_reinit(
const unsigned int face_no,
const unsigned int subface_no);
4645 template <
int dim,
int spacedim>
4648 const unsigned int q_point)
const
4654 "update_values"))));
4660 return fe_values->finite_element_output.shape_values(
4668 template <
int dim,
int spacedim>
4671 const unsigned int q_point)
const
4676 "update_gradients")));
4691 template <
int dim,
int spacedim>
4694 const unsigned int q_point)
const
4699 "update_hessians")));
4713 template <
int dim,
int spacedim>
4716 const unsigned int q_point)
const
4721 "update_3rd_derivatives")));
4736 template <
int dim,
int spacedim>
4739 const unsigned int q_point)
const
4755 .single_nonzero_component_index] =
4756 fe_values->finite_element_output.shape_values(snc, q_point);
4757 return return_value;
4762 for (
unsigned int d = 0;
d < dim; ++
d)
4764 .is_nonzero_shape_function_component[
d])
4765 return_value[
d] =
fe_values->finite_element_output.shape_values(
4768 return return_value;
4774 template <
int dim,
int spacedim>
4777 const unsigned int q_point)
const
4782 "update_gradients")));
4793 .single_nonzero_component_index] =
4794 fe_values->finite_element_output.shape_gradients[snc][q_point];
4795 return return_value;
4800 for (
unsigned int d = 0;
d < dim; ++
d)
4802 .is_nonzero_shape_function_component[
d])
4804 fe_values->finite_element_output.shape_gradients
4807 return return_value;
4813 template <
int dim,
int spacedim>
4816 const unsigned int q_point)
const
4822 "update_gradients")));
4828 return divergence_type();
4832 .single_nonzero_component_index];
4835 divergence_type return_value = 0;
4836 for (
unsigned int d = 0;
d < dim; ++
d)
4838 .is_nonzero_shape_function_component[
d])
4840 fe_values->finite_element_output.shape_gradients
4843 return return_value;
4849 template <
int dim,
int spacedim>
4852 const unsigned int q_point)
const
4859 "update_gradients")));
4874 "Computing the curl in 1d is not a useful operation"));
4882 curl_type return_value;
4886 .single_nonzero_component_index == 0)
4889 .shape_gradients[snc][q_point][1];
4891 return_value[0] =
fe_values->finite_element_output
4892 .shape_gradients[snc][q_point][0];
4894 return return_value;
4899 curl_type return_value;
4901 return_value[0] = 0.0;
4904 .is_nonzero_shape_function_component[0])
4908 .row_index[0]][q_point][1];
4911 .is_nonzero_shape_function_component[1])
4915 .row_index[1]][q_point][0];
4917 return return_value;
4925 curl_type return_value;
4928 .single_nonzero_component_index)
4932 return_value[0] = 0;
4933 return_value[1] =
fe_values->finite_element_output
4934 .shape_gradients[snc][q_point][2];
4937 .shape_gradients[snc][q_point][1];
4938 return return_value;
4945 .shape_gradients[snc][q_point][2];
4946 return_value[1] = 0;
4947 return_value[2] =
fe_values->finite_element_output
4948 .shape_gradients[snc][q_point][0];
4949 return return_value;
4954 return_value[0] =
fe_values->finite_element_output
4955 .shape_gradients[snc][q_point][1];
4958 .shape_gradients[snc][q_point][0];
4959 return_value[2] = 0;
4960 return return_value;
4967 curl_type return_value;
4969 for (
unsigned int i = 0; i < dim; ++i)
4970 return_value[i] = 0.0;
4973 .is_nonzero_shape_function_component[0])
4978 .row_index[0]][q_point][2];
4982 .row_index[0]][q_point][1];
4986 .is_nonzero_shape_function_component[1])
4991 .row_index[1]][q_point][2];
4995 .row_index[1]][q_point][0];
4999 .is_nonzero_shape_function_component[2])
5004 .row_index[2]][q_point][1];
5008 .row_index[2]][q_point][0];
5011 return return_value;
5022 template <
int dim,
int spacedim>
5025 const unsigned int q_point)
const
5031 "update_hessians")));
5042 .single_nonzero_component_index] =
5043 fe_values->finite_element_output.shape_hessians[snc][q_point];
5044 return return_value;
5049 for (
unsigned int d = 0;
d < dim; ++
d)
5051 .is_nonzero_shape_function_component[
d])
5053 fe_values->finite_element_output.shape_hessians
5056 return return_value;
5062 template <
int dim,
int spacedim>
5065 const unsigned int q_point)
const
5071 "update_3rd_derivatives")));
5082 .single_nonzero_component_index] =
5083 fe_values->finite_element_output.shape_3rd_derivatives[snc][q_point];
5084 return return_value;
5089 for (
unsigned int d = 0;
d < dim; ++
d)
5091 .is_nonzero_shape_function_component[
d])
5093 fe_values->finite_element_output.shape_3rd_derivatives
5096 return return_value;
5108 inline ::SymmetricTensor<2, 1>
5109 symmetrize_single_row(
const unsigned int n, const ::Tensor<1, 1> &t)
5119 inline ::SymmetricTensor<2, 2>
5120 symmetrize_single_row(
const unsigned int n, const ::Tensor<1, 2> &t)
5126 return {{t[0], 0, t[1] / 2}};
5130 return {{0, t[1], t[0] / 2}};
5142 inline ::SymmetricTensor<2, 3>
5143 symmetrize_single_row(
const unsigned int n, const ::Tensor<1, 3> &t)
5149 return {{t[0], 0, 0, t[1] / 2, t[2] / 2, 0}};
5153 return {{0, t[1], 0, t[0] / 2, 0, t[2] / 2}};
5157 return {{0, 0, t[2], 0, t[0] / 2, t[1] / 2}};
5170 template <
int dim,
int spacedim>
5173 const unsigned int q_point)
const
5178 "update_gradients")));
5184 return symmetric_gradient_type();
5186 return internal::symmetrize_single_row(
5188 fe_values->finite_element_output.shape_gradients[snc][q_point]);
5192 for (
unsigned int d = 0;
d < dim; ++
d)
5194 .is_nonzero_shape_function_component[
d])
5196 fe_values->finite_element_output.shape_gradients
5205 template <
int dim,
int spacedim>
5208 const unsigned int q_point)
const
5228 const unsigned int comp =
5230 return_value[value_type::unrolled_to_component_indices(comp)] =
5231 fe_values->finite_element_output.shape_values(snc, q_point);
5232 return return_value;
5237 for (
unsigned int d = 0;
d < value_type::n_independent_components; ++
d)
5239 .is_nonzero_shape_function_component[
d])
5240 return_value[value_type::unrolled_to_component_indices(
d)] =
5241 fe_values->finite_element_output.shape_values(
5243 return return_value;
5249 template <
int dim,
int spacedim>
5252 const unsigned int shape_function,
5253 const unsigned int q_point)
const
5258 "update_gradients")));
5266 return divergence_type();
5289 const unsigned int comp =
5291 const unsigned int ii =
5292 value_type::unrolled_to_component_indices(comp)[0];
5293 const unsigned int jj =
5294 value_type::unrolled_to_component_indices(comp)[1];
5307 const ::Tensor<1, spacedim> &phi_grad =
5308 fe_values->finite_element_output.shape_gradients[snc][q_point];
5310 divergence_type return_value;
5311 return_value[ii] = phi_grad[jj];
5314 return_value[jj] = phi_grad[ii];
5316 return return_value;
5321 divergence_type return_value;
5322 return return_value;
5328 template <
int dim,
int spacedim>
5331 const unsigned int q_point)
const
5351 const unsigned int comp =
5355 return_value[indices] =
5356 fe_values->finite_element_output.shape_values(snc, q_point);
5357 return return_value;
5362 for (
unsigned int d = 0;
d < dim * dim; ++
d)
5364 .is_nonzero_shape_function_component[
d])
5368 return_value[indices] =
5369 fe_values->finite_element_output.shape_values(
5372 return return_value;
5378 template <
int dim,
int spacedim>
5381 const unsigned int q_point)
const
5386 "update_gradients")));
5394 return divergence_type();
5408 const unsigned int comp =
5412 const unsigned int ii = indices[0];
5413 const unsigned int jj = indices[1];
5415 const ::Tensor<1, spacedim> &phi_grad =
5416 fe_values->finite_element_output.shape_gradients[snc][q_point];
5418 divergence_type return_value;
5420 return_value[ii] = phi_grad[jj];
5422 return return_value;
5427 divergence_type return_value;
5428 return return_value;
5434 template <
int dim,
int spacedim>
5437 const unsigned int q_point)
const
5442 "update_gradients")));
5464 const unsigned int comp =
5468 const unsigned int ii = indices[0];
5469 const unsigned int jj = indices[1];
5471 const ::Tensor<1, spacedim> &phi_grad =
5472 fe_values->finite_element_output.shape_gradients[snc][q_point];
5475 return_value[ii][jj] = phi_grad;
5477 return return_value;
5483 return return_value;
5495 template <
int dim,
int spacedim>
5502 , dof_handler(&cell->get_dof_handler())
5503 , level_dof_access(lda)
5508 template <
int dim,
int spacedim>
5515 return fe_values_views_cache.scalars[
scalar.component];
5520 template <
int dim,
int spacedim>
5526 fe_values_views_cache.vectors.size());
5533 template <
int dim,
int spacedim>
5540 fe_values_views_cache.symmetric_second_order_tensors.size(),
5543 fe_values_views_cache.symmetric_second_order_tensors.size()));
5545 return fe_values_views_cache
5551 template <
int dim,
int spacedim>
5557 fe_values_views_cache.second_order_tensors.size());
5559 return fe_values_views_cache
5565 template <
int dim,
int spacedim>
5566 inline const double &
5568 const unsigned int q_point)
const
5573 Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
5574 Assert(present_cell.is_initialized(), ExcNotReinited());
5577 if (fe->is_primitive())
5578 return this->finite_element_output.shape_values(i, q_point);
5589 const unsigned int row =
5590 this->finite_element_output
5591 .shape_function_to_row_table[i * fe->n_components() +
5592 fe->system_to_component_index(i).first];
5593 return this->finite_element_output.shape_values(row, q_point);
5599 template <
int dim,
int spacedim>
5602 const unsigned int i,
5603 const unsigned int q_point,
5604 const unsigned int component)
const
5610 Assert(present_cell.is_initialized(), ExcNotReinited());
5615 if (fe->get_nonzero_components(i)[component] ==
false)
5621 const unsigned int row =
5622 this->finite_element_output
5623 .shape_function_to_row_table[i * fe->n_components() + component];
5624 return this->finite_element_output.shape_values(row, q_point);
5629 template <
int dim,
int spacedim>
5632 const unsigned int q_point)
const
5637 Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
5638 Assert(present_cell.is_initialized(), ExcNotReinited());
5641 if (fe->is_primitive())
5642 return this->finite_element_output.shape_gradients[i][q_point];
5653 const unsigned int row =
5654 this->finite_element_output
5655 .shape_function_to_row_table[i * fe->n_components() +
5656 fe->system_to_component_index(i).first];
5657 return this->finite_element_output.shape_gradients[row][q_point];
5663 template <
int dim,
int spacedim>
5666 const unsigned int i,
5667 const unsigned int q_point,
5668 const unsigned int component)
const
5674 Assert(present_cell.is_initialized(), ExcNotReinited());
5678 if (fe->get_nonzero_components(i)[component] ==
false)
5684 const unsigned int row =
5685 this->finite_element_output
5686 .shape_function_to_row_table[i * fe->n_components() + component];
5687 return this->finite_element_output.shape_gradients[row][q_point];
5692 template <
int dim,
int spacedim>
5695 const unsigned int q_point)
const
5700 Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
5701 Assert(present_cell.is_initialized(), ExcNotReinited());
5704 if (fe->is_primitive())
5705 return this->finite_element_output.shape_hessians[i][q_point];
5716 const unsigned int row =
5717 this->finite_element_output
5718 .shape_function_to_row_table[i * fe->n_components() +
5719 fe->system_to_component_index(i).first];
5720 return this->finite_element_output.shape_hessians[row][q_point];
5726 template <
int dim,
int spacedim>
5729 const unsigned int i,
5730 const unsigned int q_point,
5731 const unsigned int component)
const
5737 Assert(present_cell.is_initialized(), ExcNotReinited());
5741 if (fe->get_nonzero_components(i)[component] ==
false)
5747 const unsigned int row =
5748 this->finite_element_output
5749 .shape_function_to_row_table[i * fe->n_components() + component];
5750 return this->finite_element_output.shape_hessians[row][q_point];
5755 template <
int dim,
int spacedim>
5758 const unsigned int i,
5759 const unsigned int q_point)
const
5764 Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
5765 Assert(present_cell.is_initialized(), ExcNotReinited());
5768 if (fe->is_primitive())
5769 return this->finite_element_output.shape_3rd_derivatives[i][q_point];
5780 const unsigned int row =
5781 this->finite_element_output
5782 .shape_function_to_row_table[i * fe->n_components() +
5783 fe->system_to_component_index(i).first];
5784 return this->finite_element_output.shape_3rd_derivatives[row][q_point];
5790 template <
int dim,
int spacedim>
5793 const unsigned int i,
5794 const unsigned int q_point,
5795 const unsigned int component)
const
5801 Assert(present_cell.is_initialized(), ExcNotReinited());
5805 if (fe->get_nonzero_components(i)[component] ==
false)
5811 const unsigned int row =
5812 this->finite_element_output
5813 .shape_function_to_row_table[i * fe->n_components() + component];
5814 return this->finite_element_output.shape_3rd_derivatives[row][q_point];
5819 template <
int dim,
int spacedim>
5828 template <
int dim,
int spacedim>
5837 template <
int dim,
int spacedim>
5841 return this->update_flags;
5846 template <
int dim,
int spacedim>
5847 inline const std::vector<Point<spacedim>> &
5852 Assert(present_cell.is_initialized(), ExcNotReinited());
5853 return this->mapping_output.quadrature_points;
5858 template <
int dim,
int spacedim>
5859 inline const std::vector<double> &
5864 Assert(present_cell.is_initialized(), ExcNotReinited());
5865 return this->mapping_output.JxW_values;
5870 template <
int dim,
int spacedim>
5871 inline const std::vector<DerivativeForm<1, dim, spacedim>> &
5876 Assert(present_cell.is_initialized(), ExcNotReinited());
5877 return this->mapping_output.jacobians;
5882 template <
int dim,
int spacedim>
5883 inline const std::vector<DerivativeForm<2, dim, spacedim>> &
5888 Assert(present_cell.is_initialized(), ExcNotReinited());
5889 return this->mapping_output.jacobian_grads;
5894 template <
int dim,
int spacedim>
5897 const unsigned int q_point)
const
5901 Assert(present_cell.is_initialized(), ExcNotReinited());
5902 return this->mapping_output.jacobian_pushed_forward_grads[q_point];
5907 template <
int dim,
int spacedim>
5908 inline const std::vector<Tensor<3, spacedim>> &
5913 Assert(present_cell.is_initialized(), ExcNotReinited());
5914 return this->mapping_output.jacobian_pushed_forward_grads;
5919 template <
int dim,
int spacedim>
5922 const unsigned int q_point)
const
5926 Assert(present_cell.is_initialized(), ExcNotReinited());
5927 return this->mapping_output.jacobian_2nd_derivatives[q_point];
5932 template <
int dim,
int spacedim>
5933 inline const std::vector<DerivativeForm<3, dim, spacedim>> &
5938 Assert(present_cell.is_initialized(), ExcNotReinited());
5939 return this->mapping_output.jacobian_2nd_derivatives;
5944 template <
int dim,
int spacedim>
5947 const unsigned int q_point)
const
5951 "update_jacobian_pushed_forward_2nd_derivatives"));
5952 Assert(present_cell.is_initialized(), ExcNotReinited());
5953 return this->mapping_output.jacobian_pushed_forward_2nd_derivatives[q_point];
5958 template <
int dim,
int spacedim>
5959 inline const std::vector<Tensor<4, spacedim>> &
5964 "update_jacobian_pushed_forward_2nd_derivatives"));
5965 Assert(present_cell.is_initialized(), ExcNotReinited());
5966 return this->mapping_output.jacobian_pushed_forward_2nd_derivatives;
5971 template <
int dim,
int spacedim>
5974 const unsigned int q_point)
const
5978 Assert(present_cell.is_initialized(), ExcNotReinited());
5979 return this->mapping_output.jacobian_3rd_derivatives[q_point];
5984 template <
int dim,
int spacedim>
5985 inline const std::vector<DerivativeForm<4, dim, spacedim>> &
5990 Assert(present_cell.is_initialized(), ExcNotReinited());
5991 return this->mapping_output.jacobian_3rd_derivatives;
5996 template <
int dim,
int spacedim>
5999 const unsigned int q_point)
const
6003 "update_jacobian_pushed_forward_3rd_derivatives"));
6004 Assert(present_cell.is_initialized(), ExcNotReinited());
6005 return this->mapping_output.jacobian_pushed_forward_3rd_derivatives[q_point];
6010 template <
int dim,
int spacedim>
6011 inline const std::vector<Tensor<5, spacedim>> &
6016 "update_jacobian_pushed_forward_3rd_derivatives"));
6017 Assert(present_cell.is_initialized(), ExcNotReinited());
6018 return this->mapping_output.jacobian_pushed_forward_3rd_derivatives;
6023 template <
int dim,
int spacedim>
6024 inline const std::vector<DerivativeForm<1, spacedim, dim>> &
6029 Assert(present_cell.is_initialized(), ExcNotReinited());
6030 return this->mapping_output.inverse_jacobians;
6035 template <
int dim,
int spacedim>
6039 return {0
U, dofs_per_cell};
6044 template <
int dim,
int spacedim>
6047 const unsigned int start_dof_index)
const
6049 Assert(start_dof_index <= dofs_per_cell,
6051 return {start_dof_index, dofs_per_cell};
6056 template <
int dim,
int spacedim>
6059 const unsigned int end_dof_index)
const
6061 Assert(end_dof_index < dofs_per_cell,
6063 return {0
U, end_dof_index + 1};
6068 template <
int dim,
int spacedim>
6072 return {0
U, n_quadrature_points};
6077 template <
int dim,
int spacedim>
6084 Assert(present_cell.is_initialized(), ExcNotReinited());
6086 return this->mapping_output.quadrature_points[q_point];
6091 template <
int dim,
int spacedim>
6098 Assert(present_cell.is_initialized(), ExcNotReinited());
6100 return this->mapping_output.JxW_values[q_point];
6105 template <
int dim,
int spacedim>
6112 Assert(present_cell.is_initialized(), ExcNotReinited());
6114 return this->mapping_output.jacobians[q_point];
6119 template <
int dim,
int spacedim>
6126 Assert(present_cell.is_initialized(), ExcNotReinited());
6128 return this->mapping_output.jacobian_grads[q_point];
6133 template <
int dim,
int spacedim>
6140 Assert(present_cell.is_initialized(), ExcNotReinited());
6142 return this->mapping_output.inverse_jacobians[q_point];
6147 template <
int dim,
int spacedim>
6153 "update_normal_vectors")));
6155 Assert(present_cell.is_initialized(), ExcNotReinited());
6157 return this->mapping_output.normal_vectors[q_point];
6165 template <
int dim,
int spacedim>
6174 template <
int dim,
int spacedim>
6185 template <
int dim,
int spacedim>
6189 return present_face_no;
6193 template <
int dim,
int spacedim>
6197 return present_face_index;
6203 template <
int dim,
int spacedim>
6207 return quadrature[quadrature.size() == 1 ? 0 : present_face_no];
6212 template <
int dim,
int spacedim>
6221 template <
int dim,
int spacedim>
6230 template <
int dim,
int spacedim>
6237 "update_boundary_forms")));
6239 return this->mapping_output.boundary_forms[q_point];
unsigned int get_face_index() const
unsigned int present_face_no
const Tensor< 1, spacedim > & boundary_form(const unsigned int q_point) const
const Quadrature< dim - 1 > & get_quadrature() const
unsigned int present_face_index
unsigned int get_face_number() const
const hp::QCollection< dim - 1 > quadrature
void reinit(const TriaIterator< DoFCellAccessor< dim, spacedim, level_dof_access >> &cell, const unsigned int face_no)
void reinit(const TriaIterator< DoFCellAccessor< dim, spacedim, level_dof_access >> &cell, const typename Triangulation< dim, spacedim >::face_iterator &face)
const FEFaceValues< dim, spacedim > & get_present_fe_values() const
void reinit(const TriaIterator< DoFCellAccessor< dim, spacedim, level_dof_access >> &cell, const unsigned int face_no, const unsigned int subface_no)
void reinit(const TriaIterator< DoFCellAccessor< dim, spacedim, level_dof_access >> &cell, const typename Triangulation< dim, spacedim >::face_iterator &face, const typename Triangulation< dim, spacedim >::face_iterator &subface)
const FESubfaceValues< dim, spacedim > & get_present_fe_values() const
CellIteratorContainer(const TriaIterator< DoFCellAccessor< dim, spacedim, lda >> &cell)
const DoFHandler< dim, spacedim > * dof_handler
Triangulation< dim, spacedim >::cell_iterator cell
CellSimilarity::Similarity cell_similarity
const double & shape_value(const unsigned int i, const unsigned int q_point) const
Tensor< 2, spacedim > shape_hessian_component(const unsigned int i, const unsigned int q_point, const unsigned int component) const
const Tensor< 4, spacedim > & jacobian_pushed_forward_2nd_derivative(const unsigned int q_point) const
CellIteratorContainer present_cell
const Tensor< 2, spacedim > & shape_hessian(const unsigned int i, const unsigned int q_point) const
const Point< spacedim > & quadrature_point(const unsigned int q_point) const
::internal::FEValuesViews::Cache< dim, spacedim > fe_values_views_cache
FEValuesBase(const FEValuesBase &)=delete
const DerivativeForm< 1, spacedim, dim > & inverse_jacobian(const unsigned int q_point) const
boost::signals2::connection tria_listener_mesh_transform
const Tensor< 3, spacedim > & shape_3rd_derivative(const unsigned int i, const unsigned int q_point) const
const Tensor< 1, spacedim > & normal_vector(const unsigned int q_point) const
std_cxx20::ranges::iota_view< unsigned int, unsigned int > dof_indices_ending_at(const unsigned int end_dof_index) const
const FEValuesViews::Vector< dim, spacedim > & operator[](const FEValuesExtractors::Vector &vector) const
internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > mapping_output
const Tensor< 3, spacedim > & jacobian_pushed_forward_grad(const unsigned int q_point) const
const std::vector< DerivativeForm< 1, spacedim, dim > > & get_inverse_jacobians() const
const SmartPointer< const Mapping< dim, spacedim >, FEValuesBase< dim, spacedim > > mapping
UpdateFlags get_update_flags() const
const std::vector< Tensor< 3, spacedim > > & get_jacobian_pushed_forward_grads() const
const unsigned int dofs_per_cell
Tensor< 3, spacedim > shape_3rd_derivative_component(const unsigned int i, const unsigned int q_point, const unsigned int component) const
const SmartPointer< const FiniteElement< dim, spacedim >, FEValuesBase< dim, spacedim > > fe
const DerivativeForm< 1, dim, spacedim > & jacobian(const unsigned int q_point) const
const std::vector< double > & get_JxW_values() const
const FEValuesViews::Scalar< dim, spacedim > & operator[](const FEValuesExtractors::Scalar &scalar) const
const DerivativeForm< 2, dim, spacedim > & jacobian_grad(const unsigned int q_point) const
double shape_value_component(const unsigned int i, const unsigned int q_point, const unsigned int component) const
std_cxx20::ranges::iota_view< unsigned int, unsigned int > dof_indices_starting_at(const unsigned int start_dof_index) const
FEValuesBase & operator=(const FEValuesBase &)=delete
const unsigned int n_quadrature_points
const Tensor< 1, spacedim > & shape_grad(const unsigned int i, const unsigned int q_point) const
const FEValuesViews::SymmetricTensor< 2, dim, spacedim > & operator[](const FEValuesExtractors::SymmetricTensor< 2 > &tensor) const
std_cxx20::ranges::iota_view< unsigned int, unsigned int > dof_indices() const
const std::vector< DerivativeForm< 2, dim, spacedim > > & get_jacobian_grads() const
std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > mapping_data
boost::signals2::connection tria_listener_refinement
const std::vector< DerivativeForm< 3, dim, spacedim > > & get_jacobian_2nd_derivatives() const
std_cxx20::ranges::iota_view< unsigned int, unsigned int > quadrature_point_indices() const
const DerivativeForm< 4, dim, spacedim > & jacobian_3rd_derivative(const unsigned int q_point) const
const std::vector< Tensor< 4, spacedim > > & get_jacobian_pushed_forward_2nd_derivatives() const
const std::vector< DerivativeForm< 1, dim, spacedim > > & get_jacobians() const
::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > finite_element_output
const Tensor< 5, spacedim > & jacobian_pushed_forward_3rd_derivative(const unsigned int q_point) const
const DerivativeForm< 3, dim, spacedim > & jacobian_2nd_derivative(const unsigned int q_point) const
Tensor< 1, spacedim > shape_grad_component(const unsigned int i, const unsigned int q_point, const unsigned int component) const
const Mapping< dim, spacedim > & get_mapping() const
const std::vector< Tensor< 5, spacedim > > & get_jacobian_pushed_forward_3rd_derivatives() const
const FiniteElement< dim, spacedim > & get_fe() const
const std::vector< DerivativeForm< 4, dim, spacedim > > & get_jacobian_3rd_derivatives() const
const std::vector< Point< spacedim > > & get_quadrature_points() const
std::unique_ptr< typename FiniteElement< dim, spacedim >::InternalDataBase > fe_data
double JxW(const unsigned int q_point) const
const FEValuesViews::Tensor< 2, dim, spacedim > & operator[](const FEValuesExtractors::Tensor< 2 > &tensor) const
const unsigned int max_n_quadrature_points
typename ProductType< Number, hessian_type >::type solution_hessian_type
value_type value(const unsigned int shape_function, const unsigned int q_point) const
void get_function_values_from_local_dof_values(const InputVector &dof_values, std::vector< solution_value_type< typename InputVector::value_type >> &values) const
void get_function_values(const InputVector &fe_function, std::vector< solution_value_type< typename InputVector::value_type >> &values) const
const unsigned int component
void get_function_gradients_from_local_dof_values(const InputVector &dof_values, std::vector< solution_gradient_type< typename InputVector::value_type >> &gradients) const
void get_function_third_derivatives(const InputVector &fe_function, std::vector< solution_third_derivative_type< typename InputVector::value_type >> &third_derivatives) const
gradient_type gradient(const unsigned int shape_function, const unsigned int q_point) const
std::vector< ShapeFunctionData > shape_function_data
::Tensor< 1, spacedim > gradient_type
typename ProductType< Number, value_type >::type solution_laplacian_type
third_derivative_type third_derivative(const unsigned int shape_function, const unsigned int q_point) const
typename ProductType< Number, third_derivative_type >::type solution_third_derivative_type
void get_function_gradients(const InputVector &fe_function, std::vector< solution_gradient_type< typename InputVector::value_type >> &gradients) const
typename ProductType< Number, value_type >::type solution_value_type
::Tensor< 2, spacedim > hessian_type
typename ProductType< Number, gradient_type >::type solution_gradient_type
void get_function_hessians_from_local_dof_values(const InputVector &dof_values, std::vector< solution_hessian_type< typename InputVector::value_type >> &hessians) const
Scalar & operator=(const Scalar< dim, spacedim > &)=delete
const SmartPointer< const FEValuesBase< dim, spacedim > > fe_values
Scalar(Scalar< dim, spacedim > &&)=default
void get_function_laplacians_from_local_dof_values(const InputVector &dof_values, std::vector< solution_laplacian_type< typename InputVector::value_type >> &laplacians) const
Scalar(const Scalar< dim, spacedim > &)=delete
void get_function_laplacians(const InputVector &fe_function, std::vector< solution_laplacian_type< typename InputVector::value_type >> &laplacians) const
void get_function_hessians(const InputVector &fe_function, std::vector< solution_hessian_type< typename InputVector::value_type >> &hessians) const
hessian_type hessian(const unsigned int shape_function, const unsigned int q_point) const
Scalar & operator=(Scalar< dim, spacedim > &&) noexcept=default
void get_function_third_derivatives_from_local_dof_values(const InputVector &dof_values, std::vector< solution_third_derivative_type< typename InputVector::value_type >> &third_derivatives) const
::Tensor< 3, spacedim > third_derivative_type
SymmetricTensor(const SymmetricTensor< 2, dim, spacedim > &)=delete
SymmetricTensor(SymmetricTensor< 2, dim, spacedim > &&)=default
SymmetricTensor & operator=(const SymmetricTensor< 2, dim, spacedim > &)=delete
SymmetricTensor & operator=(SymmetricTensor< 2, dim, spacedim > &&) noexcept=default
typename ProductType< Number, value_type >::type solution_value_type
typename ProductType< Number, divergence_type >::type solution_divergence_type
typename ProductType< Number, value_type >::type solution_value_type
typename ProductType< Number, divergence_type >::type solution_divergence_type
const unsigned int first_tensor_component
typename ProductType< Number, gradient_type >::type solution_gradient_type
divergence_type divergence(const unsigned int shape_function, const unsigned int q_point) const
Tensor(const Tensor< 2, dim, spacedim > &)=delete
const SmartPointer< const FEValuesBase< dim, spacedim > > fe_values
value_type value(const unsigned int shape_function, const unsigned int q_point) const
gradient_type gradient(const unsigned int shape_function, const unsigned int q_point) const
Tensor(Tensor< 2, dim, spacedim > &&)=default
std::vector< ShapeFunctionData > shape_function_data
Tensor & operator=(const Tensor< 2, dim, spacedim > &)=delete
Tensor & operator=(Tensor< 2, dim, spacedim > &&)=default
Vector & operator=(Vector< dim, spacedim > &&)=default
typename ProductType< Number, third_derivative_type >::type solution_third_derivative_type
typename ProductType< Number, divergence_type >::type solution_divergence_type
Vector(const Vector< dim, spacedim > &)=delete
hessian_type hessian(const unsigned int shape_function, const unsigned int q_point) const
typename ProductType< Number, hessian_type >::type solution_hessian_type
const SmartPointer< const FEValuesBase< dim, spacedim > > fe_values
Vector & operator=(const Vector< dim, spacedim > &)=delete
typename ProductType< Number, symmetric_gradient_type >::type solution_symmetric_gradient_type
gradient_type gradient(const unsigned int shape_function, const unsigned int q_point) const
divergence_type divergence(const unsigned int shape_function, const unsigned int q_point) const
third_derivative_type third_derivative(const unsigned int shape_function, const unsigned int q_point) const
typename ProductType< Number, gradient_type >::type solution_gradient_type
typename ProductType< Number, value_type >::type solution_value_type
symmetric_gradient_type symmetric_gradient(const unsigned int shape_function, const unsigned int q_point) const
Vector(Vector< dim, spacedim > &&)=default
typename ::internal::CurlType< spacedim >::type curl_type
const unsigned int first_vector_component
typename ProductType< Number, curl_type >::type solution_curl_type
std::vector< ShapeFunctionData > shape_function_data
value_type value(const unsigned int shape_function, const unsigned int q_point) const
curl_type curl(const unsigned int shape_function, const unsigned int q_point) const
typename ProductType< Number, value_type >::type solution_laplacian_type
void reinit(const TriaIterator< DoFCellAccessor< dim, spacedim, level_dof_access >> &cell)
const Quadrature< dim > quadrature
const FEValues< dim, spacedim > & get_present_fe_values() const
const Quadrature< dim > & get_quadrature() const
static constexpr TableIndices< rank_ > unrolled_to_component_indices(const unsigned int i)
typename Tensor< rank_ - 1, dim, Number >::tensor_type value_type
#define DEAL_II_DEPRECATED
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
#define DeclException0(Exception0)
static ::ExceptionBase & ExcAccessToUninitializedField()
static ::ExceptionBase & ExcInternalError()
#define Assert(cond, exc)
static ::ExceptionBase & ExcNotImplemented()
#define AssertIndexRange(index, range)
#define DeclExceptionMsg(Exception, defaulttext)
static ::ExceptionBase & ExcFENotPrimitive()
static ::ExceptionBase & ExcIndexRange(std::size_t arg1, std::size_t arg2, std::size_t arg3)
#define DeclException1(Exception1, type1, outsequence)
static ::ExceptionBase & ExcMessage(std::string arg1)
TriaIterator< CellAccessor< dim, spacedim > > cell_iterator
@ update_jacobian_pushed_forward_2nd_derivatives
@ update_jacobian_pushed_forward_grads
@ update_hessians
Second derivatives of shape functions.
@ update_jacobian_3rd_derivatives
@ update_values
Shape function values.
@ update_jacobian_grads
Gradient of volume element.
@ update_normal_vectors
Normal vectors.
@ update_3rd_derivatives
Third derivatives of shape functions.
@ update_JxW_values
Transformed quadrature weights.
@ update_jacobians
Volume element.
@ update_inverse_jacobians
Volume element.
@ update_gradients
Shape function gradients.
@ update_quadrature_points
Transformed quadrature points.
@ update_jacobian_pushed_forward_3rd_derivatives
@ update_boundary_forms
Outer normal vector, not normalized.
@ update_jacobian_2nd_derivatives
typename ::internal::FEValuesViews::ViewType< dim, spacedim, Extractor >::type View
std::enable_if_t< std::is_fundamental< T >::value, std::size_t > memory_consumption(const T &t)
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
void reinit(MatrixBlock< MatrixType > &v, const BlockSparsityPattern &p)
boost::integer_range< IncrementableType > iota_view
unsigned int global_dof_index
typename ProductType< Number, typename Scalar< dim, spacedim >::value_type >::type laplacian_type
typename ProductType< Number, typename Scalar< dim, spacedim >::value_type >::type value_type
typename ProductType< Number, typename Scalar< dim, spacedim >::gradient_type >::type gradient_type
typename ProductType< Number, typename Scalar< dim, spacedim >::hessian_type >::type hessian_type
typename ProductType< Number, typename Scalar< dim, spacedim >::third_derivative_type >::type third_derivative_type
bool is_nonzero_shape_function_component
typename ProductType< Number, typename SymmetricTensor< 2, dim, spacedim >::divergence_type >::type divergence_type
typename ProductType< Number, typename SymmetricTensor< 2, dim, spacedim >::value_type >::type value_type
int single_nonzero_component
unsigned int single_nonzero_component_index
typename ProductType< Number, typename Tensor< 2, dim, spacedim >::gradient_type >::type gradient_type
typename ProductType< Number, typename Tensor< 2, dim, spacedim >::value_type >::type value_type
typename ProductType< Number, typename Tensor< 2, dim, spacedim >::divergence_type >::type divergence_type
unsigned int single_nonzero_component_index
int single_nonzero_component
typename ProductType< Number, typename Vector< dim, spacedim >::third_derivative_type >::type third_derivative_type
typename ProductType< Number, typename Vector< dim, spacedim >::value_type >::type laplacian_type
typename ProductType< Number, typename Vector< dim, spacedim >::hessian_type >::type hessian_type
typename ProductType< Number, typename Vector< dim, spacedim >::symmetric_gradient_type >::type symmetric_gradient_type
typename ProductType< Number, typename Vector< dim, spacedim >::gradient_type >::type gradient_type
typename ProductType< Number, typename Vector< dim, spacedim >::curl_type >::type curl_type
typename ProductType< Number, typename Vector< dim, spacedim >::divergence_type >::type divergence_type
typename ProductType< Number, typename Vector< dim, spacedim >::value_type >::type value_type
int single_nonzero_component
unsigned int single_nonzero_component_index
typename internal::ProductTypeImpl< typename std::decay< T >::type, typename std::decay< U >::type >::type type
std::vector<::FEValuesViews::Scalar< dim, spacedim > > scalars
std::vector<::FEValuesViews::Vector< dim, spacedim > > vectors
std::vector<::FEValuesViews::SymmetricTensor< 2, dim, spacedim > > symmetric_second_order_tensors
std::vector<::FEValuesViews::Tensor< 2, dim, spacedim > > second_order_tensors
constexpr DEAL_II_HOST SymmetricTensor< 2, dim, Number > symmetrize(const Tensor< 2, dim, Number > &t)