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,
2513 const unsigned int point_no)
const;
2537 const unsigned int point_no,
2567 const unsigned int quadrature_point)
const;
2587 const unsigned int point_no,
2611 const unsigned int point_no)
const;
2631 const unsigned int point_no,
2655 const unsigned int point_no)
const;
2675 const unsigned int point_no,
2718 template <
class InputVector>
2721 const InputVector & fe_function,
2722 std::vector<typename InputVector::value_type> &
values)
const;
2737 template <
class InputVector>
2740 const InputVector & fe_function,
2799 template <
class InputVector>
2802 const InputVector & fe_function,
2804 std::vector<typename InputVector::value_type> &
values)
const;
2814 template <
class InputVector>
2817 const InputVector & fe_function,
2843 template <
class InputVector>
2846 const InputVector & fe_function,
2849 const bool quadrature_points_fastest)
const;
2891 template <
class InputVector>
2894 const InputVector &fe_function,
2914 template <
class InputVector>
2917 const InputVector &fe_function,
2930 template <
class InputVector>
2933 const InputVector & fe_function,
2946 template <
class InputVector>
2949 const InputVector & fe_function,
2954 const bool quadrature_points_fastest =
false)
const;
2999 template <
class InputVector>
3002 const InputVector &fe_function,
3023 template <
class InputVector>
3026 const InputVector &fe_function,
3030 const bool quadrature_points_fastest =
false)
const;
3040 template <
class InputVector>
3043 const InputVector & fe_function,
3056 template <
class InputVector>
3059 const InputVector & fe_function,
3064 const bool quadrature_points_fastest =
false)
const;
3106 template <
class InputVector>
3109 const InputVector & fe_function,
3110 std::vector<typename InputVector::value_type> &laplacians)
const;
3131 template <
class InputVector>
3134 const InputVector & fe_function,
3145 template <
class InputVector>
3148 const InputVector & fe_function,
3150 std::vector<typename InputVector::value_type> & laplacians)
const;
3160 template <
class InputVector>
3163 const InputVector & fe_function,
3175 template <
class InputVector>
3178 const InputVector & fe_function,
3180 std::vector<std::vector<typename InputVector::value_type>> &laplacians,
3181 const bool quadrature_points_fastest =
false)
const;
3225 template <
class InputVector>
3228 const InputVector &fe_function,
3230 &third_derivatives)
const;
3250 template <
class InputVector>
3253 const InputVector &fe_function,
3256 & third_derivatives,
3257 const bool quadrature_points_fastest =
false)
const;
3267 template <
class InputVector>
3270 const InputVector & fe_function,
3273 &third_derivatives)
const;
3283 template <
class InputVector>
3286 const InputVector & fe_function,
3291 const bool quadrature_points_fastest =
false)
const;
3432 const std::vector<Point<spacedim>> &
3451 JxW(
const unsigned int quadrature_point)
const;
3456 const std::vector<double> &
3474 const std::vector<DerivativeForm<1, dim, spacedim>> &
3493 const std::vector<DerivativeForm<2, dim, spacedim>> &
3513 const std::vector<Tensor<3, spacedim>> &
3532 const std::vector<DerivativeForm<3, dim, spacedim>> &
3546 const unsigned int quadrature_point)
const;
3554 const std::vector<Tensor<4, spacedim>> &
3574 const std::vector<DerivativeForm<4, dim, spacedim>> &
3588 const unsigned int quadrature_point)
const;
3596 const std::vector<Tensor<5, spacedim>> &
3614 const std::vector<DerivativeForm<1, spacedim, dim>> &
3646 const std::vector<Tensor<1, spacedim>> &
3647 get_normal_vectors()
const;
3735 get_cell_similarity()
const;
3755 <<
"You are requesting information from an FEValues/FEFaceValues/FESubfaceValues "
3756 <<
"object for which this kind of information has not been computed. What "
3757 <<
"information these objects compute is determined by the update_* flags you "
3758 <<
"pass to the constructor. Here, the operation you are attempting requires "
3760 <<
"> flag to be set, but it was apparently not specified "
3761 <<
"upon construction.");
3769 "FEValues object is not reinit'ed to any cell");
3779 "The FiniteElement you provided to FEValues and the FiniteElement that belongs "
3780 "to the DoFHandler that provided the cell iterator do not match.");
3788 <<
"The shape function with index " << arg1
3789 <<
" is not primitive, i.e. it is vector-valued and "
3790 <<
"has more than one non-zero vector component. This "
3791 <<
"function cannot be called for these shape functions. "
3792 <<
"Maybe you want to use the same function with the "
3793 <<
"_component suffix?");
3803 "The given FiniteElement is not a primitive element but the requested operation "
3804 "only works for those. See FiniteElement::is_primitive() for more information.");
3818 "You have previously called the FEValues::reinit() function with a "
3819 "cell iterator of type Triangulation<dim,spacedim>::cell_iterator. However, "
3820 "when you do this, you cannot call some functions in the FEValues "
3821 "class, such as the get_function_values/gradients/hessians/third_derivatives "
3822 "functions. If you need these functions, then you need to call "
3823 "FEValues::reinit() with an iterator type that allows to extract "
3824 "degrees of freedom, such as DoFHandler<dim,spacedim>::cell_iterator.");
3848 is_initialized()
const;
3863 n_dofs_for_dof_handler()
const;
3869 template <
typename VectorType>
3871 get_interpolated_dof_values(
3872 const VectorType & in,
3880 get_interpolated_dof_values(
const IndexSet & in,
3921 invalidate_present_cell();
3933 maybe_invalidate_previous_present_cell(
3947 std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
3970 std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
3996 compute_update_flags(
const UpdateFlags update_flags)
const;
4011 check_cell_similarity(
4026 template <
int,
int,
int>
4028 template <
int,
int,
int>
4043 template <
int dim,
int spacedim = dim>
4051 static constexpr
unsigned int integral_dimension = dim;
4098 template <
bool level_dof_access>
4182 template <
int dim,
int spacedim = dim>
4190 static constexpr
unsigned int integral_dimension = dim - 1;
4237 const std::vector<Tensor<1, spacedim>> &
4238 get_boundary_forms()
const;
4302 template <
int dim,
int spacedim = dim>
4310 static constexpr
unsigned int dimension = dim;
4312 static constexpr
unsigned int space_dimension = spacedim;
4318 static constexpr
unsigned int integral_dimension = dim - 1;
4363 template <
bool level_dof_access>
4367 const unsigned int face_no);
4375 template <
bool level_dof_access>
4396 const unsigned int face_no);
4448 do_reinit(
const unsigned int face_no);
4468 template <
int dim,
int spacedim = dim>
4475 static constexpr
unsigned int dimension = dim;
4480 static constexpr
unsigned int space_dimension = spacedim;
4486 static constexpr
unsigned int integral_dimension = dim - 1;
4533 template <
bool level_dof_access>
4537 const unsigned int face_no,
4538 const unsigned int subface_no);
4544 template <
bool level_dof_access>
4566 const unsigned int face_no,
4567 const unsigned int subface_no);
4638 do_reinit(
const unsigned int face_no,
const unsigned int subface_no);
4649 template <
int dim,
int spacedim>
4652 const unsigned int q_point)
const
4658 "update_values"))));
4664 return fe_values->finite_element_output.shape_values(
4672 template <
int dim,
int spacedim>
4675 const unsigned int q_point)
const
4680 "update_gradients")));
4695 template <
int dim,
int spacedim>
4698 const unsigned int q_point)
const
4703 "update_hessians")));
4717 template <
int dim,
int spacedim>
4720 const unsigned int q_point)
const
4725 "update_3rd_derivatives")));
4740 template <
int dim,
int spacedim>
4743 const unsigned int q_point)
const
4759 .single_nonzero_component_index] =
4760 fe_values->finite_element_output.shape_values(snc, q_point);
4761 return return_value;
4766 for (
unsigned int d = 0;
d < dim; ++
d)
4768 .is_nonzero_shape_function_component[
d])
4769 return_value[
d] =
fe_values->finite_element_output.shape_values(
4772 return return_value;
4778 template <
int dim,
int spacedim>
4781 const unsigned int q_point)
const
4786 "update_gradients")));
4797 .single_nonzero_component_index] =
4798 fe_values->finite_element_output.shape_gradients[snc][q_point];
4799 return return_value;
4804 for (
unsigned int d = 0;
d < dim; ++
d)
4806 .is_nonzero_shape_function_component[
d])
4808 fe_values->finite_element_output.shape_gradients
4811 return return_value;
4817 template <
int dim,
int spacedim>
4820 const unsigned int q_point)
const
4826 "update_gradients")));
4832 return divergence_type();
4836 .single_nonzero_component_index];
4839 divergence_type return_value = 0;
4840 for (
unsigned int d = 0;
d < dim; ++
d)
4842 .is_nonzero_shape_function_component[
d])
4844 fe_values->finite_element_output.shape_gradients
4847 return return_value;
4853 template <
int dim,
int spacedim>
4856 const unsigned int q_point)
const
4863 "update_gradients")));
4878 "Computing the curl in 1d is not a useful operation"));
4886 curl_type return_value;
4890 .single_nonzero_component_index == 0)
4893 .shape_gradients[snc][q_point][1];
4895 return_value[0] =
fe_values->finite_element_output
4896 .shape_gradients[snc][q_point][0];
4898 return return_value;
4903 curl_type return_value;
4905 return_value[0] = 0.0;
4908 .is_nonzero_shape_function_component[0])
4912 .row_index[0]][q_point][1];
4915 .is_nonzero_shape_function_component[1])
4919 .row_index[1]][q_point][0];
4921 return return_value;
4929 curl_type return_value;
4932 .single_nonzero_component_index)
4936 return_value[0] = 0;
4937 return_value[1] =
fe_values->finite_element_output
4938 .shape_gradients[snc][q_point][2];
4941 .shape_gradients[snc][q_point][1];
4942 return return_value;
4949 .shape_gradients[snc][q_point][2];
4950 return_value[1] = 0;
4951 return_value[2] =
fe_values->finite_element_output
4952 .shape_gradients[snc][q_point][0];
4953 return return_value;
4958 return_value[0] =
fe_values->finite_element_output
4959 .shape_gradients[snc][q_point][1];
4962 .shape_gradients[snc][q_point][0];
4963 return_value[2] = 0;
4964 return return_value;
4971 curl_type return_value;
4973 for (
unsigned int i = 0; i < dim; ++i)
4974 return_value[i] = 0.0;
4977 .is_nonzero_shape_function_component[0])
4982 .row_index[0]][q_point][2];
4986 .row_index[0]][q_point][1];
4990 .is_nonzero_shape_function_component[1])
4995 .row_index[1]][q_point][2];
4999 .row_index[1]][q_point][0];
5003 .is_nonzero_shape_function_component[2])
5008 .row_index[2]][q_point][1];
5012 .row_index[2]][q_point][0];
5015 return return_value;
5026 template <
int dim,
int spacedim>
5029 const unsigned int q_point)
const
5035 "update_hessians")));
5046 .single_nonzero_component_index] =
5047 fe_values->finite_element_output.shape_hessians[snc][q_point];
5048 return return_value;
5053 for (
unsigned int d = 0;
d < dim; ++
d)
5055 .is_nonzero_shape_function_component[
d])
5057 fe_values->finite_element_output.shape_hessians
5060 return return_value;
5066 template <
int dim,
int spacedim>
5069 const unsigned int q_point)
const
5075 "update_3rd_derivatives")));
5086 .single_nonzero_component_index] =
5087 fe_values->finite_element_output.shape_3rd_derivatives[snc][q_point];
5088 return return_value;
5093 for (
unsigned int d = 0;
d < dim; ++
d)
5095 .is_nonzero_shape_function_component[
d])
5097 fe_values->finite_element_output.shape_3rd_derivatives
5100 return return_value;
5112 inline ::SymmetricTensor<2, 1>
5113 symmetrize_single_row(
const unsigned int n, const ::Tensor<1, 1> &t)
5123 inline ::SymmetricTensor<2, 2>
5124 symmetrize_single_row(
const unsigned int n, const ::Tensor<1, 2> &t)
5130 return {{t[0], 0, t[1] / 2}};
5134 return {{0, t[1], t[0] / 2}};
5146 inline ::SymmetricTensor<2, 3>
5147 symmetrize_single_row(
const unsigned int n, const ::Tensor<1, 3> &t)
5153 return {{t[0], 0, 0, t[1] / 2, t[2] / 2, 0}};
5157 return {{0, t[1], 0, t[0] / 2, 0, t[2] / 2}};
5161 return {{0, 0, t[2], 0, t[0] / 2, t[1] / 2}};
5174 template <
int dim,
int spacedim>
5177 const unsigned int q_point)
const
5182 "update_gradients")));
5188 return symmetric_gradient_type();
5190 return internal::symmetrize_single_row(
5192 fe_values->finite_element_output.shape_gradients[snc][q_point]);
5196 for (
unsigned int d = 0;
d < dim; ++
d)
5198 .is_nonzero_shape_function_component[
d])
5200 fe_values->finite_element_output.shape_gradients
5209 template <
int dim,
int spacedim>
5212 const unsigned int q_point)
const
5232 const unsigned int comp =
5234 return_value[value_type::unrolled_to_component_indices(comp)] =
5235 fe_values->finite_element_output.shape_values(snc, q_point);
5236 return return_value;
5241 for (
unsigned int d = 0;
d < value_type::n_independent_components; ++
d)
5243 .is_nonzero_shape_function_component[
d])
5244 return_value[value_type::unrolled_to_component_indices(
d)] =
5245 fe_values->finite_element_output.shape_values(
5247 return return_value;
5253 template <
int dim,
int spacedim>
5256 const unsigned int shape_function,
5257 const unsigned int q_point)
const
5262 "update_gradients")));
5270 return divergence_type();
5293 const unsigned int comp =
5295 const unsigned int ii =
5296 value_type::unrolled_to_component_indices(comp)[0];
5297 const unsigned int jj =
5298 value_type::unrolled_to_component_indices(comp)[1];
5311 const ::Tensor<1, spacedim> &phi_grad =
5312 fe_values->finite_element_output.shape_gradients[snc][q_point];
5314 divergence_type return_value;
5315 return_value[ii] = phi_grad[jj];
5318 return_value[jj] = phi_grad[ii];
5320 return return_value;
5325 divergence_type return_value;
5326 return return_value;
5332 template <
int dim,
int spacedim>
5335 const unsigned int q_point)
const
5355 const unsigned int comp =
5359 return_value[indices] =
5360 fe_values->finite_element_output.shape_values(snc, q_point);
5361 return return_value;
5366 for (
unsigned int d = 0;
d < dim * dim; ++
d)
5368 .is_nonzero_shape_function_component[
d])
5372 return_value[indices] =
5373 fe_values->finite_element_output.shape_values(
5376 return return_value;
5382 template <
int dim,
int spacedim>
5385 const unsigned int q_point)
const
5390 "update_gradients")));
5398 return divergence_type();
5412 const unsigned int comp =
5416 const unsigned int ii = indices[0];
5417 const unsigned int jj = indices[1];
5419 const ::Tensor<1, spacedim> &phi_grad =
5420 fe_values->finite_element_output.shape_gradients[snc][q_point];
5422 divergence_type return_value;
5424 return_value[ii] = phi_grad[jj];
5426 return return_value;
5431 divergence_type return_value;
5432 return return_value;
5438 template <
int dim,
int spacedim>
5441 const unsigned int q_point)
const
5446 "update_gradients")));
5468 const unsigned int comp =
5472 const unsigned int ii = indices[0];
5473 const unsigned int jj = indices[1];
5475 const ::Tensor<1, spacedim> &phi_grad =
5476 fe_values->finite_element_output.shape_gradients[snc][q_point];
5479 return_value[ii][jj] = phi_grad;
5481 return return_value;
5487 return return_value;
5499 template <
int dim,
int spacedim>
5506 , dof_handler(&cell->get_dof_handler())
5507 , level_dof_access(lda)
5512 template <
int dim,
int spacedim>
5519 return fe_values_views_cache.scalars[
scalar.component];
5524 template <
int dim,
int spacedim>
5530 fe_values_views_cache.vectors.size());
5537 template <
int dim,
int spacedim>
5544 fe_values_views_cache.symmetric_second_order_tensors.size(),
5547 fe_values_views_cache.symmetric_second_order_tensors.size()));
5549 return fe_values_views_cache
5555 template <
int dim,
int spacedim>
5561 fe_values_views_cache.second_order_tensors.size());
5563 return fe_values_views_cache
5569 template <
int dim,
int spacedim>
5570 inline const double &
5572 const unsigned int j)
const
5577 Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
5578 Assert(present_cell.is_initialized(), ExcNotReinited());
5581 if (fe->is_primitive())
5582 return this->finite_element_output.shape_values(i, j);
5593 const unsigned int row =
5594 this->finite_element_output
5595 .shape_function_to_row_table[i * fe->n_components() +
5596 fe->system_to_component_index(i).first];
5597 return this->finite_element_output.shape_values(row, j);
5603 template <
int dim,
int spacedim>
5606 const unsigned int i,
5607 const unsigned int j,
5608 const unsigned int component)
const
5614 Assert(present_cell.is_initialized(), ExcNotReinited());
5619 if (fe->get_nonzero_components(i)[component] ==
false)
5625 const unsigned int row =
5626 this->finite_element_output
5627 .shape_function_to_row_table[i * fe->n_components() + component];
5628 return this->finite_element_output.shape_values(row, j);
5633 template <
int dim,
int spacedim>
5636 const unsigned int j)
const
5641 Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
5642 Assert(present_cell.is_initialized(), ExcNotReinited());
5645 if (fe->is_primitive())
5646 return this->finite_element_output.shape_gradients[i][j];
5657 const unsigned int row =
5658 this->finite_element_output
5659 .shape_function_to_row_table[i * fe->n_components() +
5660 fe->system_to_component_index(i).first];
5661 return this->finite_element_output.shape_gradients[row][j];
5667 template <
int dim,
int spacedim>
5670 const unsigned int i,
5671 const unsigned int j,
5672 const unsigned int component)
const
5678 Assert(present_cell.is_initialized(), ExcNotReinited());
5682 if (fe->get_nonzero_components(i)[component] ==
false)
5688 const unsigned int row =
5689 this->finite_element_output
5690 .shape_function_to_row_table[i * fe->n_components() + component];
5691 return this->finite_element_output.shape_gradients[row][j];
5696 template <
int dim,
int spacedim>
5699 const unsigned int j)
const
5704 Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
5705 Assert(present_cell.is_initialized(), ExcNotReinited());
5708 if (fe->is_primitive())
5709 return this->finite_element_output.shape_hessians[i][j];
5720 const unsigned int row =
5721 this->finite_element_output
5722 .shape_function_to_row_table[i * fe->n_components() +
5723 fe->system_to_component_index(i).first];
5724 return this->finite_element_output.shape_hessians[row][j];
5730 template <
int dim,
int spacedim>
5733 const unsigned int i,
5734 const unsigned int j,
5735 const unsigned int component)
const
5741 Assert(present_cell.is_initialized(), ExcNotReinited());
5745 if (fe->get_nonzero_components(i)[component] ==
false)
5751 const unsigned int row =
5752 this->finite_element_output
5753 .shape_function_to_row_table[i * fe->n_components() + component];
5754 return this->finite_element_output.shape_hessians[row][j];
5759 template <
int dim,
int spacedim>
5762 const unsigned int j)
const
5767 Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
5768 Assert(present_cell.is_initialized(), ExcNotReinited());
5771 if (fe->is_primitive())
5772 return this->finite_element_output.shape_3rd_derivatives[i][j];
5783 const unsigned int row =
5784 this->finite_element_output
5785 .shape_function_to_row_table[i * fe->n_components() +
5786 fe->system_to_component_index(i).first];
5787 return this->finite_element_output.shape_3rd_derivatives[row][j];
5793 template <
int dim,
int spacedim>
5796 const unsigned int i,
5797 const unsigned int j,
5798 const unsigned int component)
const
5804 Assert(present_cell.is_initialized(), ExcNotReinited());
5808 if (fe->get_nonzero_components(i)[component] ==
false)
5814 const unsigned int row =
5815 this->finite_element_output
5816 .shape_function_to_row_table[i * fe->n_components() + component];
5817 return this->finite_element_output.shape_3rd_derivatives[row][j];
5822 template <
int dim,
int spacedim>
5831 template <
int dim,
int spacedim>
5840 template <
int dim,
int spacedim>
5844 return this->update_flags;
5849 template <
int dim,
int spacedim>
5850 inline const std::vector<Point<spacedim>> &
5855 Assert(present_cell.is_initialized(), ExcNotReinited());
5856 return this->mapping_output.quadrature_points;
5861 template <
int dim,
int spacedim>
5862 inline const std::vector<double> &
5867 Assert(present_cell.is_initialized(), ExcNotReinited());
5868 return this->mapping_output.JxW_values;
5873 template <
int dim,
int spacedim>
5874 inline const std::vector<DerivativeForm<1, dim, spacedim>> &
5879 Assert(present_cell.is_initialized(), ExcNotReinited());
5880 return this->mapping_output.jacobians;
5885 template <
int dim,
int spacedim>
5886 inline const std::vector<DerivativeForm<2, dim, spacedim>> &
5891 Assert(present_cell.is_initialized(), ExcNotReinited());
5892 return this->mapping_output.jacobian_grads;
5897 template <
int dim,
int spacedim>
5900 const unsigned int i)
const
5904 Assert(present_cell.is_initialized(), ExcNotReinited());
5905 return this->mapping_output.jacobian_pushed_forward_grads[i];
5910 template <
int dim,
int spacedim>
5911 inline const std::vector<Tensor<3, spacedim>> &
5916 Assert(present_cell.is_initialized(), ExcNotReinited());
5917 return this->mapping_output.jacobian_pushed_forward_grads;
5922 template <
int dim,
int spacedim>
5928 Assert(present_cell.is_initialized(), ExcNotReinited());
5929 return this->mapping_output.jacobian_2nd_derivatives[i];
5934 template <
int dim,
int spacedim>
5935 inline const std::vector<DerivativeForm<3, dim, spacedim>> &
5940 Assert(present_cell.is_initialized(), ExcNotReinited());
5941 return this->mapping_output.jacobian_2nd_derivatives;
5946 template <
int dim,
int spacedim>
5949 const unsigned int i)
const
5953 "update_jacobian_pushed_forward_2nd_derivatives"));
5954 Assert(present_cell.is_initialized(), ExcNotReinited());
5955 return this->mapping_output.jacobian_pushed_forward_2nd_derivatives[i];
5960 template <
int dim,
int spacedim>
5961 inline const std::vector<Tensor<4, spacedim>> &
5966 "update_jacobian_pushed_forward_2nd_derivatives"));
5967 Assert(present_cell.is_initialized(), ExcNotReinited());
5968 return this->mapping_output.jacobian_pushed_forward_2nd_derivatives;
5973 template <
int dim,
int spacedim>
5979 Assert(present_cell.is_initialized(), ExcNotReinited());
5980 return this->mapping_output.jacobian_3rd_derivatives[i];
5985 template <
int dim,
int spacedim>
5986 inline const std::vector<DerivativeForm<4, dim, spacedim>> &
5991 Assert(present_cell.is_initialized(), ExcNotReinited());
5992 return this->mapping_output.jacobian_3rd_derivatives;
5997 template <
int dim,
int spacedim>
6000 const unsigned int i)
const
6004 "update_jacobian_pushed_forward_3rd_derivatives"));
6005 Assert(present_cell.is_initialized(), ExcNotReinited());
6006 return this->mapping_output.jacobian_pushed_forward_3rd_derivatives[i];
6011 template <
int dim,
int spacedim>
6012 inline const std::vector<Tensor<5, spacedim>> &
6017 "update_jacobian_pushed_forward_3rd_derivatives"));
6018 Assert(present_cell.is_initialized(), ExcNotReinited());
6019 return this->mapping_output.jacobian_pushed_forward_3rd_derivatives;
6024 template <
int dim,
int spacedim>
6025 inline const std::vector<DerivativeForm<1, spacedim, dim>> &
6030 Assert(present_cell.is_initialized(), ExcNotReinited());
6031 return this->mapping_output.inverse_jacobians;
6036 template <
int dim,
int spacedim>
6040 return {0
U, dofs_per_cell};
6045 template <
int dim,
int spacedim>
6048 const unsigned int start_dof_index)
const
6050 Assert(start_dof_index <= dofs_per_cell,
6052 return {start_dof_index, dofs_per_cell};
6057 template <
int dim,
int spacedim>
6060 const unsigned int end_dof_index)
const
6062 Assert(end_dof_index < dofs_per_cell,
6064 return {0
U, end_dof_index + 1};
6069 template <
int dim,
int spacedim>
6073 return {0
U, n_quadrature_points};
6078 template <
int dim,
int spacedim>
6085 Assert(present_cell.is_initialized(), ExcNotReinited());
6087 return this->mapping_output.quadrature_points[i];
6092 template <
int dim,
int spacedim>
6099 Assert(present_cell.is_initialized(), ExcNotReinited());
6101 return this->mapping_output.JxW_values[i];
6106 template <
int dim,
int spacedim>
6113 Assert(present_cell.is_initialized(), ExcNotReinited());
6115 return this->mapping_output.jacobians[i];
6120 template <
int dim,
int spacedim>
6127 Assert(present_cell.is_initialized(), ExcNotReinited());
6129 return this->mapping_output.jacobian_grads[i];
6134 template <
int dim,
int spacedim>
6141 Assert(present_cell.is_initialized(), ExcNotReinited());
6143 return this->mapping_output.inverse_jacobians[i];
6148 template <
int dim,
int spacedim>
6154 "update_normal_vectors")));
6156 Assert(present_cell.is_initialized(), ExcNotReinited());
6158 return this->mapping_output.normal_vectors[i];
6166 template <
int dim,
int spacedim>
6175 template <
int dim,
int spacedim>
6186 template <
int dim,
int spacedim>
6190 return present_face_no;
6194 template <
int dim,
int spacedim>
6198 return present_face_index;
6204 template <
int dim,
int spacedim>
6208 return quadrature[quadrature.size() == 1 ? 0 : present_face_no];
6213 template <
int dim,
int spacedim>
6222 template <
int dim,
int spacedim>
6231 template <
int dim,
int spacedim>
6238 "update_boundary_forms")));
6240 return this->mapping_output.boundary_forms[i];
unsigned int get_face_index() const
unsigned int present_face_no
const Quadrature< dim - 1 > & get_quadrature() const
unsigned int present_face_index
const Tensor< 1, spacedim > & boundary_form(const unsigned int i) const
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
CellIteratorContainer present_cell
const double & shape_value(const unsigned int function_no, const unsigned int point_no) const
::internal::FEValuesViews::Cache< dim, spacedim > fe_values_views_cache
FEValuesBase(const FEValuesBase &)=delete
const DerivativeForm< 2, dim, spacedim > & jacobian_grad(const unsigned int quadrature_point) const
const DerivativeForm< 4, dim, spacedim > & jacobian_3rd_derivative(const unsigned int quadrature_point) const
boost::signals2::connection tria_listener_mesh_transform
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< 1, spacedim > & shape_grad(const unsigned int function_no, const unsigned int quadrature_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 DerivativeForm< 3, dim, spacedim > & jacobian_2nd_derivative(const unsigned int quadrature_point) const
const unsigned int dofs_per_cell
const SmartPointer< const FiniteElement< dim, spacedim >, FEValuesBase< dim, spacedim > > fe
const DerivativeForm< 1, dim, spacedim > & jacobian(const unsigned int quadrature_point) const
Tensor< 1, spacedim > shape_grad_component(const unsigned int function_no, const unsigned int point_no, const unsigned int component) const
const std::vector< double > & get_JxW_values() const
const FEValuesViews::Scalar< dim, spacedim > & operator[](const FEValuesExtractors::Scalar &scalar) const
Tensor< 3, spacedim > shape_3rd_derivative_component(const unsigned int function_no, const unsigned int point_no, 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< 2, spacedim > & shape_hessian(const unsigned int function_no, const unsigned int point_no) const
Tensor< 2, spacedim > shape_hessian_component(const unsigned int function_no, const unsigned int point_no, const unsigned int component) 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
const DerivativeForm< 1, spacedim, dim > & inverse_jacobian(const unsigned int quadrature_point) const
double shape_value_component(const unsigned int function_no, const unsigned int point_no, const unsigned int component) const
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 Point< spacedim > & quadrature_point(const unsigned int q) const
const Tensor< 3, spacedim > & jacobian_pushed_forward_grad(const unsigned int quadrature_point) const
double JxW(const unsigned int quadrature_point) const
const std::vector< Tensor< 4, spacedim > > & get_jacobian_pushed_forward_2nd_derivatives() const
const Tensor< 1, spacedim > & normal_vector(const unsigned int i) const
const Tensor< 5, spacedim > & jacobian_pushed_forward_3rd_derivative(const unsigned int quadrature_point) const
const std::vector< DerivativeForm< 1, dim, spacedim > > & get_jacobians() const
::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > finite_element_output
const Tensor< 4, spacedim > & jacobian_pushed_forward_2nd_derivative(const unsigned int quadrature_point) const
const Tensor< 3, spacedim > & shape_3rd_derivative(const unsigned int function_no, const unsigned int point_no) 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
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)