Loading [MathJax]/extensions/TeX/AMSsymbols.js
 deal.II version GIT relicensing-2846-g6fb608615f 2025-03-15 04:10:00+00:00
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages Concepts
fe_values_base.h
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 2023 by the deal.II authors
5//
6// This file is part of the deal.II library.
7//
8// Part of the source code is dual licensed under Apache-2.0 WITH
9// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
10// governing the source code and code contributions can be found in
11// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
12//
13// ------------------------------------------------------------------------
14
15#ifndef dealii_fe_values_base_h
16#define dealii_fe_values_base_h
17
18
19#include <deal.II/base/config.h>
20
24#include <deal.II/base/point.h>
28
31
32#include <deal.II/fe/fe.h>
36#include <deal.II/fe/mapping.h>
38
39#include <deal.II/grid/tria.h>
41
43
45
46#include <boost/signals2/connection.hpp>
47
48#include <algorithm>
49#include <memory>
50#include <optional>
51#include <type_traits>
52#include <variant>
53
55
155template <int dim, int spacedim>
157{
158public:
162 static constexpr unsigned int dimension = dim;
163
167 static constexpr unsigned int space_dimension = spacedim;
168
176 const unsigned int n_quadrature_points;
177
187 const unsigned int max_n_quadrature_points;
188
194 const unsigned int dofs_per_cell;
195
196
204 FEValuesBase(const unsigned int n_q_points,
205 const unsigned int dofs_per_cell,
209
215 operator=(const FEValuesBase &) = delete;
216
221 FEValuesBase(const FEValuesBase &) = delete;
222
226 virtual ~FEValuesBase() override;
227
248 void
250
254
276 const double &
277 shape_value(const unsigned int i, const unsigned int q_point) const;
278
299 double
300 shape_value_component(const unsigned int i,
301 const unsigned int q_point,
302 const unsigned int component) const;
303
329 const Tensor<1, spacedim> &
330 shape_grad(const unsigned int i, const unsigned int q_point) const;
331
349 shape_grad_component(const unsigned int i,
350 const unsigned int q_point,
351 const unsigned int component) const;
352
372 const Tensor<2, spacedim> &
373 shape_hessian(const unsigned int i, const unsigned int q_point) const;
374
392 shape_hessian_component(const unsigned int i,
393 const unsigned int q_point,
394 const unsigned int component) const;
395
415 const Tensor<3, spacedim> &
416 shape_3rd_derivative(const unsigned int i, const unsigned int q_point) const;
417
435 shape_3rd_derivative_component(const unsigned int i,
436 const unsigned int q_point,
437 const unsigned int component) const;
438
441
487 template <typename Number>
488 void
489 get_function_values(const ReadVector<Number> &fe_function,
490 std::vector<Number> &values) const;
491
505 template <typename Number>
506 void
507 get_function_values(const ReadVector<Number> &fe_function,
508 std::vector<Vector<Number>> &values) const;
509
566 template <typename Number>
567 void
568 get_function_values(const ReadVector<Number> &fe_function,
570 std::vector<Number> &values) const;
571
580 template <typename Number>
581 void
582 get_function_values(const ReadVector<Number> &fe_function,
584 std::vector<Vector<Number>> &values) const;
585
586
608 template <typename Number>
609 void
610 get_function_values(const ReadVector<Number> &fe_function,
612 ArrayView<std::vector<Number>> values,
613 const bool quadrature_points_fastest) const;
614
617
661 template <typename Number>
662 void
664 const ReadVector<Number> &fe_function,
665 std::vector<Tensor<1, spacedim, Number>> &gradients) const;
666
683 template <typename Number>
684 void
686 const ReadVector<Number> &fe_function,
687 std::vector<std::vector<Tensor<1, spacedim, Number>>> &gradients) const;
688
697 template <typename Number>
698 void
700 const ReadVector<Number> &fe_function,
702 std::vector<Tensor<1, spacedim, Number>> &gradients) const;
703
712 template <typename Number>
713 void
715 const ReadVector<Number> &fe_function,
717 ArrayView<std::vector<Tensor<1, spacedim, Number>>> gradients,
718 const bool quadrature_points_fastest = false) const;
719
724
763 template <typename Number>
764 void
766 const ReadVector<Number> &fe_function,
767 std::vector<Tensor<2, spacedim, Number>> &hessians) const;
768
786 template <typename Number>
787 void
789 const ReadVector<Number> &fe_function,
790 std::vector<std::vector<Tensor<2, spacedim, Number>>> &hessians,
791 const bool quadrature_points_fastest = false) const;
792
801 template <typename Number>
802 void
804 const ReadVector<Number> &fe_function,
806 std::vector<Tensor<2, spacedim, Number>> &hessians) const;
807
816 template <typename Number>
817 void
819 const ReadVector<Number> &fe_function,
821 ArrayView<std::vector<Tensor<2, spacedim, Number>>> hessians,
822 const bool quadrature_points_fastest = false) const;
823
864 template <typename Number>
865 void
867 std::vector<Number> &laplacians) const;
868
888 template <typename Number>
889 void
891 std::vector<Vector<Number>> &laplacians) const;
892
901 template <typename Number>
902 void
904 const ReadVector<Number> &fe_function,
906 std::vector<Number> &laplacians) const;
907
916 template <typename Number>
917 void
919 const ReadVector<Number> &fe_function,
921 std::vector<Vector<Number>> &laplacians) const;
922
931 template <typename Number>
932 void
934 const ReadVector<Number> &fe_function,
936 std::vector<std::vector<Number>> &laplacians,
937 const bool quadrature_points_fastest = false) const;
938
941
981 template <typename Number>
982 void
984 const ReadVector<Number> &fe_function,
985 std::vector<Tensor<3, spacedim, Number>> &third_derivatives) const;
986
1006 template <typename Number>
1007 void
1009 const ReadVector<Number> &fe_function,
1010 std::vector<std::vector<Tensor<3, spacedim, Number>>> &third_derivatives,
1011 const bool quadrature_points_fastest = false) const;
1012
1021 template <typename Number>
1022 void
1024 const ReadVector<Number> &fe_function,
1026 std::vector<Tensor<3, spacedim, Number>> &third_derivatives) const;
1027
1036 template <typename Number>
1037 void
1039 const ReadVector<Number> &fe_function,
1041 ArrayView<std::vector<Tensor<3, spacedim, Number>>> third_derivatives,
1042 const bool quadrature_points_fastest = false) const;
1046
1073
1107 dof_indices_starting_at(const unsigned int start_dof_index) const;
1108
1140 dof_indices_ending_at(const unsigned int end_dof_index) const;
1141
1145
1169
1176 const Point<spacedim> &
1177 quadrature_point(const unsigned int q_point) const;
1178
1184 const std::vector<Point<spacedim>> &
1186
1203 double
1204 JxW(const unsigned int q_point) const;
1205
1209 const std::vector<double> &
1211
1219 jacobian(const unsigned int q_point) const;
1220
1227 const std::vector<DerivativeForm<1, dim, spacedim>> &
1229
1238 jacobian_grad(const unsigned int q_point) const;
1239
1246 const std::vector<DerivativeForm<2, dim, spacedim>> &
1248
1257 const Tensor<3, spacedim> &
1258 jacobian_pushed_forward_grad(const unsigned int q_point) const;
1259
1266 const std::vector<Tensor<3, spacedim>> &
1268
1277 jacobian_2nd_derivative(const unsigned int q_point) const;
1278
1285 const std::vector<DerivativeForm<3, dim, spacedim>> &
1287
1297 const Tensor<4, spacedim> &
1298 jacobian_pushed_forward_2nd_derivative(const unsigned int q_point) const;
1299
1306 const std::vector<Tensor<4, spacedim>> &
1308
1318 jacobian_3rd_derivative(const unsigned int q_point) const;
1319
1326 const std::vector<DerivativeForm<4, dim, spacedim>> &
1328
1338 const Tensor<5, spacedim> &
1339 jacobian_pushed_forward_3rd_derivative(const unsigned int q_point) const;
1340
1347 const std::vector<Tensor<5, spacedim>> &
1349
1357 inverse_jacobian(const unsigned int q_point) const;
1358
1365 const std::vector<DerivativeForm<1, spacedim, dim>> &
1367
1387 const Tensor<1, spacedim> &
1388 normal_vector(const unsigned int q_point) const;
1389
1397 const std::vector<Tensor<1, spacedim>> &
1398 get_normal_vectors() const;
1399
1403
1415
1426
1438
1439
1450
1454
1461
1466 get_fe() const;
1467
1473
1478 get_cell() const;
1479
1486 get_cell_similarity() const;
1487
1492 std::size_t
1493 memory_consumption() const;
1505 std::string,
1506 << "You are requesting information from an FEValues/FEFaceValues/FESubfaceValues "
1507 << "object for which this kind of information has not been computed. What "
1508 << "information these objects compute is determined by the update_* flags you "
1509 << "pass to the constructor. Here, the operation you are attempting requires "
1510 << "the <" << arg1
1511 << "> flag to be set, but it was apparently not specified "
1512 << "upon construction.");
1513
1520 "FEValues object is not reinit'ed to any cell");
1521
1530 "The FiniteElement you provided to FEValues and the FiniteElement that belongs "
1531 "to the DoFHandler that provided the cell iterator do not match.");
1538 int,
1539 << "The shape function with index " << arg1
1540 << " is not primitive, i.e. it is vector-valued and "
1541 << "has more than one non-zero vector component. This "
1542 << "function cannot be called for these shape functions. "
1543 << "Maybe you want to use the same function with the "
1544 << "_component suffix?");
1545
1554 "The given FiniteElement is not a primitive element but the requested operation "
1555 "only works for those. See FiniteElement::is_primitive() for more information.");
1556
1557protected:
1575 {
1576 public:
1579 "You have previously called the FEValues::reinit() function with a "
1580 "cell iterator of type Triangulation<dim,spacedim>::cell_iterator. However, "
1581 "when you do this, you cannot call some functions in the FEValues "
1582 "class, such as the get_function_values/gradients/hessians/third_derivatives "
1583 "functions. If you need these functions, then you need to call "
1584 "FEValues::reinit() with an iterator type that allows to extract "
1585 "degrees of freedom, such as DoFHandler<dim,spacedim>::cell_iterator.");
1586
1592
1598
1604
1610
1614 bool
1615 is_initialized() const;
1616
1623 operator typename Triangulation<dim, spacedim>::cell_iterator() const;
1624
1630 n_dofs_for_dof_handler() const;
1631
1636 template <typename Number>
1637 void
1639 Vector<Number> &out) const;
1640
1641 private:
1647 std::optional<
1648 std::variant<typename Triangulation<dim, spacedim>::cell_iterator,
1652 };
1653
1660
1668 boost::signals2::connection tria_listener_refinement;
1669
1677 boost::signals2::connection tria_listener_mesh_transform;
1678
1684 void
1686
1696 void
1698 const typename Triangulation<dim, spacedim>::cell_iterator &cell);
1699
1706
1712 std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
1714
1721
1729
1735 std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
1737
1743 spacedim>
1745
1746
1751
1762
1769
1775 void
1777 const typename Triangulation<dim, spacedim>::cell_iterator &cell);
1778
1779private:
1784
1789
1790 // Make the view classes friends of this class, since they access internal
1791 // data.
1792 template <int, int>
1794 template <int, int>
1796 template <int, int, int>
1798 template <int, int, int>
1800};
1801
1802#ifndef DOXYGEN
1803
1804/*---------------------- Inline functions: FEValuesBase ---------------------*/
1805
1806template <int dim, int spacedim>
1809 const FEValuesExtractors::Scalar &scalar) const
1810{
1811 AssertIndexRange(scalar.component, fe_values_views_cache.scalars.size());
1812
1813 return fe_values_views_cache.scalars[scalar.component].value_or_initialize(
1814 [scalar, this]() {
1815 return FEValuesViews::Scalar<dim, spacedim>(*this, scalar.component);
1816 });
1817}
1818
1819
1820
1821template <int dim, int spacedim>
1824 const FEValuesExtractors::Vector &vector) const
1825{
1827 fe_values_views_cache.vectors.size());
1828
1829 return fe_values_views_cache.vectors[vector.first_vector_component]
1830 .value_or_initialize([vector, this]() {
1832 *this, vector.first_vector_component);
1833 });
1834}
1835
1836
1837
1838template <int dim, int spacedim>
1841 const FEValuesExtractors::SymmetricTensor<2> &tensor) const
1842{
1843 Assert(
1844 tensor.first_tensor_component <
1845 fe_values_views_cache.symmetric_second_order_tensors.size(),
1847 0,
1848 fe_values_views_cache.symmetric_second_order_tensors.size()));
1849
1850 return fe_values_views_cache
1852 .value_or_initialize([tensor, this]() {
1854 *this, tensor.first_tensor_component);
1855 });
1856}
1857
1858
1859
1860template <int dim, int spacedim>
1863 const FEValuesExtractors::Tensor<2> &tensor) const
1864{
1866 fe_values_views_cache.second_order_tensors.size());
1867
1868 return fe_values_views_cache
1870 .value_or_initialize([tensor, this]() {
1872 *this, tensor.first_tensor_component);
1873 });
1874}
1875
1876
1877
1878template <int dim, int spacedim>
1879inline const double &
1880FEValuesBase<dim, spacedim>::shape_value(const unsigned int i,
1881 const unsigned int q_point) const
1882{
1883 AssertIndexRange(i, fe->n_dofs_per_cell());
1884 Assert(this->update_flags & update_values,
1885 ExcAccessToUninitializedField("update_values"));
1886 Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
1887 Assert(present_cell.is_initialized(), ExcNotReinited());
1888 // if the entire FE is primitive,
1889 // then we can take a short-cut:
1890 if (fe->is_primitive())
1891 return this->finite_element_output.shape_values(i, q_point);
1892 else
1893 {
1894 // otherwise, use the mapping
1895 // between shape function
1896 // numbers and rows. note that
1897 // by the assertions above, we
1898 // know that this particular
1899 // shape function is primitive,
1900 // so we can call
1901 // system_to_component_index
1902 const unsigned int row =
1903 this->finite_element_output
1904 .shape_function_to_row_table[i * fe->n_components() +
1905 fe->system_to_component_index(i).first];
1906 return this->finite_element_output.shape_values(row, q_point);
1907 }
1908}
1909
1910
1911
1912template <int dim, int spacedim>
1913inline double
1915 const unsigned int i,
1916 const unsigned int q_point,
1917 const unsigned int component) const
1918{
1919 AssertIndexRange(i, fe->n_dofs_per_cell());
1920 Assert(this->update_flags & update_values,
1921 ExcAccessToUninitializedField("update_values"));
1922 AssertIndexRange(component, fe->n_components());
1923 Assert(present_cell.is_initialized(), ExcNotReinited());
1924
1925 // check whether the shape function
1926 // is non-zero at all within
1927 // this component:
1928 if (fe->get_nonzero_components(i)[component] == false)
1929 return 0;
1930
1931 // look up the right row in the
1932 // table and take the data from
1933 // there
1934 const unsigned int row =
1935 this->finite_element_output
1936 .shape_function_to_row_table[i * fe->n_components() + component];
1937 return this->finite_element_output.shape_values(row, q_point);
1938}
1939
1940
1941
1942template <int dim, int spacedim>
1943inline const Tensor<1, spacedim> &
1944FEValuesBase<dim, spacedim>::shape_grad(const unsigned int i,
1945 const unsigned int q_point) const
1946{
1947 AssertIndexRange(i, fe->n_dofs_per_cell());
1948 Assert(this->update_flags & update_gradients,
1949 ExcAccessToUninitializedField("update_gradients"));
1950 Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
1951 Assert(present_cell.is_initialized(), ExcNotReinited());
1952 // if the entire FE is primitive,
1953 // then we can take a short-cut:
1954 if (fe->is_primitive())
1955 return this->finite_element_output.shape_gradients[i][q_point];
1956 else
1957 {
1958 // otherwise, use the mapping
1959 // between shape function
1960 // numbers and rows. note that
1961 // by the assertions above, we
1962 // know that this particular
1963 // shape function is primitive,
1964 // so we can call
1965 // system_to_component_index
1966 const unsigned int row =
1967 this->finite_element_output
1968 .shape_function_to_row_table[i * fe->n_components() +
1969 fe->system_to_component_index(i).first];
1970 return this->finite_element_output.shape_gradients[row][q_point];
1971 }
1972}
1973
1974
1975
1976template <int dim, int spacedim>
1979 const unsigned int i,
1980 const unsigned int q_point,
1981 const unsigned int component) const
1982{
1983 AssertIndexRange(i, fe->n_dofs_per_cell());
1984 Assert(this->update_flags & update_gradients,
1985 ExcAccessToUninitializedField("update_gradients"));
1986 AssertIndexRange(component, fe->n_components());
1987 Assert(present_cell.is_initialized(), ExcNotReinited());
1988 // check whether the shape function
1989 // is non-zero at all within
1990 // this component:
1991 if (fe->get_nonzero_components(i)[component] == false)
1992 return Tensor<1, spacedim>();
1993
1994 // look up the right row in the
1995 // table and take the data from
1996 // there
1997 const unsigned int row =
1998 this->finite_element_output
1999 .shape_function_to_row_table[i * fe->n_components() + component];
2000 return this->finite_element_output.shape_gradients[row][q_point];
2001}
2002
2003
2004
2005template <int dim, int spacedim>
2006inline const Tensor<2, spacedim> &
2008 const unsigned int q_point) const
2009{
2010 AssertIndexRange(i, fe->n_dofs_per_cell());
2011 Assert(this->update_flags & update_hessians,
2012 ExcAccessToUninitializedField("update_hessians"));
2013 Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
2014 Assert(present_cell.is_initialized(), ExcNotReinited());
2015 // if the entire FE is primitive,
2016 // then we can take a short-cut:
2017 if (fe->is_primitive())
2018 return this->finite_element_output.shape_hessians[i][q_point];
2019 else
2020 {
2021 // otherwise, use the mapping
2022 // between shape function
2023 // numbers and rows. note that
2024 // by the assertions above, we
2025 // know that this particular
2026 // shape function is primitive,
2027 // so we can call
2028 // system_to_component_index
2029 const unsigned int row =
2030 this->finite_element_output
2031 .shape_function_to_row_table[i * fe->n_components() +
2032 fe->system_to_component_index(i).first];
2033 return this->finite_element_output.shape_hessians[row][q_point];
2034 }
2035}
2036
2037
2038
2039template <int dim, int spacedim>
2042 const unsigned int i,
2043 const unsigned int q_point,
2044 const unsigned int component) const
2045{
2046 AssertIndexRange(i, fe->n_dofs_per_cell());
2047 Assert(this->update_flags & update_hessians,
2048 ExcAccessToUninitializedField("update_hessians"));
2049 AssertIndexRange(component, fe->n_components());
2050 Assert(present_cell.is_initialized(), ExcNotReinited());
2051 // check whether the shape function
2052 // is non-zero at all within
2053 // this component:
2054 if (fe->get_nonzero_components(i)[component] == false)
2055 return Tensor<2, spacedim>();
2056
2057 // look up the right row in the
2058 // table and take the data from
2059 // there
2060 const unsigned int row =
2061 this->finite_element_output
2062 .shape_function_to_row_table[i * fe->n_components() + component];
2063 return this->finite_element_output.shape_hessians[row][q_point];
2064}
2065
2066
2067
2068template <int dim, int spacedim>
2069inline const Tensor<3, spacedim> &
2071 const unsigned int i,
2072 const unsigned int q_point) const
2073{
2074 AssertIndexRange(i, fe->n_dofs_per_cell());
2075 Assert(this->update_flags & update_3rd_derivatives,
2076 ExcAccessToUninitializedField("update_3rd_derivatives"));
2077 Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
2078 Assert(present_cell.is_initialized(), ExcNotReinited());
2079 // if the entire FE is primitive,
2080 // then we can take a short-cut:
2081 if (fe->is_primitive())
2082 return this->finite_element_output.shape_3rd_derivatives[i][q_point];
2083 else
2084 {
2085 // otherwise, use the mapping
2086 // between shape function
2087 // numbers and rows. note that
2088 // by the assertions above, we
2089 // know that this particular
2090 // shape function is primitive,
2091 // so we can call
2092 // system_to_component_index
2093 const unsigned int row =
2094 this->finite_element_output
2095 .shape_function_to_row_table[i * fe->n_components() +
2096 fe->system_to_component_index(i).first];
2097 return this->finite_element_output.shape_3rd_derivatives[row][q_point];
2098 }
2099}
2100
2101
2102
2103template <int dim, int spacedim>
2106 const unsigned int i,
2107 const unsigned int q_point,
2108 const unsigned int component) const
2109{
2110 AssertIndexRange(i, fe->n_dofs_per_cell());
2111 Assert(this->update_flags & update_3rd_derivatives,
2112 ExcAccessToUninitializedField("update_3rd_derivatives"));
2113 AssertIndexRange(component, fe->n_components());
2114 Assert(present_cell.is_initialized(), ExcNotReinited());
2115 // check whether the shape function
2116 // is non-zero at all within
2117 // this component:
2118 if (fe->get_nonzero_components(i)[component] == false)
2119 return Tensor<3, spacedim>();
2120
2121 // look up the right row in the
2122 // table and take the data from
2123 // there
2124 const unsigned int row =
2125 this->finite_element_output
2126 .shape_function_to_row_table[i * fe->n_components() + component];
2127 return this->finite_element_output.shape_3rd_derivatives[row][q_point];
2128}
2129
2130
2131
2132template <int dim, int spacedim>
2133inline const FiniteElement<dim, spacedim> &
2135{
2136 return *fe;
2137}
2138
2139
2140
2141template <int dim, int spacedim>
2142inline const Mapping<dim, spacedim> &
2144{
2145 return *mapping;
2146}
2147
2148
2149
2150template <int dim, int spacedim>
2151inline UpdateFlags
2153{
2154 return this->update_flags;
2155}
2156
2157
2158
2159template <int dim, int spacedim>
2160inline const std::vector<Point<spacedim>> &
2162{
2163 Assert(this->update_flags & update_quadrature_points,
2164 ExcAccessToUninitializedField("update_quadrature_points"));
2165 Assert(present_cell.is_initialized(), ExcNotReinited());
2166 return this->mapping_output.quadrature_points;
2167}
2168
2169
2170
2171template <int dim, int spacedim>
2172inline const std::vector<double> &
2174{
2175 Assert(this->update_flags & update_JxW_values,
2176 ExcAccessToUninitializedField("update_JxW_values"));
2177 Assert(present_cell.is_initialized(), ExcNotReinited());
2178 return this->mapping_output.JxW_values;
2179}
2180
2181
2182
2183template <int dim, int spacedim>
2184inline const std::vector<DerivativeForm<1, dim, spacedim>> &
2186{
2187 Assert(this->update_flags & update_jacobians,
2188 ExcAccessToUninitializedField("update_jacobians"));
2189 Assert(present_cell.is_initialized(), ExcNotReinited());
2190 return this->mapping_output.jacobians;
2191}
2192
2193
2194
2195template <int dim, int spacedim>
2196inline const std::vector<DerivativeForm<2, dim, spacedim>> &
2198{
2199 Assert(this->update_flags & update_jacobian_grads,
2200 ExcAccessToUninitializedField("update_jacobians_grads"));
2201 Assert(present_cell.is_initialized(), ExcNotReinited());
2202 return this->mapping_output.jacobian_grads;
2203}
2204
2205
2206
2207template <int dim, int spacedim>
2208inline const Tensor<3, spacedim> &
2210 const unsigned int q_point) const
2211{
2212 Assert(this->update_flags & update_jacobian_pushed_forward_grads,
2213 ExcAccessToUninitializedField("update_jacobian_pushed_forward_grads"));
2214 Assert(present_cell.is_initialized(), ExcNotReinited());
2215 return this->mapping_output.jacobian_pushed_forward_grads[q_point];
2216}
2217
2218
2219
2220template <int dim, int spacedim>
2221inline const std::vector<Tensor<3, spacedim>> &
2223{
2224 Assert(this->update_flags & update_jacobian_pushed_forward_grads,
2225 ExcAccessToUninitializedField("update_jacobian_pushed_forward_grads"));
2226 Assert(present_cell.is_initialized(), ExcNotReinited());
2227 return this->mapping_output.jacobian_pushed_forward_grads;
2228}
2229
2230
2231
2232template <int dim, int spacedim>
2235 const unsigned int q_point) const
2236{
2237 Assert(this->update_flags & update_jacobian_2nd_derivatives,
2238 ExcAccessToUninitializedField("update_jacobian_2nd_derivatives"));
2239 Assert(present_cell.is_initialized(), ExcNotReinited());
2240 return this->mapping_output.jacobian_2nd_derivatives[q_point];
2241}
2242
2243
2244
2245template <int dim, int spacedim>
2246inline const std::vector<DerivativeForm<3, dim, spacedim>> &
2248{
2249 Assert(this->update_flags & update_jacobian_2nd_derivatives,
2250 ExcAccessToUninitializedField("update_jacobian_2nd_derivatives"));
2251 Assert(present_cell.is_initialized(), ExcNotReinited());
2252 return this->mapping_output.jacobian_2nd_derivatives;
2253}
2254
2255
2256
2257template <int dim, int spacedim>
2258inline const Tensor<4, spacedim> &
2260 const unsigned int q_point) const
2261{
2263 ExcAccessToUninitializedField(
2264 "update_jacobian_pushed_forward_2nd_derivatives"));
2265 Assert(present_cell.is_initialized(), ExcNotReinited());
2266 return this->mapping_output.jacobian_pushed_forward_2nd_derivatives[q_point];
2267}
2268
2269
2270
2271template <int dim, int spacedim>
2272inline const std::vector<Tensor<4, spacedim>> &
2274{
2276 ExcAccessToUninitializedField(
2277 "update_jacobian_pushed_forward_2nd_derivatives"));
2278 Assert(present_cell.is_initialized(), ExcNotReinited());
2279 return this->mapping_output.jacobian_pushed_forward_2nd_derivatives;
2280}
2281
2282
2283
2284template <int dim, int spacedim>
2287 const unsigned int q_point) const
2288{
2289 Assert(this->update_flags & update_jacobian_3rd_derivatives,
2290 ExcAccessToUninitializedField("update_jacobian_3rd_derivatives"));
2291 Assert(present_cell.is_initialized(), ExcNotReinited());
2292 return this->mapping_output.jacobian_3rd_derivatives[q_point];
2293}
2294
2295
2296
2297template <int dim, int spacedim>
2298inline const std::vector<DerivativeForm<4, dim, spacedim>> &
2300{
2301 Assert(this->update_flags & update_jacobian_3rd_derivatives,
2302 ExcAccessToUninitializedField("update_jacobian_3rd_derivatives"));
2303 Assert(present_cell.is_initialized(), ExcNotReinited());
2304 return this->mapping_output.jacobian_3rd_derivatives;
2305}
2306
2307
2308
2309template <int dim, int spacedim>
2310inline const Tensor<5, spacedim> &
2312 const unsigned int q_point) const
2313{
2315 ExcAccessToUninitializedField(
2316 "update_jacobian_pushed_forward_3rd_derivatives"));
2317 Assert(present_cell.is_initialized(), ExcNotReinited());
2318 return this->mapping_output.jacobian_pushed_forward_3rd_derivatives[q_point];
2319}
2320
2321
2322
2323template <int dim, int spacedim>
2324inline const std::vector<Tensor<5, spacedim>> &
2326{
2328 ExcAccessToUninitializedField(
2329 "update_jacobian_pushed_forward_3rd_derivatives"));
2330 Assert(present_cell.is_initialized(), ExcNotReinited());
2331 return this->mapping_output.jacobian_pushed_forward_3rd_derivatives;
2332}
2333
2334
2335
2336template <int dim, int spacedim>
2337inline const std::vector<DerivativeForm<1, spacedim, dim>> &
2339{
2340 Assert(this->update_flags & update_inverse_jacobians,
2341 ExcAccessToUninitializedField("update_inverse_jacobians"));
2342 Assert(present_cell.is_initialized(), ExcNotReinited());
2343 return this->mapping_output.inverse_jacobians;
2344}
2345
2346
2347
2348template <int dim, int spacedim>
2351{
2353 0U, dofs_per_cell);
2354}
2355
2356
2357
2358template <int dim, int spacedim>
2361 const unsigned int start_dof_index) const
2362{
2363 Assert(start_dof_index <= dofs_per_cell,
2364 ExcIndexRange(start_dof_index, 0, dofs_per_cell + 1));
2366 start_dof_index, dofs_per_cell);
2367}
2368
2369
2370
2371template <int dim, int spacedim>
2374 const unsigned int end_dof_index) const
2375{
2376 Assert(end_dof_index < dofs_per_cell,
2377 ExcIndexRange(end_dof_index, 0, dofs_per_cell));
2379 0U, end_dof_index + 1);
2380}
2381
2382
2383
2384template <int dim, int spacedim>
2387{
2389 0U, n_quadrature_points);
2390}
2391
2392
2393
2394template <int dim, int spacedim>
2395inline const Point<spacedim> &
2396FEValuesBase<dim, spacedim>::quadrature_point(const unsigned int q_point) const
2397{
2398 Assert(this->update_flags & update_quadrature_points,
2399 ExcAccessToUninitializedField("update_quadrature_points"));
2400 AssertIndexRange(q_point, this->mapping_output.quadrature_points.size());
2401 Assert(present_cell.is_initialized(), ExcNotReinited());
2402
2403 return this->mapping_output.quadrature_points[q_point];
2404}
2405
2406
2407
2408template <int dim, int spacedim>
2409inline double
2410FEValuesBase<dim, spacedim>::JxW(const unsigned int q_point) const
2411{
2412 Assert(this->update_flags & update_JxW_values,
2413 ExcAccessToUninitializedField("update_JxW_values"));
2414 AssertIndexRange(q_point, this->mapping_output.JxW_values.size());
2415 Assert(present_cell.is_initialized(), ExcNotReinited());
2416
2417 return this->mapping_output.JxW_values[q_point];
2418}
2419
2420
2421
2422template <int dim, int spacedim>
2424FEValuesBase<dim, spacedim>::jacobian(const unsigned int q_point) const
2425{
2426 Assert(this->update_flags & update_jacobians,
2427 ExcAccessToUninitializedField("update_jacobians"));
2428 AssertIndexRange(q_point, this->mapping_output.jacobians.size());
2429 Assert(present_cell.is_initialized(), ExcNotReinited());
2430
2431 return this->mapping_output.jacobians[q_point];
2432}
2433
2434
2435
2436template <int dim, int spacedim>
2438FEValuesBase<dim, spacedim>::jacobian_grad(const unsigned int q_point) const
2439{
2440 Assert(this->update_flags & update_jacobian_grads,
2441 ExcAccessToUninitializedField("update_jacobians_grads"));
2442 AssertIndexRange(q_point, this->mapping_output.jacobian_grads.size());
2443 Assert(present_cell.is_initialized(), ExcNotReinited());
2444
2445 return this->mapping_output.jacobian_grads[q_point];
2446}
2447
2448
2449
2450template <int dim, int spacedim>
2452FEValuesBase<dim, spacedim>::inverse_jacobian(const unsigned int q_point) const
2453{
2454 Assert(this->update_flags & update_inverse_jacobians,
2455 ExcAccessToUninitializedField("update_inverse_jacobians"));
2456 AssertIndexRange(q_point, this->mapping_output.inverse_jacobians.size());
2457 Assert(present_cell.is_initialized(), ExcNotReinited());
2458
2459 return this->mapping_output.inverse_jacobians[q_point];
2460}
2461
2462
2463
2464template <int dim, int spacedim>
2465inline const Tensor<1, spacedim> &
2466FEValuesBase<dim, spacedim>::normal_vector(const unsigned int q_point) const
2467{
2468 Assert(this->update_flags & update_normal_vectors,
2470 "update_normal_vectors")));
2471 AssertIndexRange(q_point, this->mapping_output.normal_vectors.size());
2472 Assert(present_cell.is_initialized(), ExcNotReinited());
2473
2474 return this->mapping_output.normal_vectors[q_point];
2475}
2476
2477#endif
2478
2480
2481#endif
typename LevelSelector::cell_iterator level_cell_iterator
types::global_dof_index n_dofs_for_dof_handler() const
std::optional< std::variant< typename Triangulation< dim, spacedim >::cell_iterator, typename DoFHandler< dim, spacedim >::cell_iterator, typename DoFHandler< dim, spacedim >::level_cell_iterator > > cell
void get_interpolated_dof_values(const ReadVector< Number > &in, Vector< Number > &out) const
CellSimilarity::Similarity cell_similarity
const Tensor< 5, spacedim > & jacobian_pushed_forward_3rd_derivative(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
const std::vector< double > & get_JxW_values() const
const DerivativeForm< 1, dim, spacedim > & jacobian(const unsigned int q_point) const
::internal::FEValuesViews::Cache< dim, spacedim > fe_values_views_cache
FEValuesBase(const FEValuesBase &)=delete
const FEValuesViews::Scalar< dim, spacedim > & operator[](const FEValuesExtractors::Scalar &scalar) const
Triangulation< dim, spacedim >::cell_iterator get_cell() const
boost::signals2::connection tria_listener_mesh_transform
void get_function_values(const ReadVector< Number > &fe_function, std::vector< Number > &values) const
const std::vector< Point< spacedim > > & get_quadrature_points() const
const std::vector< DerivativeForm< 1, dim, spacedim > > & get_jacobians() const
const DerivativeForm< 2, dim, spacedim > & jacobian_grad(const unsigned int q_point) const
internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > mapping_output
virtual ~FEValuesBase() override
const std::vector< DerivativeForm< 3, dim, spacedim > > & get_jacobian_2nd_derivatives() const
UpdateFlags get_update_flags() const
const FEValuesViews::Tensor< 2, dim, spacedim > & operator[](const FEValuesExtractors::Tensor< 2 > &tensor) 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
CellIteratorWrapper present_cell
const Tensor< 1, spacedim > & normal_vector(const unsigned int q_point) const
UpdateFlags update_flags
static constexpr unsigned int space_dimension
Tensor< 2, spacedim > shape_hessian_component(const unsigned int i, const unsigned int q_point, const unsigned int component) const
void get_function_hessians(const ReadVector< Number > &fe_function, std::vector< Tensor< 2, spacedim, Number > > &hessians) 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 Point< spacedim > & quadrature_point(const unsigned int q_point) const
std_cxx20::ranges::iota_view< unsigned int, unsigned int > dof_indices_starting_at(const unsigned int start_dof_index) const
double shape_value_component(const unsigned int i, const unsigned int q_point, const unsigned int component) const
const Mapping< dim, spacedim > & get_mapping() const
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< 3, spacedim > & jacobian_pushed_forward_grad(const unsigned int q_point) const
std_cxx20::ranges::iota_view< unsigned int, unsigned int > dof_indices() const
std_cxx20::ranges::iota_view< unsigned int, unsigned int > quadrature_point_indices() const
void get_function_third_derivatives(const ReadVector< Number > &fe_function, std::vector< Tensor< 3, spacedim, Number > > &third_derivatives) const
const Tensor< 2, spacedim > & shape_hessian(const unsigned int i, const unsigned int q_point) const
const std::vector< Tensor< 5, spacedim > > & get_jacobian_pushed_forward_3rd_derivatives() const
const Tensor< 4, spacedim > & jacobian_pushed_forward_2nd_derivative(const unsigned int q_point) const
const Tensor< 3, spacedim > & shape_3rd_derivative(const unsigned int i, const unsigned int q_point) 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
const std::vector< DerivativeForm< 2, dim, spacedim > > & get_jacobian_grads() const
boost::signals2::connection tria_listener_refinement
const std::vector< Tensor< 3, spacedim > > & get_jacobian_pushed_forward_grads() const
const DerivativeForm< 4, dim, spacedim > & jacobian_3rd_derivative(const unsigned int q_point) const
const DerivativeForm< 1, spacedim, dim > & inverse_jacobian(const unsigned int q_point) const
void invalidate_present_cell()
const std::vector< Tensor< 4, spacedim > > & get_jacobian_pushed_forward_2nd_derivatives() const
const DerivativeForm< 3, dim, spacedim > & jacobian_2nd_derivative(const unsigned int q_point) const
const ObserverPointer< const Mapping< dim, spacedim >, FEValuesBase< dim, spacedim > > mapping
Tensor< 1, spacedim > shape_grad_component(const unsigned int i, const unsigned int q_point, const unsigned int component) const
bool check_for_cell_similarity_allowed
::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > finite_element_output
const FEValuesViews::SymmetricTensor< 2, dim, spacedim > & operator[](const FEValuesExtractors::SymmetricTensor< 2 > &tensor) const
const ObserverPointer< const FiniteElement< dim, spacedim >, FEValuesBase< dim, spacedim > > fe
const Tensor< 1, spacedim > & shape_grad(const unsigned int i, const unsigned int q_point) const
const std::vector< DerivativeForm< 4, dim, spacedim > > & get_jacobian_3rd_derivatives() const
std::unique_ptr< typename FiniteElement< dim, spacedim >::InternalDataBase > fe_data
const FiniteElement< dim, spacedim > & get_fe() const
double JxW(const unsigned int q_point) const
const double & shape_value(const unsigned int i, const unsigned int q_point) const
void maybe_invalidate_previous_present_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell)
FEValuesBase & operator=(const FEValuesBase &)=delete
const std::vector< DerivativeForm< 1, spacedim, dim > > & get_inverse_jacobians() const
const unsigned int max_n_quadrature_points
static constexpr unsigned int dimension
Abstract base class for mapping classes.
Definition mapping.h:320
Definition point.h:113
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:40
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:41
static ::ExceptionBase & ExcAccessToUninitializedField(std::string arg1)
#define Assert(cond, exc)
static ::ExceptionBase & ExcNotReinited()
static ::ExceptionBase & ExcNeedsDoFHandler()
static ::ExceptionBase & ExcShapeFunctionNotPrimitive(int arg1)
#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 & ExcFEDontMatch()
static ::ExceptionBase & ExcFENotPrimitive()
typename ActiveSelector::cell_iterator cell_iterator
TriaIterator< CellAccessor< dim, spacedim > > cell_iterator
Definition tria.h:1557
UpdateFlags
@ 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
Definition iota_view.h:45
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