Reference documentation for deal.II version GIT relicensing-487-ge9eb5ab491 2024-04-25 07:20:02+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\}}\)
Loading...
Searching...
No Matches
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
23#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 <algorithm>
47#include <memory>
48#include <optional>
49#include <type_traits>
50#include <variant>
51
53
153template <int dim, int spacedim>
155{
156public:
160 static constexpr unsigned int dimension = dim;
161
165 static constexpr unsigned int space_dimension = spacedim;
166
174 const unsigned int n_quadrature_points;
175
185 const unsigned int max_n_quadrature_points;
186
192 const unsigned int dofs_per_cell;
193
194
202 FEValuesBase(const unsigned int n_q_points,
203 const unsigned int dofs_per_cell,
207
213 operator=(const FEValuesBase &) = delete;
214
219 FEValuesBase(const FEValuesBase &) = delete;
220
224 virtual ~FEValuesBase() override;
225
246 void
248
252
274 const double &
275 shape_value(const unsigned int i, const unsigned int q_point) const;
276
297 double
298 shape_value_component(const unsigned int i,
299 const unsigned int q_point,
300 const unsigned int component) const;
301
327 const Tensor<1, spacedim> &
328 shape_grad(const unsigned int i, const unsigned int q_point) const;
329
347 shape_grad_component(const unsigned int i,
348 const unsigned int q_point,
349 const unsigned int component) const;
350
370 const Tensor<2, spacedim> &
371 shape_hessian(const unsigned int i, const unsigned int q_point) const;
372
390 shape_hessian_component(const unsigned int i,
391 const unsigned int q_point,
392 const unsigned int component) const;
393
413 const Tensor<3, spacedim> &
414 shape_3rd_derivative(const unsigned int i, const unsigned int q_point) const;
415
433 shape_3rd_derivative_component(const unsigned int i,
434 const unsigned int q_point,
435 const unsigned int component) const;
436
439
485 template <typename Number>
486 void
487 get_function_values(const ReadVector<Number> &fe_function,
488 std::vector<Number> &values) const;
489
503 template <typename Number>
504 void
505 get_function_values(const ReadVector<Number> &fe_function,
506 std::vector<Vector<Number>> &values) const;
507
564 template <typename Number>
565 void
566 get_function_values(const ReadVector<Number> &fe_function,
568 std::vector<Number> &values) const;
569
578 template <typename Number>
579 void
580 get_function_values(const ReadVector<Number> &fe_function,
582 std::vector<Vector<Number>> &values) const;
583
584
606 template <typename Number>
607 void
608 get_function_values(const ReadVector<Number> &fe_function,
610 ArrayView<std::vector<Number>> values,
611 const bool quadrature_points_fastest) const;
612
615
659 template <typename Number>
660 void
662 const ReadVector<Number> &fe_function,
663 std::vector<Tensor<1, spacedim, Number>> &gradients) const;
664
681 template <typename Number>
682 void
684 const ReadVector<Number> &fe_function,
685 std::vector<std::vector<Tensor<1, spacedim, Number>>> &gradients) const;
686
695 template <typename Number>
696 void
698 const ReadVector<Number> &fe_function,
700 std::vector<Tensor<1, spacedim, Number>> &gradients) const;
701
710 template <typename Number>
711 void
713 const ReadVector<Number> &fe_function,
715 ArrayView<std::vector<Tensor<1, spacedim, Number>>> gradients,
716 const bool quadrature_points_fastest = false) const;
717
722
761 template <typename Number>
762 void
764 const ReadVector<Number> &fe_function,
765 std::vector<Tensor<2, spacedim, Number>> &hessians) const;
766
784 template <typename Number>
785 void
787 const ReadVector<Number> &fe_function,
788 std::vector<std::vector<Tensor<2, spacedim, Number>>> &hessians,
789 const bool quadrature_points_fastest = false) const;
790
799 template <typename Number>
800 void
802 const ReadVector<Number> &fe_function,
804 std::vector<Tensor<2, spacedim, Number>> &hessians) const;
805
814 template <typename Number>
815 void
817 const ReadVector<Number> &fe_function,
819 ArrayView<std::vector<Tensor<2, spacedim, Number>>> hessians,
820 const bool quadrature_points_fastest = false) const;
821
862 template <typename Number>
863 void
865 std::vector<Number> &laplacians) const;
866
886 template <typename Number>
887 void
889 std::vector<Vector<Number>> &laplacians) const;
890
899 template <typename Number>
900 void
902 const ReadVector<Number> &fe_function,
904 std::vector<Number> &laplacians) const;
905
914 template <typename Number>
915 void
917 const ReadVector<Number> &fe_function,
919 std::vector<Vector<Number>> &laplacians) const;
920
929 template <typename Number>
930 void
932 const ReadVector<Number> &fe_function,
934 std::vector<std::vector<Number>> &laplacians,
935 const bool quadrature_points_fastest = false) const;
936
939
979 template <typename Number>
980 void
982 const ReadVector<Number> &fe_function,
983 std::vector<Tensor<3, spacedim, Number>> &third_derivatives) const;
984
1004 template <typename Number>
1005 void
1007 const ReadVector<Number> &fe_function,
1008 std::vector<std::vector<Tensor<3, spacedim, Number>>> &third_derivatives,
1009 const bool quadrature_points_fastest = false) const;
1010
1019 template <typename Number>
1020 void
1022 const ReadVector<Number> &fe_function,
1024 std::vector<Tensor<3, spacedim, Number>> &third_derivatives) const;
1025
1034 template <typename Number>
1035 void
1037 const ReadVector<Number> &fe_function,
1039 ArrayView<std::vector<Tensor<3, spacedim, Number>>> third_derivatives,
1040 const bool quadrature_points_fastest = false) const;
1044
1071
1105 dof_indices_starting_at(const unsigned int start_dof_index) const;
1106
1138 dof_indices_ending_at(const unsigned int end_dof_index) const;
1139
1143
1167
1174 const Point<spacedim> &
1175 quadrature_point(const unsigned int q_point) const;
1176
1182 const std::vector<Point<spacedim>> &
1184
1201 double
1202 JxW(const unsigned int q_point) const;
1203
1207 const std::vector<double> &
1209
1217 jacobian(const unsigned int q_point) const;
1218
1225 const std::vector<DerivativeForm<1, dim, spacedim>> &
1227
1236 jacobian_grad(const unsigned int q_point) const;
1237
1244 const std::vector<DerivativeForm<2, dim, spacedim>> &
1246
1255 const Tensor<3, spacedim> &
1256 jacobian_pushed_forward_grad(const unsigned int q_point) const;
1257
1264 const std::vector<Tensor<3, spacedim>> &
1266
1275 jacobian_2nd_derivative(const unsigned int q_point) const;
1276
1283 const std::vector<DerivativeForm<3, dim, spacedim>> &
1285
1295 const Tensor<4, spacedim> &
1296 jacobian_pushed_forward_2nd_derivative(const unsigned int q_point) const;
1297
1304 const std::vector<Tensor<4, spacedim>> &
1306
1316 jacobian_3rd_derivative(const unsigned int q_point) const;
1317
1324 const std::vector<DerivativeForm<4, dim, spacedim>> &
1326
1336 const Tensor<5, spacedim> &
1337 jacobian_pushed_forward_3rd_derivative(const unsigned int q_point) const;
1338
1345 const std::vector<Tensor<5, spacedim>> &
1347
1355 inverse_jacobian(const unsigned int q_point) const;
1356
1363 const std::vector<DerivativeForm<1, spacedim, dim>> &
1365
1385 const Tensor<1, spacedim> &
1386 normal_vector(const unsigned int q_point) const;
1387
1395 const std::vector<Tensor<1, spacedim>> &
1396 get_normal_vectors() const;
1397
1401
1413
1424
1436
1437
1448
1452
1459
1464 get_fe() const;
1465
1471
1476 get_cell() const;
1477
1484 get_cell_similarity() const;
1485
1490 std::size_t
1491 memory_consumption() const;
1503 std::string,
1504 << "You are requesting information from an FEValues/FEFaceValues/FESubfaceValues "
1505 << "object for which this kind of information has not been computed. What "
1506 << "information these objects compute is determined by the update_* flags you "
1507 << "pass to the constructor. Here, the operation you are attempting requires "
1508 << "the <" << arg1
1509 << "> flag to be set, but it was apparently not specified "
1510 << "upon construction.");
1511
1518 "FEValues object is not reinit'ed to any cell");
1519
1528 "The FiniteElement you provided to FEValues and the FiniteElement that belongs "
1529 "to the DoFHandler that provided the cell iterator do not match.");
1536 int,
1537 << "The shape function with index " << arg1
1538 << " is not primitive, i.e. it is vector-valued and "
1539 << "has more than one non-zero vector component. This "
1540 << "function cannot be called for these shape functions. "
1541 << "Maybe you want to use the same function with the "
1542 << "_component suffix?");
1543
1552 "The given FiniteElement is not a primitive element but the requested operation "
1553 "only works for those. See FiniteElement::is_primitive() for more information.");
1554
1555protected:
1573 {
1574 public:
1577 "You have previously called the FEValues::reinit() function with a "
1578 "cell iterator of type Triangulation<dim,spacedim>::cell_iterator. However, "
1579 "when you do this, you cannot call some functions in the FEValues "
1580 "class, such as the get_function_values/gradients/hessians/third_derivatives "
1581 "functions. If you need these functions, then you need to call "
1582 "FEValues::reinit() with an iterator type that allows to extract "
1583 "degrees of freedom, such as DoFHandler<dim,spacedim>::cell_iterator.");
1584
1590
1596
1602
1608
1612 bool
1613 is_initialized() const;
1614
1621 operator typename Triangulation<dim, spacedim>::cell_iterator() const;
1622
1628 n_dofs_for_dof_handler() const;
1629
1634 template <typename Number>
1635 void
1637 Vector<Number> &out) const;
1638
1639 private:
1645 std::optional<
1646 std::variant<typename Triangulation<dim, spacedim>::cell_iterator,
1650 };
1651
1658
1666 boost::signals2::connection tria_listener_refinement;
1667
1675 boost::signals2::connection tria_listener_mesh_transform;
1676
1682 void
1684
1694 void
1696 const typename Triangulation<dim, spacedim>::cell_iterator &cell);
1697
1703
1709 std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
1711
1718
1726
1732 std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
1734
1740 spacedim>
1742
1743
1748
1759
1766
1772 void
1774 const typename Triangulation<dim, spacedim>::cell_iterator &cell);
1775
1776private:
1781
1786
1787 // Make the view classes friends of this class, since they access internal
1788 // data.
1789 template <int, int>
1791 template <int, int>
1793 template <int, int, int>
1795 template <int, int, int>
1797};
1798
1799#ifndef DOXYGEN
1800
1801/*---------------------- Inline functions: FEValuesBase ---------------------*/
1802
1803template <int dim, int spacedim>
1806 const FEValuesExtractors::Scalar &scalar) const
1807{
1808 AssertIndexRange(scalar.component, fe_values_views_cache.scalars.size());
1809
1810 return fe_values_views_cache.scalars[scalar.component].value_or_initialize(
1811 [scalar, this]() {
1812 return FEValuesViews::Scalar<dim, spacedim>(*this, scalar.component);
1813 });
1814}
1815
1816
1817
1818template <int dim, int spacedim>
1821 const FEValuesExtractors::Vector &vector) const
1822{
1824 fe_values_views_cache.vectors.size());
1825
1826 return fe_values_views_cache.vectors[vector.first_vector_component]
1827 .value_or_initialize([vector, this]() {
1829 *this, vector.first_vector_component);
1830 });
1831}
1832
1833
1834
1835template <int dim, int spacedim>
1838 const FEValuesExtractors::SymmetricTensor<2> &tensor) const
1839{
1840 Assert(
1841 tensor.first_tensor_component <
1842 fe_values_views_cache.symmetric_second_order_tensors.size(),
1844 0,
1845 fe_values_views_cache.symmetric_second_order_tensors.size()));
1846
1847 return fe_values_views_cache
1849 .value_or_initialize([tensor, this]() {
1851 *this, tensor.first_tensor_component);
1852 });
1853}
1854
1855
1856
1857template <int dim, int spacedim>
1860 const FEValuesExtractors::Tensor<2> &tensor) const
1861{
1863 fe_values_views_cache.second_order_tensors.size());
1864
1865 return fe_values_views_cache
1867 .value_or_initialize([tensor, this]() {
1869 *this, tensor.first_tensor_component);
1870 });
1871}
1872
1873
1874
1875template <int dim, int spacedim>
1876inline const double &
1877FEValuesBase<dim, spacedim>::shape_value(const unsigned int i,
1878 const unsigned int q_point) const
1879{
1880 AssertIndexRange(i, fe->n_dofs_per_cell());
1881 Assert(this->update_flags & update_values,
1882 ExcAccessToUninitializedField("update_values"));
1883 Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
1884 Assert(present_cell.is_initialized(), ExcNotReinited());
1885 // if the entire FE is primitive,
1886 // then we can take a short-cut:
1887 if (fe->is_primitive())
1888 return this->finite_element_output.shape_values(i, q_point);
1889 else
1890 {
1891 // otherwise, use the mapping
1892 // between shape function
1893 // numbers and rows. note that
1894 // by the assertions above, we
1895 // know that this particular
1896 // shape function is primitive,
1897 // so we can call
1898 // system_to_component_index
1899 const unsigned int row =
1900 this->finite_element_output
1901 .shape_function_to_row_table[i * fe->n_components() +
1902 fe->system_to_component_index(i).first];
1903 return this->finite_element_output.shape_values(row, q_point);
1904 }
1905}
1906
1907
1908
1909template <int dim, int spacedim>
1910inline double
1912 const unsigned int i,
1913 const unsigned int q_point,
1914 const unsigned int component) const
1915{
1916 AssertIndexRange(i, fe->n_dofs_per_cell());
1917 Assert(this->update_flags & update_values,
1918 ExcAccessToUninitializedField("update_values"));
1919 AssertIndexRange(component, fe->n_components());
1920 Assert(present_cell.is_initialized(), ExcNotReinited());
1921
1922 // check whether the shape function
1923 // is non-zero at all within
1924 // this component:
1925 if (fe->get_nonzero_components(i)[component] == false)
1926 return 0;
1927
1928 // look up the right row in the
1929 // table and take the data from
1930 // there
1931 const unsigned int row =
1932 this->finite_element_output
1933 .shape_function_to_row_table[i * fe->n_components() + component];
1934 return this->finite_element_output.shape_values(row, q_point);
1935}
1936
1937
1938
1939template <int dim, int spacedim>
1940inline const Tensor<1, spacedim> &
1941FEValuesBase<dim, spacedim>::shape_grad(const unsigned int i,
1942 const unsigned int q_point) const
1943{
1944 AssertIndexRange(i, fe->n_dofs_per_cell());
1945 Assert(this->update_flags & update_gradients,
1946 ExcAccessToUninitializedField("update_gradients"));
1947 Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
1948 Assert(present_cell.is_initialized(), ExcNotReinited());
1949 // if the entire FE is primitive,
1950 // then we can take a short-cut:
1951 if (fe->is_primitive())
1952 return this->finite_element_output.shape_gradients[i][q_point];
1953 else
1954 {
1955 // otherwise, use the mapping
1956 // between shape function
1957 // numbers and rows. note that
1958 // by the assertions above, we
1959 // know that this particular
1960 // shape function is primitive,
1961 // so we can call
1962 // system_to_component_index
1963 const unsigned int row =
1964 this->finite_element_output
1965 .shape_function_to_row_table[i * fe->n_components() +
1966 fe->system_to_component_index(i).first];
1967 return this->finite_element_output.shape_gradients[row][q_point];
1968 }
1969}
1970
1971
1972
1973template <int dim, int spacedim>
1976 const unsigned int i,
1977 const unsigned int q_point,
1978 const unsigned int component) const
1979{
1980 AssertIndexRange(i, fe->n_dofs_per_cell());
1981 Assert(this->update_flags & update_gradients,
1982 ExcAccessToUninitializedField("update_gradients"));
1983 AssertIndexRange(component, fe->n_components());
1984 Assert(present_cell.is_initialized(), ExcNotReinited());
1985 // check whether the shape function
1986 // is non-zero at all within
1987 // this component:
1988 if (fe->get_nonzero_components(i)[component] == false)
1989 return Tensor<1, spacedim>();
1990
1991 // look up the right row in the
1992 // table and take the data from
1993 // there
1994 const unsigned int row =
1995 this->finite_element_output
1996 .shape_function_to_row_table[i * fe->n_components() + component];
1997 return this->finite_element_output.shape_gradients[row][q_point];
1998}
1999
2000
2001
2002template <int dim, int spacedim>
2003inline const Tensor<2, spacedim> &
2005 const unsigned int q_point) const
2006{
2007 AssertIndexRange(i, fe->n_dofs_per_cell());
2008 Assert(this->update_flags & update_hessians,
2009 ExcAccessToUninitializedField("update_hessians"));
2010 Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
2011 Assert(present_cell.is_initialized(), ExcNotReinited());
2012 // if the entire FE is primitive,
2013 // then we can take a short-cut:
2014 if (fe->is_primitive())
2015 return this->finite_element_output.shape_hessians[i][q_point];
2016 else
2017 {
2018 // otherwise, use the mapping
2019 // between shape function
2020 // numbers and rows. note that
2021 // by the assertions above, we
2022 // know that this particular
2023 // shape function is primitive,
2024 // so we can call
2025 // system_to_component_index
2026 const unsigned int row =
2027 this->finite_element_output
2028 .shape_function_to_row_table[i * fe->n_components() +
2029 fe->system_to_component_index(i).first];
2030 return this->finite_element_output.shape_hessians[row][q_point];
2031 }
2032}
2033
2034
2035
2036template <int dim, int spacedim>
2039 const unsigned int i,
2040 const unsigned int q_point,
2041 const unsigned int component) const
2042{
2043 AssertIndexRange(i, fe->n_dofs_per_cell());
2044 Assert(this->update_flags & update_hessians,
2045 ExcAccessToUninitializedField("update_hessians"));
2046 AssertIndexRange(component, fe->n_components());
2047 Assert(present_cell.is_initialized(), ExcNotReinited());
2048 // check whether the shape function
2049 // is non-zero at all within
2050 // this component:
2051 if (fe->get_nonzero_components(i)[component] == false)
2052 return Tensor<2, spacedim>();
2053
2054 // look up the right row in the
2055 // table and take the data from
2056 // there
2057 const unsigned int row =
2058 this->finite_element_output
2059 .shape_function_to_row_table[i * fe->n_components() + component];
2060 return this->finite_element_output.shape_hessians[row][q_point];
2061}
2062
2063
2064
2065template <int dim, int spacedim>
2066inline const Tensor<3, spacedim> &
2068 const unsigned int i,
2069 const unsigned int q_point) const
2070{
2071 AssertIndexRange(i, fe->n_dofs_per_cell());
2072 Assert(this->update_flags & update_3rd_derivatives,
2073 ExcAccessToUninitializedField("update_3rd_derivatives"));
2074 Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
2075 Assert(present_cell.is_initialized(), ExcNotReinited());
2076 // if the entire FE is primitive,
2077 // then we can take a short-cut:
2078 if (fe->is_primitive())
2079 return this->finite_element_output.shape_3rd_derivatives[i][q_point];
2080 else
2081 {
2082 // otherwise, use the mapping
2083 // between shape function
2084 // numbers and rows. note that
2085 // by the assertions above, we
2086 // know that this particular
2087 // shape function is primitive,
2088 // so we can call
2089 // system_to_component_index
2090 const unsigned int row =
2091 this->finite_element_output
2092 .shape_function_to_row_table[i * fe->n_components() +
2093 fe->system_to_component_index(i).first];
2094 return this->finite_element_output.shape_3rd_derivatives[row][q_point];
2095 }
2096}
2097
2098
2099
2100template <int dim, int spacedim>
2103 const unsigned int i,
2104 const unsigned int q_point,
2105 const unsigned int component) const
2106{
2107 AssertIndexRange(i, fe->n_dofs_per_cell());
2108 Assert(this->update_flags & update_3rd_derivatives,
2109 ExcAccessToUninitializedField("update_3rd_derivatives"));
2110 AssertIndexRange(component, fe->n_components());
2111 Assert(present_cell.is_initialized(), ExcNotReinited());
2112 // check whether the shape function
2113 // is non-zero at all within
2114 // this component:
2115 if (fe->get_nonzero_components(i)[component] == false)
2116 return Tensor<3, spacedim>();
2117
2118 // look up the right row in the
2119 // table and take the data from
2120 // there
2121 const unsigned int row =
2122 this->finite_element_output
2123 .shape_function_to_row_table[i * fe->n_components() + component];
2124 return this->finite_element_output.shape_3rd_derivatives[row][q_point];
2125}
2126
2127
2128
2129template <int dim, int spacedim>
2130inline const FiniteElement<dim, spacedim> &
2132{
2133 return *fe;
2134}
2135
2136
2137
2138template <int dim, int spacedim>
2139inline const Mapping<dim, spacedim> &
2141{
2142 return *mapping;
2143}
2144
2145
2146
2147template <int dim, int spacedim>
2148inline UpdateFlags
2150{
2151 return this->update_flags;
2152}
2153
2154
2155
2156template <int dim, int spacedim>
2157inline const std::vector<Point<spacedim>> &
2159{
2160 Assert(this->update_flags & update_quadrature_points,
2161 ExcAccessToUninitializedField("update_quadrature_points"));
2162 Assert(present_cell.is_initialized(), ExcNotReinited());
2163 return this->mapping_output.quadrature_points;
2164}
2165
2166
2167
2168template <int dim, int spacedim>
2169inline const std::vector<double> &
2171{
2172 Assert(this->update_flags & update_JxW_values,
2173 ExcAccessToUninitializedField("update_JxW_values"));
2174 Assert(present_cell.is_initialized(), ExcNotReinited());
2175 return this->mapping_output.JxW_values;
2176}
2177
2178
2179
2180template <int dim, int spacedim>
2181inline const std::vector<DerivativeForm<1, dim, spacedim>> &
2183{
2184 Assert(this->update_flags & update_jacobians,
2185 ExcAccessToUninitializedField("update_jacobians"));
2186 Assert(present_cell.is_initialized(), ExcNotReinited());
2187 return this->mapping_output.jacobians;
2188}
2189
2190
2191
2192template <int dim, int spacedim>
2193inline const std::vector<DerivativeForm<2, dim, spacedim>> &
2195{
2196 Assert(this->update_flags & update_jacobian_grads,
2197 ExcAccessToUninitializedField("update_jacobians_grads"));
2198 Assert(present_cell.is_initialized(), ExcNotReinited());
2199 return this->mapping_output.jacobian_grads;
2200}
2201
2202
2203
2204template <int dim, int spacedim>
2205inline const Tensor<3, spacedim> &
2207 const unsigned int q_point) const
2208{
2209 Assert(this->update_flags & update_jacobian_pushed_forward_grads,
2210 ExcAccessToUninitializedField("update_jacobian_pushed_forward_grads"));
2211 Assert(present_cell.is_initialized(), ExcNotReinited());
2212 return this->mapping_output.jacobian_pushed_forward_grads[q_point];
2213}
2214
2215
2216
2217template <int dim, int spacedim>
2218inline const std::vector<Tensor<3, spacedim>> &
2220{
2221 Assert(this->update_flags & update_jacobian_pushed_forward_grads,
2222 ExcAccessToUninitializedField("update_jacobian_pushed_forward_grads"));
2223 Assert(present_cell.is_initialized(), ExcNotReinited());
2224 return this->mapping_output.jacobian_pushed_forward_grads;
2225}
2226
2227
2228
2229template <int dim, int spacedim>
2232 const unsigned int q_point) const
2233{
2234 Assert(this->update_flags & update_jacobian_2nd_derivatives,
2235 ExcAccessToUninitializedField("update_jacobian_2nd_derivatives"));
2236 Assert(present_cell.is_initialized(), ExcNotReinited());
2237 return this->mapping_output.jacobian_2nd_derivatives[q_point];
2238}
2239
2240
2241
2242template <int dim, int spacedim>
2243inline const std::vector<DerivativeForm<3, dim, spacedim>> &
2245{
2246 Assert(this->update_flags & update_jacobian_2nd_derivatives,
2247 ExcAccessToUninitializedField("update_jacobian_2nd_derivatives"));
2248 Assert(present_cell.is_initialized(), ExcNotReinited());
2249 return this->mapping_output.jacobian_2nd_derivatives;
2250}
2251
2252
2253
2254template <int dim, int spacedim>
2255inline const Tensor<4, spacedim> &
2257 const unsigned int q_point) const
2258{
2260 ExcAccessToUninitializedField(
2261 "update_jacobian_pushed_forward_2nd_derivatives"));
2262 Assert(present_cell.is_initialized(), ExcNotReinited());
2263 return this->mapping_output.jacobian_pushed_forward_2nd_derivatives[q_point];
2264}
2265
2266
2267
2268template <int dim, int spacedim>
2269inline const std::vector<Tensor<4, spacedim>> &
2271{
2273 ExcAccessToUninitializedField(
2274 "update_jacobian_pushed_forward_2nd_derivatives"));
2275 Assert(present_cell.is_initialized(), ExcNotReinited());
2276 return this->mapping_output.jacobian_pushed_forward_2nd_derivatives;
2277}
2278
2279
2280
2281template <int dim, int spacedim>
2284 const unsigned int q_point) const
2285{
2286 Assert(this->update_flags & update_jacobian_3rd_derivatives,
2287 ExcAccessToUninitializedField("update_jacobian_3rd_derivatives"));
2288 Assert(present_cell.is_initialized(), ExcNotReinited());
2289 return this->mapping_output.jacobian_3rd_derivatives[q_point];
2290}
2291
2292
2293
2294template <int dim, int spacedim>
2295inline const std::vector<DerivativeForm<4, dim, spacedim>> &
2297{
2298 Assert(this->update_flags & update_jacobian_3rd_derivatives,
2299 ExcAccessToUninitializedField("update_jacobian_3rd_derivatives"));
2300 Assert(present_cell.is_initialized(), ExcNotReinited());
2301 return this->mapping_output.jacobian_3rd_derivatives;
2302}
2303
2304
2305
2306template <int dim, int spacedim>
2307inline const Tensor<5, spacedim> &
2309 const unsigned int q_point) const
2310{
2312 ExcAccessToUninitializedField(
2313 "update_jacobian_pushed_forward_3rd_derivatives"));
2314 Assert(present_cell.is_initialized(), ExcNotReinited());
2315 return this->mapping_output.jacobian_pushed_forward_3rd_derivatives[q_point];
2316}
2317
2318
2319
2320template <int dim, int spacedim>
2321inline const std::vector<Tensor<5, spacedim>> &
2323{
2325 ExcAccessToUninitializedField(
2326 "update_jacobian_pushed_forward_3rd_derivatives"));
2327 Assert(present_cell.is_initialized(), ExcNotReinited());
2328 return this->mapping_output.jacobian_pushed_forward_3rd_derivatives;
2329}
2330
2331
2332
2333template <int dim, int spacedim>
2334inline const std::vector<DerivativeForm<1, spacedim, dim>> &
2336{
2337 Assert(this->update_flags & update_inverse_jacobians,
2338 ExcAccessToUninitializedField("update_inverse_jacobians"));
2339 Assert(present_cell.is_initialized(), ExcNotReinited());
2340 return this->mapping_output.inverse_jacobians;
2341}
2342
2343
2344
2345template <int dim, int spacedim>
2348{
2349 return {0U, dofs_per_cell};
2350}
2351
2352
2353
2354template <int dim, int spacedim>
2357 const unsigned int start_dof_index) const
2358{
2359 Assert(start_dof_index <= dofs_per_cell,
2360 ExcIndexRange(start_dof_index, 0, dofs_per_cell + 1));
2361 return {start_dof_index, dofs_per_cell};
2362}
2363
2364
2365
2366template <int dim, int spacedim>
2369 const unsigned int end_dof_index) const
2370{
2371 Assert(end_dof_index < dofs_per_cell,
2372 ExcIndexRange(end_dof_index, 0, dofs_per_cell));
2373 return {0U, end_dof_index + 1};
2374}
2375
2376
2377
2378template <int dim, int spacedim>
2381{
2382 return {0U, n_quadrature_points};
2383}
2384
2385
2386
2387template <int dim, int spacedim>
2388inline const Point<spacedim> &
2389FEValuesBase<dim, spacedim>::quadrature_point(const unsigned int q_point) const
2390{
2391 Assert(this->update_flags & update_quadrature_points,
2392 ExcAccessToUninitializedField("update_quadrature_points"));
2393 AssertIndexRange(q_point, this->mapping_output.quadrature_points.size());
2394 Assert(present_cell.is_initialized(), ExcNotReinited());
2395
2396 return this->mapping_output.quadrature_points[q_point];
2397}
2398
2399
2400
2401template <int dim, int spacedim>
2402inline double
2403FEValuesBase<dim, spacedim>::JxW(const unsigned int q_point) const
2404{
2405 Assert(this->update_flags & update_JxW_values,
2406 ExcAccessToUninitializedField("update_JxW_values"));
2407 AssertIndexRange(q_point, this->mapping_output.JxW_values.size());
2408 Assert(present_cell.is_initialized(), ExcNotReinited());
2409
2410 return this->mapping_output.JxW_values[q_point];
2411}
2412
2413
2414
2415template <int dim, int spacedim>
2417FEValuesBase<dim, spacedim>::jacobian(const unsigned int q_point) const
2418{
2419 Assert(this->update_flags & update_jacobians,
2420 ExcAccessToUninitializedField("update_jacobians"));
2421 AssertIndexRange(q_point, this->mapping_output.jacobians.size());
2422 Assert(present_cell.is_initialized(), ExcNotReinited());
2423
2424 return this->mapping_output.jacobians[q_point];
2425}
2426
2427
2428
2429template <int dim, int spacedim>
2431FEValuesBase<dim, spacedim>::jacobian_grad(const unsigned int q_point) const
2432{
2433 Assert(this->update_flags & update_jacobian_grads,
2434 ExcAccessToUninitializedField("update_jacobians_grads"));
2435 AssertIndexRange(q_point, this->mapping_output.jacobian_grads.size());
2436 Assert(present_cell.is_initialized(), ExcNotReinited());
2437
2438 return this->mapping_output.jacobian_grads[q_point];
2439}
2440
2441
2442
2443template <int dim, int spacedim>
2445FEValuesBase<dim, spacedim>::inverse_jacobian(const unsigned int q_point) const
2446{
2447 Assert(this->update_flags & update_inverse_jacobians,
2448 ExcAccessToUninitializedField("update_inverse_jacobians"));
2449 AssertIndexRange(q_point, this->mapping_output.inverse_jacobians.size());
2450 Assert(present_cell.is_initialized(), ExcNotReinited());
2451
2452 return this->mapping_output.inverse_jacobians[q_point];
2453}
2454
2455
2456
2457template <int dim, int spacedim>
2458inline const Tensor<1, spacedim> &
2459FEValuesBase<dim, spacedim>::normal_vector(const unsigned int q_point) const
2460{
2461 Assert(this->update_flags & update_normal_vectors,
2463 "update_normal_vectors")));
2464 AssertIndexRange(q_point, this->mapping_output.normal_vectors.size());
2465 Assert(present_cell.is_initialized(), ExcNotReinited());
2466
2467 return this->mapping_output.normal_vectors[q_point];
2468}
2469
2470#endif
2471
2473
2474#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 SmartPointer< const Mapping< dim, spacedim >, FEValuesBase< dim, spacedim > > mapping
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
const SmartPointer< const FiniteElement< dim, spacedim >, FEValuesBase< dim, spacedim > > fe
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
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 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:318
Definition point.h:111
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:502
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:503
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)
Definition exceptions.h:494
static ::ExceptionBase & ExcIndexRange(std::size_t arg1, std::size_t arg2, std::size_t arg3)
#define DeclException1(Exception1, type1, outsequence)
Definition exceptions.h:516
static ::ExceptionBase & ExcFEDontMatch()
static ::ExceptionBase & ExcFENotPrimitive()
typename ActiveSelector::cell_iterator cell_iterator
TriaIterator< CellAccessor< dim, spacedim > > cell_iterator
Definition tria.h:1550
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
spacedim & mapping
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