16 #ifndef dealii_fe_values_base_h
17 #define dealii_fe_values_base_h
49 #include <type_traits>
152 template <
int dim,
int spacedim>
274 shape_value(
const unsigned int i,
const unsigned int q_point)
const;
298 const unsigned int q_point,
299 const unsigned int component)
const;
327 shape_grad(
const unsigned int i,
const unsigned int q_point)
const;
347 const unsigned int q_point,
348 const unsigned int component)
const;
390 const unsigned int q_point,
391 const unsigned int component)
const;
433 const unsigned int q_point,
434 const unsigned int component)
const;
484 template <
typename Number>
487 std::vector<Number> &
values)
const;
502 template <
typename Number>
563 template <
typename Number>
567 std::vector<Number> &
values)
const;
577 template <
typename Number>
605 template <
typename Number>
610 const bool quadrature_points_fastest)
const;
658 template <
typename Number>
680 template <
typename Number>
694 template <
typename Number>
709 template <
typename Number>
715 const bool quadrature_points_fastest =
false)
const;
760 template <
typename Number>
783 template <
typename Number>
788 const bool quadrature_points_fastest =
false)
const;
798 template <
typename Number>
813 template <
typename Number>
819 const bool quadrature_points_fastest =
false)
const;
861 template <
typename Number>
864 std::vector<Number> &laplacians)
const;
885 template <
typename Number>
898 template <
typename Number>
903 std::vector<Number> &laplacians)
const;
913 template <
typename Number>
928 template <
typename Number>
933 std::vector<std::vector<Number>> &laplacians,
934 const bool quadrature_points_fastest =
false)
const;
978 template <
typename Number>
1002 template <
typename Number>
1007 const bool quadrature_points_fastest =
false)
const;
1017 template <
typename Number>
1032 template <
typename Number>
1038 const bool quadrature_points_fastest =
false)
const;
1180 const std::vector<Point<spacedim>> &
1200 JxW(
const unsigned int q_point)
const;
1205 const std::vector<double> &
1223 const std::vector<DerivativeForm<1, dim, spacedim>> &
1242 const std::vector<DerivativeForm<2, dim, spacedim>> &
1262 const std::vector<Tensor<3, spacedim>> &
1281 const std::vector<DerivativeForm<3, dim, spacedim>> &
1302 const std::vector<Tensor<4, spacedim>> &
1322 const std::vector<DerivativeForm<4, dim, spacedim>> &
1343 const std::vector<Tensor<5, spacedim>> &
1361 const std::vector<DerivativeForm<1, spacedim, dim>> &
1393 const std::vector<Tensor<1, spacedim>> &
1502 <<
"You are requesting information from an FEValues/FEFaceValues/FESubfaceValues "
1503 <<
"object for which this kind of information has not been computed. What "
1504 <<
"information these objects compute is determined by the update_* flags you "
1505 <<
"pass to the constructor. Here, the operation you are attempting requires "
1507 <<
"> flag to be set, but it was apparently not specified "
1508 <<
"upon construction.");
1516 "FEValues object is not reinit'ed to any cell");
1526 "The FiniteElement you provided to FEValues and the FiniteElement that belongs "
1527 "to the DoFHandler that provided the cell iterator do not match.");
1535 <<
"The shape function with index " << arg1
1536 <<
" is not primitive, i.e. it is vector-valued and "
1537 <<
"has more than one non-zero vector component. This "
1538 <<
"function cannot be called for these shape functions. "
1539 <<
"Maybe you want to use the same function with the "
1540 <<
"_component suffix?");
1550 "The given FiniteElement is not a primitive element but the requested operation "
1551 "only works for those. See FiniteElement::is_primitive() for more information.");
1565 "You have previously called the FEValues::reinit() function with a "
1566 "cell iterator of type Triangulation<dim,spacedim>::cell_iterator. However, "
1567 "when you do this, you cannot call some functions in the FEValues "
1568 "class, such as the get_function_values/gradients/hessians/third_derivatives "
1569 "functions. If you need these functions, then you need to call "
1570 "FEValues::reinit() with an iterator type that allows to extract "
1571 "degrees of freedom, such as DoFHandler<dim,spacedim>::cell_iterator.");
1616 template <
typename Number>
1693 std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
1716 std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
1777 template <
int,
int,
int>
1779 template <
int,
int,
int>
1787 template <
int dim,
int spacedim>
1794 , dof_handler(&cell->get_dof_handler())
1795 , level_dof_access(lda)
1800 template <
int dim,
int spacedim>
1807 return fe_values_views_cache.
scalars[
scalar.component].value_or_initialize(
1815 template <
int dim,
int spacedim>
1821 fe_values_views_cache.
vectors.size());
1824 .value_or_initialize([vector,
this]() {
1832 template <
int dim,
int spacedim>
1844 return fe_values_views_cache
1846 .value_or_initialize([tensor,
this]() {
1854 template <
int dim,
int spacedim>
1862 return fe_values_views_cache
1864 .value_or_initialize([tensor,
this]() {
1872 template <
int dim,
int spacedim>
1873 inline const double &
1875 const unsigned int q_point)
const
1880 Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
1884 if (fe->is_primitive())
1885 return this->finite_element_output.shape_values(i, q_point);
1896 const unsigned int row =
1897 this->finite_element_output
1898 .shape_function_to_row_table[i * fe->n_components() +
1899 fe->system_to_component_index(i).first];
1900 return this->finite_element_output.shape_values(row, q_point);
1906 template <
int dim,
int spacedim>
1909 const unsigned int i,
1910 const unsigned int q_point,
1911 const unsigned int component)
const
1922 if (fe->get_nonzero_components(i)[component] ==
false)
1928 const unsigned int row =
1929 this->finite_element_output
1930 .shape_function_to_row_table[i * fe->n_components() + component];
1931 return this->finite_element_output.shape_values(row, q_point);
1936 template <
int dim,
int spacedim>
1939 const unsigned int q_point)
const
1944 Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
1948 if (fe->is_primitive())
1949 return this->finite_element_output.shape_gradients[i][q_point];
1960 const unsigned int row =
1961 this->finite_element_output
1962 .shape_function_to_row_table[i * fe->n_components() +
1963 fe->system_to_component_index(i).first];
1964 return this->finite_element_output.shape_gradients[row][q_point];
1970 template <
int dim,
int spacedim>
1973 const unsigned int i,
1974 const unsigned int q_point,
1975 const unsigned int component)
const
1985 if (fe->get_nonzero_components(i)[component] ==
false)
1991 const unsigned int row =
1992 this->finite_element_output
1993 .shape_function_to_row_table[i * fe->n_components() + component];
1994 return this->finite_element_output.shape_gradients[row][q_point];
1999 template <
int dim,
int spacedim>
2002 const unsigned int q_point)
const
2007 Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
2011 if (fe->is_primitive())
2012 return this->finite_element_output.shape_hessians[i][q_point];
2023 const unsigned int row =
2024 this->finite_element_output
2025 .shape_function_to_row_table[i * fe->n_components() +
2026 fe->system_to_component_index(i).first];
2027 return this->finite_element_output.shape_hessians[row][q_point];
2033 template <
int dim,
int spacedim>
2036 const unsigned int i,
2037 const unsigned int q_point,
2038 const unsigned int component)
const
2048 if (fe->get_nonzero_components(i)[component] ==
false)
2054 const unsigned int row =
2055 this->finite_element_output
2056 .shape_function_to_row_table[i * fe->n_components() + component];
2057 return this->finite_element_output.shape_hessians[row][q_point];
2062 template <
int dim,
int spacedim>
2065 const unsigned int i,
2066 const unsigned int q_point)
const
2071 Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
2075 if (fe->is_primitive())
2076 return this->finite_element_output.shape_3rd_derivatives[i][q_point];
2087 const unsigned int row =
2088 this->finite_element_output
2089 .shape_function_to_row_table[i * fe->n_components() +
2090 fe->system_to_component_index(i).first];
2091 return this->finite_element_output.shape_3rd_derivatives[row][q_point];
2097 template <
int dim,
int spacedim>
2100 const unsigned int i,
2101 const unsigned int q_point,
2102 const unsigned int component)
const
2112 if (fe->get_nonzero_components(i)[component] ==
false)
2118 const unsigned int row =
2119 this->finite_element_output
2120 .shape_function_to_row_table[i * fe->n_components() + component];
2121 return this->finite_element_output.shape_3rd_derivatives[row][q_point];
2126 template <
int dim,
int spacedim>
2135 template <
int dim,
int spacedim>
2144 template <
int dim,
int spacedim>
2148 return this->update_flags;
2153 template <
int dim,
int spacedim>
2154 inline const std::vector<Point<spacedim>> &
2160 return this->mapping_output.quadrature_points;
2165 template <
int dim,
int spacedim>
2166 inline const std::vector<double> &
2172 return this->mapping_output.JxW_values;
2177 template <
int dim,
int spacedim>
2178 inline const std::vector<DerivativeForm<1, dim, spacedim>> &
2184 return this->mapping_output.jacobians;
2189 template <
int dim,
int spacedim>
2190 inline const std::vector<DerivativeForm<2, dim, spacedim>> &
2196 return this->mapping_output.jacobian_grads;
2201 template <
int dim,
int spacedim>
2204 const unsigned int q_point)
const
2209 return this->mapping_output.jacobian_pushed_forward_grads[q_point];
2214 template <
int dim,
int spacedim>
2215 inline const std::vector<Tensor<3, spacedim>> &
2221 return this->mapping_output.jacobian_pushed_forward_grads;
2226 template <
int dim,
int spacedim>
2229 const unsigned int q_point)
const
2234 return this->mapping_output.jacobian_2nd_derivatives[q_point];
2239 template <
int dim,
int spacedim>
2240 inline const std::vector<DerivativeForm<3, dim, spacedim>> &
2246 return this->mapping_output.jacobian_2nd_derivatives;
2251 template <
int dim,
int spacedim>
2254 const unsigned int q_point)
const
2258 "update_jacobian_pushed_forward_2nd_derivatives"));
2260 return this->mapping_output.jacobian_pushed_forward_2nd_derivatives[q_point];
2265 template <
int dim,
int spacedim>
2266 inline const std::vector<Tensor<4, spacedim>> &
2271 "update_jacobian_pushed_forward_2nd_derivatives"));
2273 return this->mapping_output.jacobian_pushed_forward_2nd_derivatives;
2278 template <
int dim,
int spacedim>
2281 const unsigned int q_point)
const
2286 return this->mapping_output.jacobian_3rd_derivatives[q_point];
2291 template <
int dim,
int spacedim>
2292 inline const std::vector<DerivativeForm<4, dim, spacedim>> &
2298 return this->mapping_output.jacobian_3rd_derivatives;
2303 template <
int dim,
int spacedim>
2306 const unsigned int q_point)
const
2310 "update_jacobian_pushed_forward_3rd_derivatives"));
2312 return this->mapping_output.jacobian_pushed_forward_3rd_derivatives[q_point];
2317 template <
int dim,
int spacedim>
2318 inline const std::vector<Tensor<5, spacedim>> &
2323 "update_jacobian_pushed_forward_3rd_derivatives"));
2325 return this->mapping_output.jacobian_pushed_forward_3rd_derivatives;
2330 template <
int dim,
int spacedim>
2331 inline const std::vector<DerivativeForm<1, spacedim, dim>> &
2337 return this->mapping_output.inverse_jacobians;
2342 template <
int dim,
int spacedim>
2346 return {0
U, dofs_per_cell};
2351 template <
int dim,
int spacedim>
2354 const unsigned int start_dof_index)
const
2356 Assert(start_dof_index <= dofs_per_cell,
2358 return {start_dof_index, dofs_per_cell};
2363 template <
int dim,
int spacedim>
2366 const unsigned int end_dof_index)
const
2368 Assert(end_dof_index < dofs_per_cell,
2370 return {0
U, end_dof_index + 1};
2375 template <
int dim,
int spacedim>
2379 return {0
U, n_quadrature_points};
2384 template <
int dim,
int spacedim>
2393 return this->mapping_output.quadrature_points[q_point];
2398 template <
int dim,
int spacedim>
2407 return this->mapping_output.JxW_values[q_point];
2412 template <
int dim,
int spacedim>
2421 return this->mapping_output.jacobians[q_point];
2426 template <
int dim,
int spacedim>
2435 return this->mapping_output.jacobian_grads[q_point];
2440 template <
int dim,
int spacedim>
2449 return this->mapping_output.inverse_jacobians[q_point];
2454 template <
int dim,
int spacedim>
2460 "update_normal_vectors")));
2464 return this->mapping_output.normal_vectors[q_point];
CellIteratorContainer(const TriaIterator< DoFCellAccessor< dim, spacedim, lda >> &cell)
const DoFHandler< dim, spacedim > * dof_handler
bool is_initialized() const
types::global_dof_index n_dofs_for_dof_handler() const
Triangulation< dim, spacedim >::cell_iterator cell
void get_interpolated_dof_values(const ReadVector< Number > &in, Vector< Number > &out) const
CellSimilarity::Similarity cell_similarity
void get_function_hessians(const ReadVector< Number > &fe_function, std::vector< Tensor< 2, spacedim, Number >> &hessians) const
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
Triangulation< dim, spacedim >::cell_iterator get_cell() const
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
void get_function_values(const ReadVector< Number > &fe_function, std::vector< Number > &values) const
std_cxx20::ranges::iota_view< unsigned int, unsigned int > dof_indices_ending_at(const unsigned int end_dof_index) const
FEValuesBase(const unsigned int n_q_points, const unsigned int dofs_per_cell, const UpdateFlags update_flags, const Mapping< dim, spacedim > &mapping, const FiniteElement< dim, spacedim > &fe)
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
void get_function_third_derivatives(const ReadVector< Number > &fe_function, std::vector< Tensor< 3, spacedim, Number >> &third_derivatives) const
virtual ~FEValuesBase() override
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
void check_cell_similarity(const typename Triangulation< dim, spacedim >::cell_iterator &cell)
Tensor< 3, spacedim > shape_3rd_derivative_component(const unsigned int i, const unsigned int q_point, const unsigned int component) const
static constexpr unsigned int space_dimension
const SmartPointer< const FiniteElement< dim, spacedim >, FEValuesBase< dim, spacedim > > fe
const DerivativeForm< 1, dim, spacedim > & jacobian(const unsigned int q_point) const
void always_allow_check_for_cell_similarity(const bool allow)
void get_function_laplacians(const ReadVector< Number > &fe_function, std::vector< Number > &laplacians) 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
CellSimilarity::Similarity get_cell_similarity() const
const std::vector< Tensor< 1, spacedim > > & get_normal_vectors() const
std::size_t memory_consumption() const
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
UpdateFlags compute_update_flags(const UpdateFlags update_flags) const
std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > mapping_data
void get_function_gradients(const ReadVector< Number > &fe_function, std::vector< Tensor< 1, spacedim, Number >> &gradients) 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
void invalidate_present_cell()
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
bool check_for_cell_similarity_allowed
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
void maybe_invalidate_previous_present_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell)
const FEValuesViews::Tensor< 2, dim, spacedim > & operator[](const FEValuesExtractors::Tensor< 2 > &tensor) const
const unsigned int max_n_quadrature_points
static constexpr unsigned int dimension
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
static ::ExceptionBase & ExcAccessToUninitializedField()
static ::ExceptionBase & ExcNeedsDoFHandler()
static ::ExceptionBase & ExcNotReinited()
static ::ExceptionBase & ExcShapeFunctionNotPrimitive(int arg1)
#define Assert(cond, exc)
static ::ExceptionBase & ExcFENotPrimitive()
static ::ExceptionBase & ExcFEDontMatch()
#define AssertIndexRange(index, range)
#define DeclExceptionMsg(Exception, defaulttext)
static ::ExceptionBase & ExcIndexRange(std::size_t arg1, std::size_t arg2, std::size_t arg3)
#define DeclException1(Exception1, type1, outsequence)
static ::ExceptionBase & ExcAccessToUninitializedField(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_jacobian_2nd_derivatives
boost::integer_range< IncrementableType > iota_view
unsigned int global_dof_index
std::vector< Lazy<::FEValuesViews::Vector< dim, spacedim > > > vectors
std::vector< Lazy<::FEValuesViews::SymmetricTensor< 2, dim, spacedim > > > symmetric_second_order_tensors
std::vector< Lazy<::FEValuesViews::Tensor< 2, dim, spacedim > > > second_order_tensors
std::vector< Lazy<::FEValuesViews::Scalar< dim, spacedim > > > scalars