Reference documentation for deal.II version GIT 35969cdc9b 2023-12-09 01:10: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\}}\)
fe_values_base.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1998 - 2023 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
16 #ifndef dealii_fe_values_base_h
17 #define dealii_fe_values_base_h
18 
19 
20 #include <deal.II/base/config.h>
21 
24 #include <deal.II/base/point.h>
29 
32 
33 #include <deal.II/fe/fe.h>
37 #include <deal.II/fe/mapping.h>
39 
40 #include <deal.II/grid/tria.h>
42 
44 
46 
47 #include <algorithm>
48 #include <memory>
49 #include <type_traits>
50 
52 
152 template <int dim, int spacedim>
153 class FEValuesBase : public Subscriptor
154 {
155 public:
159  static constexpr unsigned int dimension = dim;
160 
164  static constexpr unsigned int space_dimension = spacedim;
165 
173  const unsigned int n_quadrature_points;
174 
184  const unsigned int max_n_quadrature_points;
185 
191  const unsigned int dofs_per_cell;
192 
193 
201  FEValuesBase(const unsigned int n_q_points,
202  const unsigned int dofs_per_cell,
206 
211  FEValuesBase &
212  operator=(const FEValuesBase &) = delete;
213 
218  FEValuesBase(const FEValuesBase &) = delete;
219 
223  virtual ~FEValuesBase() override;
224 
245  void
246  always_allow_check_for_cell_similarity(const bool allow);
247 
251 
273  const double &
274  shape_value(const unsigned int i, const unsigned int q_point) const;
275 
296  double
297  shape_value_component(const unsigned int i,
298  const unsigned int q_point,
299  const unsigned int component) const;
300 
326  const Tensor<1, spacedim> &
327  shape_grad(const unsigned int i, const unsigned int q_point) const;
328 
346  shape_grad_component(const unsigned int i,
347  const unsigned int q_point,
348  const unsigned int component) const;
349 
369  const Tensor<2, spacedim> &
370  shape_hessian(const unsigned int i, const unsigned int q_point) const;
371 
389  shape_hessian_component(const unsigned int i,
390  const unsigned int q_point,
391  const unsigned int component) const;
392 
412  const Tensor<3, spacedim> &
413  shape_3rd_derivative(const unsigned int i, const unsigned int q_point) const;
414 
432  shape_3rd_derivative_component(const unsigned int i,
433  const unsigned int q_point,
434  const unsigned int component) const;
435 
438 
484  template <typename Number>
485  void
486  get_function_values(const ReadVector<Number> &fe_function,
487  std::vector<Number> &values) const;
488 
502  template <typename Number>
503  void
504  get_function_values(const ReadVector<Number> &fe_function,
505  std::vector<Vector<Number>> &values) const;
506 
563  template <typename Number>
564  void
565  get_function_values(const ReadVector<Number> &fe_function,
567  std::vector<Number> &values) const;
568 
577  template <typename Number>
578  void
579  get_function_values(const ReadVector<Number> &fe_function,
581  std::vector<Vector<Number>> &values) const;
582 
583 
605  template <typename Number>
606  void
607  get_function_values(const ReadVector<Number> &fe_function,
609  ArrayView<std::vector<Number>> values,
610  const bool quadrature_points_fastest) const;
611 
614 
658  template <typename Number>
659  void
661  const ReadVector<Number> &fe_function,
662  std::vector<Tensor<1, spacedim, Number>> &gradients) const;
663 
680  template <typename Number>
681  void
683  const ReadVector<Number> &fe_function,
684  std::vector<std::vector<Tensor<1, spacedim, Number>>> &gradients) const;
685 
694  template <typename Number>
695  void
697  const ReadVector<Number> &fe_function,
699  std::vector<Tensor<1, spacedim, Number>> &gradients) const;
700 
709  template <typename Number>
710  void
712  const ReadVector<Number> &fe_function,
715  const bool quadrature_points_fastest = false) const;
716 
721 
760  template <typename Number>
761  void
763  const ReadVector<Number> &fe_function,
764  std::vector<Tensor<2, spacedim, Number>> &hessians) const;
765 
783  template <typename Number>
784  void
786  const ReadVector<Number> &fe_function,
787  std::vector<std::vector<Tensor<2, spacedim, Number>>> &hessians,
788  const bool quadrature_points_fastest = false) const;
789 
798  template <typename Number>
799  void
801  const ReadVector<Number> &fe_function,
803  std::vector<Tensor<2, spacedim, Number>> &hessians) const;
804 
813  template <typename Number>
814  void
816  const ReadVector<Number> &fe_function,
819  const bool quadrature_points_fastest = false) const;
820 
861  template <typename Number>
862  void
863  get_function_laplacians(const ReadVector<Number> &fe_function,
864  std::vector<Number> &laplacians) const;
865 
885  template <typename Number>
886  void
887  get_function_laplacians(const ReadVector<Number> &fe_function,
888  std::vector<Vector<Number>> &laplacians) const;
889 
898  template <typename Number>
899  void
901  const ReadVector<Number> &fe_function,
903  std::vector<Number> &laplacians) const;
904 
913  template <typename Number>
914  void
916  const ReadVector<Number> &fe_function,
918  std::vector<Vector<Number>> &laplacians) const;
919 
928  template <typename Number>
929  void
931  const ReadVector<Number> &fe_function,
933  std::vector<std::vector<Number>> &laplacians,
934  const bool quadrature_points_fastest = false) const;
935 
938 
978  template <typename Number>
979  void
981  const ReadVector<Number> &fe_function,
982  std::vector<Tensor<3, spacedim, Number>> &third_derivatives) const;
983 
1002  template <typename Number>
1003  void
1005  const ReadVector<Number> &fe_function,
1006  std::vector<std::vector<Tensor<3, spacedim, Number>>> &third_derivatives,
1007  const bool quadrature_points_fastest = false) const;
1008 
1017  template <typename Number>
1018  void
1020  const ReadVector<Number> &fe_function,
1022  std::vector<Tensor<3, spacedim, Number>> &third_derivatives) const;
1023 
1032  template <typename Number>
1033  void
1035  const ReadVector<Number> &fe_function,
1037  ArrayView<std::vector<Tensor<3, spacedim, Number>>> third_derivatives,
1038  const bool quadrature_points_fastest = false) const;
1042 
1068  dof_indices() const;
1069 
1103  dof_indices_starting_at(const unsigned int start_dof_index) const;
1104 
1136  dof_indices_ending_at(const unsigned int end_dof_index) const;
1137 
1141 
1165 
1172  const Point<spacedim> &
1173  quadrature_point(const unsigned int q_point) const;
1174 
1180  const std::vector<Point<spacedim>> &
1182 
1199  double
1200  JxW(const unsigned int q_point) const;
1201 
1205  const std::vector<double> &
1207 
1215  jacobian(const unsigned int q_point) const;
1216 
1223  const std::vector<DerivativeForm<1, dim, spacedim>> &
1224  get_jacobians() const;
1225 
1234  jacobian_grad(const unsigned int q_point) const;
1235 
1242  const std::vector<DerivativeForm<2, dim, spacedim>> &
1244 
1253  const Tensor<3, spacedim> &
1254  jacobian_pushed_forward_grad(const unsigned int q_point) const;
1255 
1262  const std::vector<Tensor<3, spacedim>> &
1264 
1273  jacobian_2nd_derivative(const unsigned int q_point) const;
1274 
1281  const std::vector<DerivativeForm<3, dim, spacedim>> &
1283 
1293  const Tensor<4, spacedim> &
1294  jacobian_pushed_forward_2nd_derivative(const unsigned int q_point) const;
1295 
1302  const std::vector<Tensor<4, spacedim>> &
1304 
1314  jacobian_3rd_derivative(const unsigned int q_point) const;
1315 
1322  const std::vector<DerivativeForm<4, dim, spacedim>> &
1324 
1334  const Tensor<5, spacedim> &
1335  jacobian_pushed_forward_3rd_derivative(const unsigned int q_point) const;
1336 
1343  const std::vector<Tensor<5, spacedim>> &
1345 
1353  inverse_jacobian(const unsigned int q_point) const;
1354 
1361  const std::vector<DerivativeForm<1, spacedim, dim>> &
1363 
1383  const Tensor<1, spacedim> &
1384  normal_vector(const unsigned int q_point) const;
1385 
1393  const std::vector<Tensor<1, spacedim>> &
1394  get_normal_vectors() const;
1395 
1399 
1411 
1422 
1434 
1435 
1446 
1450 
1455  const Mapping<dim, spacedim> &
1456  get_mapping() const;
1457 
1462  get_fe() const;
1463 
1467  UpdateFlags
1469 
1474  get_cell() const;
1475 
1482  get_cell_similarity() const;
1483 
1488  std::size_t
1489  memory_consumption() const;
1501  std::string,
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 "
1506  << "the <" << arg1
1507  << "> flag to be set, but it was apparently not specified "
1508  << "upon construction.");
1509 
1516  "FEValues object is not reinit'ed to any cell");
1517 
1526  "The FiniteElement you provided to FEValues and the FiniteElement that belongs "
1527  "to the DoFHandler that provided the cell iterator do not match.");
1534  int,
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?");
1541 
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.");
1552 
1553 protected:
1561  {
1562  public:
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.");
1572 
1577 
1581  template <bool lda>
1584 
1590 
1594  bool
1595  is_initialized() const;
1596 
1603  operator typename Triangulation<dim, spacedim>::cell_iterator() const;
1604 
1610  n_dofs_for_dof_handler() const;
1611 
1616  template <typename Number>
1617  void
1619  Vector<Number> &out) const;
1620 
1625  void
1627  Vector<IndexSet::value_type> &out) const;
1628 
1629  private:
1634  };
1635 
1642 
1650  boost::signals2::connection tria_listener_refinement;
1651 
1659  boost::signals2::connection tria_listener_mesh_transform;
1660 
1666  void
1668 
1678  void
1680  const typename Triangulation<dim, spacedim>::cell_iterator &cell);
1681 
1687 
1693  std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
1695 
1702 
1710 
1716  std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
1718 
1724  spacedim>
1726 
1727 
1732 
1741  UpdateFlags
1743 
1750 
1756  void
1758  const typename Triangulation<dim, spacedim>::cell_iterator &cell);
1759 
1760 private:
1765 
1770 
1771  // Make the view classes friends of this class, since they access internal
1772  // data.
1773  template <int, int>
1775  template <int, int>
1777  template <int, int, int>
1779  template <int, int, int>
1781 };
1782 
1783 #ifndef DOXYGEN
1784 
1785 /*---------------------- Inline functions: FEValuesBase ---------------------*/
1786 
1787 template <int dim, int spacedim>
1788 template <bool lda>
1792  : initialized(true)
1793  , cell(cell)
1794  , dof_handler(&cell->get_dof_handler())
1795  , level_dof_access(lda)
1796 {}
1797 
1798 
1799 
1800 template <int dim, int spacedim>
1803  const FEValuesExtractors::Scalar &scalar) const
1804 {
1805  AssertIndexRange(scalar.component, fe_values_views_cache.scalars.size());
1806 
1807  return fe_values_views_cache.scalars[scalar.component].value_or_initialize(
1808  [scalar, this]() {
1809  return FEValuesViews::Scalar<dim, spacedim>(*this, scalar.component);
1810  });
1811 }
1812 
1813 
1814 
1815 template <int dim, int spacedim>
1818  const FEValuesExtractors::Vector &vector) const
1819 {
1821  fe_values_views_cache.vectors.size());
1822 
1823  return fe_values_views_cache.vectors[vector.first_vector_component]
1824  .value_or_initialize([vector, this]() {
1826  *this, vector.first_vector_component);
1827  });
1828 }
1829 
1830 
1831 
1832 template <int dim, int spacedim>
1835  const FEValuesExtractors::SymmetricTensor<2> &tensor) const
1836 {
1837  Assert(
1838  tensor.first_tensor_component <
1839  fe_values_views_cache.symmetric_second_order_tensors.size(),
1841  0,
1842  fe_values_views_cache.symmetric_second_order_tensors.size()));
1843 
1844  return fe_values_views_cache
1846  .value_or_initialize([tensor, this]() {
1848  *this, tensor.first_tensor_component);
1849  });
1850 }
1851 
1852 
1853 
1854 template <int dim, int spacedim>
1857  const FEValuesExtractors::Tensor<2> &tensor) const
1858 {
1860  fe_values_views_cache.second_order_tensors.size());
1861 
1862  return fe_values_views_cache
1864  .value_or_initialize([tensor, this]() {
1866  *this, tensor.first_tensor_component);
1867  });
1868 }
1869 
1870 
1871 
1872 template <int dim, int spacedim>
1873 inline const double &
1874 FEValuesBase<dim, spacedim>::shape_value(const unsigned int i,
1875  const unsigned int q_point) const
1876 {
1877  AssertIndexRange(i, fe->n_dofs_per_cell());
1878  Assert(this->update_flags & update_values,
1879  ExcAccessToUninitializedField("update_values"));
1880  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
1881  Assert(present_cell.is_initialized(), ExcNotReinited());
1882  // if the entire FE is primitive,
1883  // then we can take a short-cut:
1884  if (fe->is_primitive())
1885  return this->finite_element_output.shape_values(i, q_point);
1886  else
1887  {
1888  // otherwise, use the mapping
1889  // between shape function
1890  // numbers and rows. note that
1891  // by the assertions above, we
1892  // know that this particular
1893  // shape function is primitive,
1894  // so we can call
1895  // system_to_component_index
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);
1901  }
1902 }
1903 
1904 
1905 
1906 template <int dim, int spacedim>
1907 inline double
1909  const unsigned int i,
1910  const unsigned int q_point,
1911  const unsigned int component) const
1912 {
1913  AssertIndexRange(i, fe->n_dofs_per_cell());
1914  Assert(this->update_flags & update_values,
1915  ExcAccessToUninitializedField("update_values"));
1916  AssertIndexRange(component, fe->n_components());
1917  Assert(present_cell.is_initialized(), ExcNotReinited());
1918 
1919  // check whether the shape function
1920  // is non-zero at all within
1921  // this component:
1922  if (fe->get_nonzero_components(i)[component] == false)
1923  return 0;
1924 
1925  // look up the right row in the
1926  // table and take the data from
1927  // there
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);
1932 }
1933 
1934 
1935 
1936 template <int dim, int spacedim>
1937 inline const Tensor<1, spacedim> &
1938 FEValuesBase<dim, spacedim>::shape_grad(const unsigned int i,
1939  const unsigned int q_point) const
1940 {
1941  AssertIndexRange(i, fe->n_dofs_per_cell());
1942  Assert(this->update_flags & update_gradients,
1943  ExcAccessToUninitializedField("update_gradients"));
1944  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
1945  Assert(present_cell.is_initialized(), ExcNotReinited());
1946  // if the entire FE is primitive,
1947  // then we can take a short-cut:
1948  if (fe->is_primitive())
1949  return this->finite_element_output.shape_gradients[i][q_point];
1950  else
1951  {
1952  // otherwise, use the mapping
1953  // between shape function
1954  // numbers and rows. note that
1955  // by the assertions above, we
1956  // know that this particular
1957  // shape function is primitive,
1958  // so we can call
1959  // system_to_component_index
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];
1965  }
1966 }
1967 
1968 
1969 
1970 template <int dim, int spacedim>
1971 inline Tensor<1, spacedim>
1973  const unsigned int i,
1974  const unsigned int q_point,
1975  const unsigned int component) const
1976 {
1977  AssertIndexRange(i, fe->n_dofs_per_cell());
1978  Assert(this->update_flags & update_gradients,
1979  ExcAccessToUninitializedField("update_gradients"));
1980  AssertIndexRange(component, fe->n_components());
1981  Assert(present_cell.is_initialized(), ExcNotReinited());
1982  // check whether the shape function
1983  // is non-zero at all within
1984  // this component:
1985  if (fe->get_nonzero_components(i)[component] == false)
1986  return Tensor<1, spacedim>();
1987 
1988  // look up the right row in the
1989  // table and take the data from
1990  // there
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];
1995 }
1996 
1997 
1998 
1999 template <int dim, int spacedim>
2000 inline const Tensor<2, spacedim> &
2001 FEValuesBase<dim, spacedim>::shape_hessian(const unsigned int i,
2002  const unsigned int q_point) const
2003 {
2004  AssertIndexRange(i, fe->n_dofs_per_cell());
2005  Assert(this->update_flags & update_hessians,
2006  ExcAccessToUninitializedField("update_hessians"));
2007  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
2008  Assert(present_cell.is_initialized(), ExcNotReinited());
2009  // if the entire FE is primitive,
2010  // then we can take a short-cut:
2011  if (fe->is_primitive())
2012  return this->finite_element_output.shape_hessians[i][q_point];
2013  else
2014  {
2015  // otherwise, use the mapping
2016  // between shape function
2017  // numbers and rows. note that
2018  // by the assertions above, we
2019  // know that this particular
2020  // shape function is primitive,
2021  // so we can call
2022  // system_to_component_index
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];
2028  }
2029 }
2030 
2031 
2032 
2033 template <int dim, int spacedim>
2034 inline Tensor<2, spacedim>
2036  const unsigned int i,
2037  const unsigned int q_point,
2038  const unsigned int component) const
2039 {
2040  AssertIndexRange(i, fe->n_dofs_per_cell());
2041  Assert(this->update_flags & update_hessians,
2042  ExcAccessToUninitializedField("update_hessians"));
2043  AssertIndexRange(component, fe->n_components());
2044  Assert(present_cell.is_initialized(), ExcNotReinited());
2045  // check whether the shape function
2046  // is non-zero at all within
2047  // this component:
2048  if (fe->get_nonzero_components(i)[component] == false)
2049  return Tensor<2, spacedim>();
2050 
2051  // look up the right row in the
2052  // table and take the data from
2053  // there
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];
2058 }
2059 
2060 
2061 
2062 template <int dim, int spacedim>
2063 inline const Tensor<3, spacedim> &
2065  const unsigned int i,
2066  const unsigned int q_point) const
2067 {
2068  AssertIndexRange(i, fe->n_dofs_per_cell());
2069  Assert(this->update_flags & update_3rd_derivatives,
2070  ExcAccessToUninitializedField("update_3rd_derivatives"));
2071  Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
2072  Assert(present_cell.is_initialized(), ExcNotReinited());
2073  // if the entire FE is primitive,
2074  // then we can take a short-cut:
2075  if (fe->is_primitive())
2076  return this->finite_element_output.shape_3rd_derivatives[i][q_point];
2077  else
2078  {
2079  // otherwise, use the mapping
2080  // between shape function
2081  // numbers and rows. note that
2082  // by the assertions above, we
2083  // know that this particular
2084  // shape function is primitive,
2085  // so we can call
2086  // system_to_component_index
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];
2092  }
2093 }
2094 
2095 
2096 
2097 template <int dim, int spacedim>
2098 inline Tensor<3, spacedim>
2100  const unsigned int i,
2101  const unsigned int q_point,
2102  const unsigned int component) const
2103 {
2104  AssertIndexRange(i, fe->n_dofs_per_cell());
2105  Assert(this->update_flags & update_3rd_derivatives,
2106  ExcAccessToUninitializedField("update_3rd_derivatives"));
2107  AssertIndexRange(component, fe->n_components());
2108  Assert(present_cell.is_initialized(), ExcNotReinited());
2109  // check whether the shape function
2110  // is non-zero at all within
2111  // this component:
2112  if (fe->get_nonzero_components(i)[component] == false)
2113  return Tensor<3, spacedim>();
2114 
2115  // look up the right row in the
2116  // table and take the data from
2117  // there
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];
2122 }
2123 
2124 
2125 
2126 template <int dim, int spacedim>
2127 inline const FiniteElement<dim, spacedim> &
2129 {
2130  return *fe;
2131 }
2132 
2133 
2134 
2135 template <int dim, int spacedim>
2136 inline const Mapping<dim, spacedim> &
2138 {
2139  return *mapping;
2140 }
2141 
2142 
2143 
2144 template <int dim, int spacedim>
2145 inline UpdateFlags
2147 {
2148  return this->update_flags;
2149 }
2150 
2151 
2152 
2153 template <int dim, int spacedim>
2154 inline const std::vector<Point<spacedim>> &
2156 {
2157  Assert(this->update_flags & update_quadrature_points,
2158  ExcAccessToUninitializedField("update_quadrature_points"));
2159  Assert(present_cell.is_initialized(), ExcNotReinited());
2160  return this->mapping_output.quadrature_points;
2161 }
2162 
2163 
2164 
2165 template <int dim, int spacedim>
2166 inline const std::vector<double> &
2168 {
2169  Assert(this->update_flags & update_JxW_values,
2170  ExcAccessToUninitializedField("update_JxW_values"));
2171  Assert(present_cell.is_initialized(), ExcNotReinited());
2172  return this->mapping_output.JxW_values;
2173 }
2174 
2175 
2176 
2177 template <int dim, int spacedim>
2178 inline const std::vector<DerivativeForm<1, dim, spacedim>> &
2180 {
2181  Assert(this->update_flags & update_jacobians,
2182  ExcAccessToUninitializedField("update_jacobians"));
2183  Assert(present_cell.is_initialized(), ExcNotReinited());
2184  return this->mapping_output.jacobians;
2185 }
2186 
2187 
2188 
2189 template <int dim, int spacedim>
2190 inline const std::vector<DerivativeForm<2, dim, spacedim>> &
2192 {
2193  Assert(this->update_flags & update_jacobian_grads,
2194  ExcAccessToUninitializedField("update_jacobians_grads"));
2195  Assert(present_cell.is_initialized(), ExcNotReinited());
2196  return this->mapping_output.jacobian_grads;
2197 }
2198 
2199 
2200 
2201 template <int dim, int spacedim>
2202 inline const Tensor<3, spacedim> &
2204  const unsigned int q_point) const
2205 {
2206  Assert(this->update_flags & update_jacobian_pushed_forward_grads,
2207  ExcAccessToUninitializedField("update_jacobian_pushed_forward_grads"));
2208  Assert(present_cell.is_initialized(), ExcNotReinited());
2209  return this->mapping_output.jacobian_pushed_forward_grads[q_point];
2210 }
2211 
2212 
2213 
2214 template <int dim, int spacedim>
2215 inline const std::vector<Tensor<3, spacedim>> &
2217 {
2218  Assert(this->update_flags & update_jacobian_pushed_forward_grads,
2219  ExcAccessToUninitializedField("update_jacobian_pushed_forward_grads"));
2220  Assert(present_cell.is_initialized(), ExcNotReinited());
2221  return this->mapping_output.jacobian_pushed_forward_grads;
2222 }
2223 
2224 
2225 
2226 template <int dim, int spacedim>
2227 inline const DerivativeForm<3, dim, spacedim> &
2229  const unsigned int q_point) const
2230 {
2231  Assert(this->update_flags & update_jacobian_2nd_derivatives,
2232  ExcAccessToUninitializedField("update_jacobian_2nd_derivatives"));
2233  Assert(present_cell.is_initialized(), ExcNotReinited());
2234  return this->mapping_output.jacobian_2nd_derivatives[q_point];
2235 }
2236 
2237 
2238 
2239 template <int dim, int spacedim>
2240 inline const std::vector<DerivativeForm<3, dim, spacedim>> &
2242 {
2243  Assert(this->update_flags & update_jacobian_2nd_derivatives,
2244  ExcAccessToUninitializedField("update_jacobian_2nd_derivatives"));
2245  Assert(present_cell.is_initialized(), ExcNotReinited());
2246  return this->mapping_output.jacobian_2nd_derivatives;
2247 }
2248 
2249 
2250 
2251 template <int dim, int spacedim>
2252 inline const Tensor<4, spacedim> &
2254  const unsigned int q_point) const
2255 {
2258  "update_jacobian_pushed_forward_2nd_derivatives"));
2259  Assert(present_cell.is_initialized(), ExcNotReinited());
2260  return this->mapping_output.jacobian_pushed_forward_2nd_derivatives[q_point];
2261 }
2262 
2263 
2264 
2265 template <int dim, int spacedim>
2266 inline const std::vector<Tensor<4, spacedim>> &
2268 {
2271  "update_jacobian_pushed_forward_2nd_derivatives"));
2272  Assert(present_cell.is_initialized(), ExcNotReinited());
2273  return this->mapping_output.jacobian_pushed_forward_2nd_derivatives;
2274 }
2275 
2276 
2277 
2278 template <int dim, int spacedim>
2279 inline const DerivativeForm<4, dim, spacedim> &
2281  const unsigned int q_point) const
2282 {
2283  Assert(this->update_flags & update_jacobian_3rd_derivatives,
2284  ExcAccessToUninitializedField("update_jacobian_3rd_derivatives"));
2285  Assert(present_cell.is_initialized(), ExcNotReinited());
2286  return this->mapping_output.jacobian_3rd_derivatives[q_point];
2287 }
2288 
2289 
2290 
2291 template <int dim, int spacedim>
2292 inline const std::vector<DerivativeForm<4, dim, spacedim>> &
2294 {
2295  Assert(this->update_flags & update_jacobian_3rd_derivatives,
2296  ExcAccessToUninitializedField("update_jacobian_3rd_derivatives"));
2297  Assert(present_cell.is_initialized(), ExcNotReinited());
2298  return this->mapping_output.jacobian_3rd_derivatives;
2299 }
2300 
2301 
2302 
2303 template <int dim, int spacedim>
2304 inline const Tensor<5, spacedim> &
2306  const unsigned int q_point) const
2307 {
2310  "update_jacobian_pushed_forward_3rd_derivatives"));
2311  Assert(present_cell.is_initialized(), ExcNotReinited());
2312  return this->mapping_output.jacobian_pushed_forward_3rd_derivatives[q_point];
2313 }
2314 
2315 
2316 
2317 template <int dim, int spacedim>
2318 inline const std::vector<Tensor<5, spacedim>> &
2320 {
2323  "update_jacobian_pushed_forward_3rd_derivatives"));
2324  Assert(present_cell.is_initialized(), ExcNotReinited());
2325  return this->mapping_output.jacobian_pushed_forward_3rd_derivatives;
2326 }
2327 
2328 
2329 
2330 template <int dim, int spacedim>
2331 inline const std::vector<DerivativeForm<1, spacedim, dim>> &
2333 {
2334  Assert(this->update_flags & update_inverse_jacobians,
2335  ExcAccessToUninitializedField("update_inverse_jacobians"));
2336  Assert(present_cell.is_initialized(), ExcNotReinited());
2337  return this->mapping_output.inverse_jacobians;
2338 }
2339 
2340 
2341 
2342 template <int dim, int spacedim>
2345 {
2346  return {0U, dofs_per_cell};
2347 }
2348 
2349 
2350 
2351 template <int dim, int spacedim>
2354  const unsigned int start_dof_index) const
2355 {
2356  Assert(start_dof_index <= dofs_per_cell,
2357  ExcIndexRange(start_dof_index, 0, dofs_per_cell + 1));
2358  return {start_dof_index, dofs_per_cell};
2359 }
2360 
2361 
2362 
2363 template <int dim, int spacedim>
2366  const unsigned int end_dof_index) const
2367 {
2368  Assert(end_dof_index < dofs_per_cell,
2369  ExcIndexRange(end_dof_index, 0, dofs_per_cell));
2370  return {0U, end_dof_index + 1};
2371 }
2372 
2373 
2374 
2375 template <int dim, int spacedim>
2378 {
2379  return {0U, n_quadrature_points};
2380 }
2381 
2382 
2383 
2384 template <int dim, int spacedim>
2385 inline const Point<spacedim> &
2386 FEValuesBase<dim, spacedim>::quadrature_point(const unsigned int q_point) const
2387 {
2388  Assert(this->update_flags & update_quadrature_points,
2389  ExcAccessToUninitializedField("update_quadrature_points"));
2390  AssertIndexRange(q_point, this->mapping_output.quadrature_points.size());
2391  Assert(present_cell.is_initialized(), ExcNotReinited());
2392 
2393  return this->mapping_output.quadrature_points[q_point];
2394 }
2395 
2396 
2397 
2398 template <int dim, int spacedim>
2399 inline double
2400 FEValuesBase<dim, spacedim>::JxW(const unsigned int q_point) const
2401 {
2402  Assert(this->update_flags & update_JxW_values,
2403  ExcAccessToUninitializedField("update_JxW_values"));
2404  AssertIndexRange(q_point, this->mapping_output.JxW_values.size());
2405  Assert(present_cell.is_initialized(), ExcNotReinited());
2406 
2407  return this->mapping_output.JxW_values[q_point];
2408 }
2409 
2410 
2411 
2412 template <int dim, int spacedim>
2413 inline const DerivativeForm<1, dim, spacedim> &
2414 FEValuesBase<dim, spacedim>::jacobian(const unsigned int q_point) const
2415 {
2416  Assert(this->update_flags & update_jacobians,
2417  ExcAccessToUninitializedField("update_jacobians"));
2418  AssertIndexRange(q_point, this->mapping_output.jacobians.size());
2419  Assert(present_cell.is_initialized(), ExcNotReinited());
2420 
2421  return this->mapping_output.jacobians[q_point];
2422 }
2423 
2424 
2425 
2426 template <int dim, int spacedim>
2427 inline const DerivativeForm<2, dim, spacedim> &
2428 FEValuesBase<dim, spacedim>::jacobian_grad(const unsigned int q_point) const
2429 {
2430  Assert(this->update_flags & update_jacobian_grads,
2431  ExcAccessToUninitializedField("update_jacobians_grads"));
2432  AssertIndexRange(q_point, this->mapping_output.jacobian_grads.size());
2433  Assert(present_cell.is_initialized(), ExcNotReinited());
2434 
2435  return this->mapping_output.jacobian_grads[q_point];
2436 }
2437 
2438 
2439 
2440 template <int dim, int spacedim>
2441 inline const DerivativeForm<1, spacedim, dim> &
2442 FEValuesBase<dim, spacedim>::inverse_jacobian(const unsigned int q_point) const
2443 {
2444  Assert(this->update_flags & update_inverse_jacobians,
2445  ExcAccessToUninitializedField("update_inverse_jacobians"));
2446  AssertIndexRange(q_point, this->mapping_output.inverse_jacobians.size());
2447  Assert(present_cell.is_initialized(), ExcNotReinited());
2448 
2449  return this->mapping_output.inverse_jacobians[q_point];
2450 }
2451 
2452 
2453 
2454 template <int dim, int spacedim>
2455 inline const Tensor<1, spacedim> &
2456 FEValuesBase<dim, spacedim>::normal_vector(const unsigned int q_point) const
2457 {
2458  Assert(this->update_flags & update_normal_vectors,
2460  "update_normal_vectors")));
2461  AssertIndexRange(q_point, this->mapping_output.normal_vectors.size());
2462  Assert(present_cell.is_initialized(), ExcNotReinited());
2463 
2464  return this->mapping_output.normal_vectors[q_point];
2465 }
2466 
2467 #endif
2468 
2470 
2471 #endif
CellIteratorContainer(const TriaIterator< DoFCellAccessor< dim, spacedim, lda >> &cell)
const DoFHandler< dim, spacedim > * dof_handler
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
UpdateFlags update_flags
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
Definition: vector.h:110
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:477
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:478
static ::ExceptionBase & ExcAccessToUninitializedField()
static ::ExceptionBase & ExcNeedsDoFHandler()
static ::ExceptionBase & ExcNotReinited()
static ::ExceptionBase & ExcShapeFunctionNotPrimitive(int arg1)
#define Assert(cond, exc)
Definition: exceptions.h:1631
static ::ExceptionBase & ExcFENotPrimitive()
static ::ExceptionBase & ExcFEDontMatch()
#define AssertIndexRange(index, range)
Definition: exceptions.h:1888
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:495
static ::ExceptionBase & ExcIndexRange(std::size_t arg1, std::size_t arg2, std::size_t arg3)
#define DeclException1(Exception1, type1, outsequence)
Definition: exceptions.h:517
static ::ExceptionBase & ExcAccessToUninitializedField(std::string arg1)
TriaIterator< CellAccessor< dim, spacedim > > cell_iterator
Definition: tria.h:1571
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
static const char U
boost::integer_range< IncrementableType > iota_view
Definition: iota_view.h:46
unsigned int global_dof_index
Definition: types.h:82
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