Reference documentation for deal.II version 8.1.0
fe_values.h
1 // ---------------------------------------------------------------------
2 // @f$Id: fe_values.h 31932 2013-12-08 02:15:54Z heister @f$
3 //
4 // Copyright (C) 1998 - 2013 by the deal.II authors
5 //
6 // This file is part of the deal.II library.
7 //
8 // The deal.II library is free software; you can use it, redistribute
9 // it, and/or modify it under the terms of the GNU Lesser General
10 // Public License as published by the Free Software Foundation; either
11 // version 2.1 of the License, or (at your option) any later version.
12 // The full text of the license can be found in the file LICENSE at
13 // the top level of the deal.II distribution.
14 //
15 // ---------------------------------------------------------------------
16 
17 #ifndef __deal2__fe_values_h
18 #define __deal2__fe_values_h
19 
20 
21 #include <deal.II/base/config.h>
23 #include <deal.II/base/subscriptor.h>
24 #include <deal.II/base/point.h>
25 #include <deal.II/base/derivative_form.h>
26 #include <deal.II/base/symmetric_tensor.h>
27 #include <deal.II/base/vector_slice.h>
28 #include <deal.II/base/quadrature.h>
29 #include <deal.II/base/table.h>
30 #include <deal.II/grid/tria.h>
31 #include <deal.II/grid/tria_iterator.h>
32 #include <deal.II/dofs/dof_handler.h>
33 #include <deal.II/dofs/dof_accessor.h>
34 #include <deal.II/hp/dof_handler.h>
35 #include <deal.II/fe/fe.h>
36 #include <deal.II/fe/fe_update_flags.h>
37 #include <deal.II/fe/fe_values_extractors.h>
38 #include <deal.II/fe/mapping.h>
39 #include <deal.II/multigrid/mg_dof_handler.h>
40 
41 #include <algorithm>
42 #include <memory>
43 
44 // dummy include in order to have the
45 // definition of PetscScalar available
46 // without including other PETSc stuff
47 #ifdef DEAL_II_WITH_PETSC
48 # include <petsc.h>
49 #endif
50 
52 
53 template <int dim> class Quadrature;
54 template <int dim, int spacedim=dim> class FEValuesBase;
55 
56 template <typename Number> class Vector;
57 template <typename Number> class BlockVector;
58 
59 
60 namespace internal
61 {
66  template <int dim>
67  struct CurlType;
68 
75  template <>
76  struct CurlType<1>
77  {
78  typedef Tensor<1,1> type;
79  };
80 
87  template <>
88  struct CurlType<2>
89  {
90  typedef Tensor<1,1> type;
91  };
92 
99  template <>
100  struct CurlType<3>
101  {
102  typedef Tensor<1,3> type;
103  };
104 }
105 
106 
107 
108 
129 namespace FEValuesViews
130 {
141  template <int dim, int spacedim=dim>
142  class Scalar
143  {
144  public:
150  typedef double value_type;
151 
157  typedef ::Tensor<1,spacedim> gradient_type;
158 
164  typedef ::Tensor<2,spacedim> hessian_type;
165 
171  {
181 
190  unsigned int row_index;
191  };
192 
196  Scalar ();
197 
203  Scalar (const FEValuesBase<dim,spacedim> &fe_values_base,
204  const unsigned int component);
205 
211 
223  value_type
224  value (const unsigned int shape_function,
225  const unsigned int q_point) const;
226 
235  gradient_type
236  gradient (const unsigned int shape_function,
237  const unsigned int q_point) const;
238 
247  hessian_type
248  hessian (const unsigned int shape_function,
249  const unsigned int q_point) const;
250 
261  template <class InputVector>
262  void get_function_values (const InputVector &fe_function,
263  std::vector<value_type> &values) const;
264 
275  template <class InputVector>
276  void get_function_gradients (const InputVector &fe_function,
277  std::vector<gradient_type> &gradients) const;
278 
289  template <class InputVector>
290  void get_function_hessians (const InputVector &fe_function,
291  std::vector<hessian_type> &hessians) const;
292 
304  template <class InputVector>
305  void get_function_laplacians (const InputVector &fe_function,
306  std::vector<value_type> &laplacians) const;
307 
308  private:
313 
318  const unsigned int component;
319 
323  std::vector<ShapeFunctionData> shape_function_data;
324  };
325 
326 
327 
355  template <int dim, int spacedim=dim>
356  class Vector
357  {
358  public:
364  typedef ::Tensor<1,spacedim> value_type;
365 
374  typedef ::Tensor<2,spacedim> gradient_type;
375 
386  typedef ::SymmetricTensor<2,spacedim> symmetric_gradient_type;
387 
393  typedef double divergence_type;
394 
401  typedef typename ::internal::CurlType<spacedim>::type curl_type;
402 
408  typedef ::Tensor<3,spacedim> hessian_type;
409 
415  {
425 
435  unsigned int row_index[spacedim];
436 
446  unsigned int single_nonzero_component_index;
447  };
448 
452  Vector ();
453 
462  Vector (const FEValuesBase<dim,spacedim> &fe_values_base,
463  const unsigned int first_vector_component);
464 
470 
485  value_type
486  value (const unsigned int shape_function,
487  const unsigned int q_point) const;
488 
500  gradient_type
501  gradient (const unsigned int shape_function,
502  const unsigned int q_point) const;
503 
517  symmetric_gradient_type
518  symmetric_gradient (const unsigned int shape_function,
519  const unsigned int q_point) const;
520 
529  divergence_type
530  divergence (const unsigned int shape_function,
531  const unsigned int q_point) const;
532 
547  curl_type
548  curl (const unsigned int shape_function,
549  const unsigned int q_point) const;
550 
559  hessian_type
560  hessian (const unsigned int shape_function,
561  const unsigned int q_point) const;
562 
573  template <class InputVector>
574  void get_function_values (const InputVector &fe_function,
575  std::vector<value_type> &values) const;
576 
587  template <class InputVector>
588  void get_function_gradients (const InputVector &fe_function,
589  std::vector<gradient_type> &gradients) const;
590 
606  template <class InputVector>
607  void
608  get_function_symmetric_gradients (const InputVector &fe_function,
609  std::vector<symmetric_gradient_type> &symmetric_gradients) const;
610 
622  template <class InputVector>
623  void get_function_divergences (const InputVector &fe_function,
624  std::vector<divergence_type> &divergences) const;
625 
637  template <class InputVector>
638  void get_function_curls (const InputVector &fe_function,
639  std::vector<curl_type> &curls) const;
640 
651  template <class InputVector>
652  void get_function_hessians (const InputVector &fe_function,
653  std::vector<hessian_type> &hessians) const;
654 
666  template <class InputVector>
667  void get_function_laplacians (const InputVector &fe_function,
668  std::vector<value_type> &laplacians) const;
669 
670  private:
675 
680  const unsigned int first_vector_component;
681 
685  std::vector<ShapeFunctionData> shape_function_data;
686  };
687 
688 
689  template <int rank, int dim, int spacedim = dim>
691 
714  template <int dim, int spacedim>
715  class SymmetricTensor<2,dim,spacedim>
716  {
717  public:
724  typedef ::SymmetricTensor<2, spacedim> value_type;
725 
735  typedef ::Tensor<1, spacedim> divergence_type;
736 
741  struct ShapeFunctionData
742  {
751  bool is_nonzero_shape_function_component[value_type::n_independent_components];
752 
762  unsigned int row_index[value_type::n_independent_components];
763 
773  unsigned int single_nonzero_component_index;
774  };
775 
779  SymmetricTensor();
780 
790  SymmetricTensor(const FEValuesBase<dim, spacedim> &fe_values_base,
791  const unsigned int first_tensor_component);
792 
798 
814  value_type
815  value (const unsigned int shape_function,
816  const unsigned int q_point) const;
817 
818 
830  divergence_type
831  divergence (const unsigned int shape_function,
832  const unsigned int q_point) const;
833 
844  template <class InputVector>
845  void get_function_values (const InputVector &fe_function,
846  std::vector<value_type> &values) const;
847 
862  template <class InputVector>
863  void get_function_divergences (const InputVector &fe_function,
864  std::vector<divergence_type> &divergences) const;
865 
866  private:
871 
876  const unsigned int first_tensor_component;
877 
881  std::vector<ShapeFunctionData> shape_function_data;
882  };
883 
884 
885  template <int rank, int dim, int spacedim = dim>
886  class Tensor;
887 
905  template <int dim, int spacedim>
906  class Tensor<2,dim,spacedim>
907  {
908  public:
909 
914  typedef ::Tensor<2, spacedim> value_type;
915 
919  typedef ::Tensor<1, spacedim> divergence_type;
920 
925  struct ShapeFunctionData
926  {
935  bool is_nonzero_shape_function_component[value_type::n_independent_components];
936 
946  unsigned int row_index[value_type::n_independent_components];
947 
957  unsigned int single_nonzero_component_index;
958  };
959 
963  Tensor();
964 
965 
975  Tensor(const FEValuesBase<dim, spacedim> &fe_values_base,
976  const unsigned int first_tensor_component);
977 
978 
983  Tensor &operator=(const Tensor<2, dim, spacedim> &);
984 
1000  value_type
1001  value (const unsigned int shape_function,
1002  const unsigned int q_point) const;
1003 
1015  divergence_type
1016  divergence (const unsigned int shape_function,
1017  const unsigned int q_point) const;
1018 
1029  template <class InputVector>
1030  void get_function_values (const InputVector &fe_function,
1031  std::vector<value_type> &values) const;
1032 
1033 
1048  template <class InputVector>
1049  void get_function_divergences (const InputVector &fe_function,
1050  std::vector<divergence_type> &divergences) const;
1051 
1052  private:
1057 
1062  const unsigned int first_tensor_component;
1063 
1067  std::vector<ShapeFunctionData> shape_function_data;
1068  };
1069 
1070 }
1071 
1072 
1073 namespace internal
1074 {
1075  namespace FEValuesViews
1076  {
1084  template <int dim, int spacedim>
1085  struct Cache
1086  {
1091  std::vector<::FEValuesViews::Scalar<dim,spacedim> > scalars;
1092  std::vector<::FEValuesViews::Vector<dim,spacedim> > vectors;
1093  std::vector<::FEValuesViews::SymmetricTensor<2,dim,spacedim> >
1094  symmetric_second_order_tensors;
1095  std::vector<::FEValuesViews::Tensor<2,dim,spacedim> >
1096  second_order_tensors;
1097 
1101  Cache (const FEValuesBase<dim,spacedim> &fe_values);
1102  };
1103  }
1104 }
1105 
1106 
1107 
1108 //TODO: Add access to mapping values to FEValuesBase
1109 //TODO: Several FEValuesBase of a system should share Mapping
1110 
1126 template <int dim, int spacedim=dim>
1128 {
1129 public:
1133  void initialize (const unsigned int n_quadrature_points,
1134  const FiniteElement<dim,spacedim> &fe,
1135  const UpdateFlags flags);
1136 
1155 
1160  typedef std::vector<std::vector<Tensor<1,spacedim> > > GradientVector;
1161 
1165  typedef std::vector<std::vector<Tensor<2,spacedim> > > HessianVector;
1166 
1171  ShapeVector shape_values;
1172 
1178  GradientVector shape_gradients;
1179 
1185  HessianVector shape_hessians;
1186 
1200  std::vector<double> JxW_values;
1201 
1205  std::vector< DerivativeForm<1,dim,spacedim> > jacobians;
1206 
1211  std::vector<DerivativeForm<2,dim,spacedim> > jacobian_grads;
1212 
1216  std::vector<DerivativeForm<1,spacedim,dim> > inverse_jacobians;
1217 
1223  std::vector<Point<spacedim> > quadrature_points;
1224 
1229  std::vector<Point<spacedim> > normal_vectors;
1230 
1235  std::vector<Tensor<1,spacedim> > boundary_forms;
1236 
1272  std::vector<unsigned int> shape_function_to_row_table;
1273 
1278 };
1279 
1280 
1393 template <int dim, int spacedim>
1394 class FEValuesBase : protected FEValuesData<dim,spacedim>,
1395  public Subscriptor
1396 {
1397 public:
1401  static const unsigned int dimension = dim;
1402 
1406  static const unsigned int space_dimension = spacedim;
1407 
1411  const unsigned int n_quadrature_points;
1412 
1418  const unsigned int dofs_per_cell;
1419 
1420 
1428  FEValuesBase (const unsigned int n_q_points,
1429  const unsigned int dofs_per_cell,
1430  const UpdateFlags update_flags,
1433 
1434 
1438  ~FEValuesBase ();
1440 
1441 
1460  const double &shape_value (const unsigned int function_no,
1461  const unsigned int point_no) const;
1462 
1481  double shape_value_component (const unsigned int function_no,
1482  const unsigned int point_no,
1483  const unsigned int component) const;
1484 
1503  const Tensor<1,spacedim> &
1504  shape_grad (const unsigned int function_no,
1505  const unsigned int quadrature_point) const;
1506 
1522  shape_grad_component (const unsigned int function_no,
1523  const unsigned int point_no,
1524  const unsigned int component) const;
1525 
1543  const Tensor<2,spacedim> &
1544  shape_hessian (const unsigned int function_no,
1545  const unsigned int point_no) const;
1546 
1550  const Tensor<2,spacedim> &
1551  shape_2nd_derivative (const unsigned int function_no,
1552  const unsigned int point_no) const DEAL_II_DEPRECATED;
1553 
1554 
1570  shape_hessian_component (const unsigned int function_no,
1571  const unsigned int point_no,
1572  const unsigned int component) const;
1573 
1578  shape_2nd_derivative_component (const unsigned int function_no,
1579  const unsigned int point_no,
1580  const unsigned int component) const DEAL_II_DEPRECATED;
1581 
1582 
1584 
1586 
1619  template <class InputVector, typename number>
1620  void get_function_values (const InputVector &fe_function,
1621  std::vector<number> &values) const;
1622 
1634  template <class InputVector, typename number>
1635  void get_function_values (const InputVector &fe_function,
1636  std::vector<Vector<number> > &values) const;
1637 
1654  template <class InputVector, typename number>
1655  void get_function_values (const InputVector &fe_function,
1656  const VectorSlice<const std::vector<types::global_dof_index> > &indices,
1657  std::vector<number> &values) const;
1658 
1678  template <class InputVector, typename number>
1679  void get_function_values (const InputVector &fe_function,
1680  const VectorSlice<const std::vector<types::global_dof_index> > &indices,
1681  std::vector<Vector<number> > &values) const;
1682 
1683 
1712  template <class InputVector>
1713  void get_function_values (const InputVector &fe_function,
1714  const VectorSlice<const std::vector<types::global_dof_index> > &indices,
1715  VectorSlice<std::vector<std::vector<double> > > values,
1716  const bool quadrature_points_fastest) const;
1717 
1719 
1721 
1754  template <class InputVector>
1755  void get_function_gradients (const InputVector &fe_function,
1756  std::vector<Tensor<1,spacedim> > &gradients) const;
1757 
1772  template <class InputVector>
1773  void get_function_gradients (const InputVector &fe_function,
1774  std::vector<std::vector<Tensor<1,spacedim> > > &gradients) const;
1775 
1780  template <class InputVector>
1781  void get_function_gradients (const InputVector &fe_function,
1782  const VectorSlice<const std::vector<types::global_dof_index> > &indices,
1783  std::vector<Tensor<1,spacedim> > &gradients) const;
1784 
1789  template <class InputVector>
1790  void get_function_gradients (const InputVector &fe_function,
1791  const VectorSlice<const std::vector<types::global_dof_index> > &indices,
1792  VectorSlice<std::vector<std::vector<Tensor<1,spacedim> > > > gradients,
1793  bool quadrature_points_fastest = false) const;
1794 
1798  template <class InputVector>
1799  void get_function_grads (const InputVector &fe_function,
1800  std::vector<Tensor<1,spacedim> > &gradients) const DEAL_II_DEPRECATED;
1801 
1805  template <class InputVector>
1806  void get_function_grads (const InputVector &fe_function,
1807  std::vector<std::vector<Tensor<1,spacedim> > > &gradients) const DEAL_II_DEPRECATED;
1808 
1812  template <class InputVector>
1813  void get_function_grads (const InputVector &fe_function,
1814  const VectorSlice<const std::vector<types::global_dof_index> > &indices,
1815  std::vector<Tensor<1,spacedim> > &gradients) const DEAL_II_DEPRECATED;
1816 
1820  template <class InputVector>
1821  void get_function_grads (const InputVector &fe_function,
1822  const VectorSlice<const std::vector<types::global_dof_index> > &indices,
1823  std::vector<std::vector<Tensor<1,spacedim> > > &gradients,
1824  bool quadrature_points_fastest = false) const DEAL_II_DEPRECATED;
1825 
1827 
1829 
1863  template <class InputVector>
1864  void
1865  get_function_hessians (const InputVector &fe_function,
1866  std::vector<Tensor<2,spacedim> > &hessians) const;
1867 
1883  template <class InputVector>
1884  void
1885  get_function_hessians (const InputVector &fe_function,
1886  std::vector<std::vector<Tensor<2,spacedim> > > &hessians,
1887  bool quadrature_points_fastest = false) const;
1888 
1893  template <class InputVector>
1894  void get_function_hessians (
1895  const InputVector &fe_function,
1896  const VectorSlice<const std::vector<types::global_dof_index> > &indices,
1897  std::vector<Tensor<2,spacedim> > &hessians) const;
1898 
1903  template <class InputVector>
1904  void get_function_hessians (
1905  const InputVector &fe_function,
1906  const VectorSlice<const std::vector<types::global_dof_index> > &indices,
1907  VectorSlice<std::vector<std::vector<Tensor<2,spacedim> > > > hessians,
1908  bool quadrature_points_fastest = false) const;
1909 
1913  template <class InputVector>
1914  void
1915  get_function_2nd_derivatives (const InputVector &,
1916  std::vector<Tensor<2,spacedim> > &) const DEAL_II_DEPRECATED;
1917 
1921  template <class InputVector>
1922  void
1923  get_function_2nd_derivatives (const InputVector &,
1924  std::vector<std::vector<Tensor<2,spacedim> > > &,
1925  bool = false) const DEAL_II_DEPRECATED;
1926 
1964  template <class InputVector, typename number>
1965  void
1966  get_function_laplacians (const InputVector &fe_function,
1967  std::vector<number> &laplacians) const;
1968 
1986  template <class InputVector, typename number>
1987  void
1988  get_function_laplacians (const InputVector &fe_function,
1989  std::vector<Vector<number> > &laplacians) const;
1990 
1995  template <class InputVector, typename number>
1997  const InputVector &fe_function,
1998  const VectorSlice<const std::vector<types::global_dof_index> > &indices,
1999  std::vector<number> &laplacians) const;
2000 
2005  template <class InputVector, typename number>
2007  const InputVector &fe_function,
2008  const VectorSlice<const std::vector<types::global_dof_index> > &indices,
2009  std::vector<Vector<number> > &laplacians) const;
2010 
2015  template <class InputVector, typename number>
2017  const InputVector &fe_function,
2018  const VectorSlice<const std::vector<types::global_dof_index> > &indices,
2019  std::vector<std::vector<number> > &laplacians,
2020  bool quadrature_points_fastest = false) const;
2022 
2024 
2025 
2029  const Point<spacedim> &quadrature_point (const unsigned int i) const;
2030 
2034  const std::vector<Point<spacedim> > &get_quadrature_points () const;
2035 
2049  double JxW (const unsigned int quadrature_point) const;
2050 
2054  const std::vector<double> &get_JxW_values () const;
2055 
2060  const DerivativeForm<1,dim,spacedim> &jacobian (const unsigned int quadrature_point) const;
2061 
2065  const std::vector<DerivativeForm<1,dim,spacedim> > &get_jacobians () const;
2066 
2072  const DerivativeForm<2,dim,spacedim> &jacobian_grad (const unsigned int quadrature_point) const;
2073 
2077  const std::vector<DerivativeForm<2,dim,spacedim> > &get_jacobian_grads () const;
2078 
2083  const DerivativeForm<1,spacedim,dim> &inverse_jacobian (const unsigned int quadrature_point) const;
2084 
2088  const std::vector<DerivativeForm<1,spacedim,dim> > &get_inverse_jacobians () const;
2098  const Point<spacedim> &normal_vector (const unsigned int i) const;
2099 
2105  const std::vector<Point<spacedim> > &get_normal_vectors () const;
2106 
2111  void transform (std::vector<Tensor<1,spacedim> > &transformed,
2112  const std::vector<Tensor<1,dim> > &original,
2113  MappingType mapping) const;
2114 
2121  const Point<spacedim> &cell_normal_vector (const unsigned int i) const DEAL_II_DEPRECATED;
2122 
2128  const std::vector<Point<spacedim> > &get_cell_normal_vectors () const DEAL_II_DEPRECATED;
2129 
2131 
2133 
2134 
2141  const FEValuesViews::Scalar<dim,spacedim> &
2142  operator[] (const FEValuesExtractors::Scalar &scalar) const;
2143 
2151  const FEValuesViews::Vector<dim,spacedim> &
2152  operator[] (const FEValuesExtractors::Vector &vector) const;
2153 
2161  const FEValuesViews::SymmetricTensor<2,dim,spacedim> &
2162  operator[] (const FEValuesExtractors::SymmetricTensor<2> &tensor) const;
2163 
2164 
2172  const FEValuesViews::Tensor<2,dim,spacedim> &
2173  operator[] (const FEValuesExtractors::Tensor<2> &tensor) const;
2174 
2176 
2178 
2179 
2183  const Mapping<dim,spacedim> &get_mapping () const;
2184 
2188  const FiniteElement<dim,spacedim> &get_fe () const;
2189 
2193  UpdateFlags get_update_flags () const;
2194 
2198  const typename Triangulation<dim,spacedim>::cell_iterator get_cell () const;
2199 
2205  CellSimilarity::Similarity get_cell_similarity () const;
2206 
2211  std::size_t memory_consumption () const;
2213 
2214 
2221  DeclException1 (ExcAccessToUninitializedField,
2222  char *,
2223  << ("You are requesting information from an FEValues/FEFaceValues/FESubfaceValues "
2224  "object for which this kind of information has not been computed. What "
2225  "information these objects compute is determined by the update_* flags you "
2226  "pass to the constructor. Here, the operation you are attempting requires "
2227  "the <")
2228  << arg1
2229  << "> flag to be set, but it was apparently not specified upon construction.");
2235  DeclException0 (ExcCannotInitializeField);
2241  DeclException0 (ExcInvalidUpdateFlag);
2247  DeclException0 (ExcFEDontMatch);
2253  DeclException1 (ExcShapeFunctionNotPrimitive,
2254  int,
2255  << "The shape function with index " << arg1
2256  << " is not primitive, i.e. it is vector-valued and "
2257  << "has more than one non-zero vector component. This "
2258  << "function cannot be called for these shape functions. "
2259  << "Maybe you want to use the same function with the "
2260  << "_component suffix?");
2266  DeclException0 (ExcFENotPrimitive);
2267 
2268 protected:
2297  class CellIteratorBase;
2298 
2303  template <typename CI> class CellIterator;
2304  class TriaCellIterator;
2305 
2311  std::auto_ptr<const CellIteratorBase> present_cell;
2312 
2320  boost::signals2::connection tria_listener;
2321 
2327  void invalidate_present_cell ();
2328 
2338  void
2339  maybe_invalidate_previous_present_cell (const typename Triangulation<dim,spacedim>::cell_iterator &cell);
2340 
2344  const SmartPointer<const Mapping<dim,spacedim>,FEValuesBase<dim,spacedim> > mapping;
2345 
2349  const SmartPointer<const FiniteElement<dim,spacedim>,FEValuesBase<dim,spacedim> > fe;
2350 
2351 
2355  SmartPointer<typename Mapping<dim,spacedim>::InternalDataBase,FEValuesBase<dim,spacedim> > mapping_data;
2356 
2360  SmartPointer<typename Mapping<dim,spacedim>::InternalDataBase,FEValuesBase<dim,spacedim> > fe_data;
2361 
2371 
2378 
2384  void
2385  check_cell_similarity (const typename Triangulation<dim,spacedim>::cell_iterator &cell);
2386 
2387 private:
2392  FEValuesBase (const FEValuesBase &);
2393 
2398  FEValuesBase &operator= (const FEValuesBase &);
2399 
2404 
2409  template <int, int> friend class FEValuesViews::Scalar;
2410  template <int, int> friend class FEValuesViews::Vector;
2411  template <int, int, int> friend class FEValuesViews::SymmetricTensor;
2412  template <int, int, int> friend class FEValuesViews::Tensor;
2413 };
2414 
2415 
2416 
2427 template <int dim, int spacedim=dim>
2428 class FEValues : public FEValuesBase<dim,spacedim>
2429 {
2430 public:
2435  static const unsigned int integral_dimension = dim;
2436 
2441  FEValues (const Mapping<dim,spacedim> &mapping,
2442  const FiniteElement<dim,spacedim> &fe,
2443  const Quadrature<dim> &quadrature,
2444  const UpdateFlags update_flags);
2445 
2450  const Quadrature<dim> &quadrature,
2451  const UpdateFlags update_flags);
2452 
2459  template <class DH, bool level_dof_access>
2460  void reinit (const TriaIterator<DoFCellAccessor<DH,level_dof_access> > cell);
2461 
2474  void reinit (const typename Triangulation<dim,spacedim>::cell_iterator &cell);
2475 
2480  const Quadrature<dim> &get_quadrature () const;
2481 
2486  std::size_t memory_consumption () const;
2487 
2501  const FEValues<dim,spacedim> &get_present_fe_values () const;
2502 
2503 private:
2508 
2512  void initialize (const UpdateFlags update_flags);
2513 
2520  void do_reinit ();
2521 };
2522 
2523 
2534 template <int dim, int spacedim=dim>
2535 class FEFaceValuesBase : public FEValuesBase<dim,spacedim>
2536 {
2537 public:
2542  static const unsigned int integral_dimension = dim-1;
2543 
2555  FEFaceValuesBase (const unsigned int n_q_points,
2556  const unsigned int dofs_per_cell,
2557  const UpdateFlags update_flags,
2560  const Quadrature<dim-1>& quadrature);
2561 
2566  const Tensor<1,spacedim> &boundary_form (const unsigned int i) const;
2567 
2572  const std::vector<Tensor<1,spacedim> > &get_boundary_forms () const;
2573 
2578  unsigned int get_face_index() const;
2579 
2584  const Quadrature<dim-1> & get_quadrature () const;
2585 
2590  std::size_t memory_consumption () const;
2591 
2592 protected:
2593 
2598  unsigned int present_face_index;
2599 
2603  const Quadrature<dim-1> quadrature;
2604 };
2605 
2606 
2607 
2622 template <int dim, int spacedim=dim>
2623 class FEFaceValues : public FEFaceValuesBase<dim,spacedim>
2624 {
2625 public:
2630  static const unsigned int dimension = dim;
2631 
2632  static const unsigned int space_dimension = spacedim;
2633 
2638  static const unsigned int integral_dimension = dim-1;
2639 
2646  const Quadrature<dim-1> &quadrature,
2647  const UpdateFlags update_flags);
2648 
2653  const Quadrature<dim-1> &quadrature,
2654  const UpdateFlags update_flags);
2655 
2660  template <class DH, bool level_dof_access>
2661  void reinit (const TriaIterator<DoFCellAccessor<DH,level_dof_access> > cell,
2662  const unsigned int face_no);
2663 
2676  void reinit (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
2677  const unsigned int face_no);
2678 
2692  const FEFaceValues<dim,spacedim> &get_present_fe_values () const;
2693 private:
2694 
2698  void initialize (const UpdateFlags update_flags);
2699 
2706  void do_reinit (const unsigned int face_no);
2707 };
2708 
2709 
2727 template <int dim, int spacedim=dim>
2728 class FESubfaceValues : public FEFaceValuesBase<dim,spacedim>
2729 {
2730 public:
2734  static const unsigned int dimension = dim;
2735 
2739  static const unsigned int space_dimension = spacedim;
2740 
2745  static const unsigned int integral_dimension = dim-1;
2746 
2753  const Quadrature<dim-1> &face_quadrature,
2754  const UpdateFlags update_flags);
2755 
2760  const Quadrature<dim-1> &face_quadrature,
2761  const UpdateFlags update_flags);
2762 
2769  template <class DH, bool level_dof_access>
2770  void reinit (const TriaIterator<DoFCellAccessor<DH,level_dof_access> > cell,
2771  const unsigned int face_no,
2772  const unsigned int subface_no);
2773 
2786  void reinit (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
2787  const unsigned int face_no,
2788  const unsigned int subface_no);
2789 
2803  const FESubfaceValues<dim,spacedim> &get_present_fe_values () const;
2804 
2810  DeclException0 (ExcReinitCalledWithBoundaryFace);
2811 
2817  DeclException0 (ExcFaceHasNoSubfaces);
2818 
2819 private:
2820 
2824  void initialize (const UpdateFlags update_flags);
2825 
2832  void do_reinit (const unsigned int face_no,
2833  const unsigned int subface_no);
2834 };
2835 
2836 
2837 #ifndef DOXYGEN
2838 
2839 
2840 /*------------------------ Inline functions: namespace FEValuesViews --------*/
2841 
2842 namespace FEValuesViews
2843 {
2844  template <int dim, int spacedim>
2845  inline
2846  typename Scalar<dim,spacedim>::value_type
2847  Scalar<dim,spacedim>::value (const unsigned int shape_function,
2848  const unsigned int q_point) const
2849  {
2850  typedef FEValuesBase<dim,spacedim> FVB;
2851  Assert (shape_function < fe_values.fe->dofs_per_cell,
2852  ExcIndexRange (shape_function, 0, fe_values.fe->dofs_per_cell));
2853  Assert (fe_values.update_flags & update_values,
2854  typename FVB::ExcAccessToUninitializedField("update_values"));
2855 
2856  // an adaptation of the FEValuesBase::shape_value_component function
2857  // except that here we know the component as fixed and we have
2858  // pre-computed and cached a bunch of information. see the comments there
2859  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
2860  return fe_values.shape_values(shape_function_data[shape_function]
2861  .row_index,
2862  q_point);
2863  else
2864  return 0;
2865  }
2866 
2867 
2868 
2869 
2870  template <int dim, int spacedim>
2871  inline
2872  typename Scalar<dim,spacedim>::gradient_type
2873  Scalar<dim,spacedim>::gradient (const unsigned int shape_function,
2874  const unsigned int q_point) const
2875  {
2876  typedef FEValuesBase<dim,spacedim> FVB;
2877  Assert (shape_function < fe_values.fe->dofs_per_cell,
2878  ExcIndexRange (shape_function, 0, fe_values.fe->dofs_per_cell));
2879  Assert (fe_values.update_flags & update_gradients,
2880  typename FVB::ExcAccessToUninitializedField("update_gradients"));
2881 
2882  // an adaptation of the
2883  // FEValuesBase::shape_grad_component
2884  // function except that here we know the
2885  // component as fixed and we have
2886  // pre-computed and cached a bunch of
2887  // information. see the comments there
2888  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
2889  return fe_values.shape_gradients[shape_function_data[shape_function]
2890  .row_index][q_point];
2891  else
2892  return gradient_type();
2893  }
2894 
2895 
2896 
2897  template <int dim, int spacedim>
2898  inline
2899  typename Scalar<dim,spacedim>::hessian_type
2900  Scalar<dim,spacedim>::hessian (const unsigned int shape_function,
2901  const unsigned int q_point) const
2902  {
2903  typedef FEValuesBase<dim,spacedim> FVB;
2904  Assert (shape_function < fe_values.fe->dofs_per_cell,
2905  ExcIndexRange (shape_function, 0, fe_values.fe->dofs_per_cell));
2906  Assert (fe_values.update_flags & update_hessians,
2907  typename FVB::ExcAccessToUninitializedField("update_hessians"));
2908 
2909  // an adaptation of the
2910  // FEValuesBase::shape_grad_component
2911  // function except that here we know the
2912  // component as fixed and we have
2913  // pre-computed and cached a bunch of
2914  // information. see the comments there
2915  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
2916  return fe_values.shape_hessians[shape_function_data[shape_function].row_index][q_point];
2917  else
2918  return hessian_type();
2919  }
2920 
2921 
2922 
2923  template <int dim, int spacedim>
2924  inline
2926  Vector<dim,spacedim>::value (const unsigned int shape_function,
2927  const unsigned int q_point) const
2928  {
2929  typedef FEValuesBase<dim,spacedim> FVB;
2930  Assert (shape_function < fe_values.fe->dofs_per_cell,
2931  ExcIndexRange (shape_function, 0, fe_values.fe->dofs_per_cell));
2932  Assert (fe_values.update_flags & update_values,
2933  typename FVB::ExcAccessToUninitializedField("update_values"));
2934 
2935  // same as for the scalar case except
2936  // that we have one more index
2937  const int snc = shape_function_data[shape_function].single_nonzero_component;
2938  if (snc == -2)
2939  return value_type();
2940  else if (snc != -1)
2941  {
2942  value_type return_value;
2943  return_value[shape_function_data[shape_function].single_nonzero_component_index]
2944  = fe_values.shape_values(snc,q_point);
2945  return return_value;
2946  }
2947  else
2948  {
2949  value_type return_value;
2950  for (unsigned int d=0; d<dim; ++d)
2951  if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
2952  return_value[d]
2953  = fe_values.shape_values(shape_function_data[shape_function].row_index[d],q_point);
2954 
2955  return return_value;
2956  }
2957  }
2958 
2959 
2960 
2961  template <int dim, int spacedim>
2962  inline
2964  Vector<dim,spacedim>::gradient (const unsigned int shape_function,
2965  const unsigned int q_point) const
2966  {
2967  typedef FEValuesBase<dim,spacedim> FVB;
2968  Assert (shape_function < fe_values.fe->dofs_per_cell,
2969  ExcIndexRange (shape_function, 0, fe_values.fe->dofs_per_cell));
2970  Assert (fe_values.update_flags & update_gradients,
2971  typename FVB::ExcAccessToUninitializedField("update_gradients"));
2972 
2973  // same as for the scalar case except
2974  // that we have one more index
2975  const int snc = shape_function_data[shape_function].single_nonzero_component;
2976  if (snc == -2)
2977  return gradient_type();
2978  else if (snc != -1)
2979  {
2980  gradient_type return_value;
2981  return_value[shape_function_data[shape_function].single_nonzero_component_index]
2982  = fe_values.shape_gradients[snc][q_point];
2983  return return_value;
2984  }
2985  else
2986  {
2987  gradient_type return_value;
2988  for (unsigned int d=0; d<dim; ++d)
2989  if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
2990  return_value[d]
2991  = fe_values.shape_gradients[shape_function_data[shape_function].row_index[d]][q_point];
2992 
2993  return return_value;
2994  }
2995  }
2996 
2997 
2998 
2999  template <int dim, int spacedim>
3000  inline
3002  Vector<dim,spacedim>::divergence (const unsigned int shape_function,
3003  const unsigned int q_point) const
3004  {
3005  // this function works like in
3006  // the case above
3007  typedef FEValuesBase<dim,spacedim> FVB;
3008  Assert (shape_function < fe_values.fe->dofs_per_cell,
3009  ExcIndexRange (shape_function, 0, fe_values.fe->dofs_per_cell));
3010  Assert (fe_values.update_flags & update_gradients,
3011  typename FVB::ExcAccessToUninitializedField("update_gradients"));
3012 
3013  // same as for the scalar case except
3014  // that we have one more index
3015  const int snc = shape_function_data[shape_function].single_nonzero_component;
3016  if (snc == -2)
3017  return divergence_type();
3018  else if (snc != -1)
3019  return
3020  fe_values.shape_gradients[snc][q_point][shape_function_data[shape_function].single_nonzero_component_index];
3021  else
3022  {
3023  divergence_type return_value = 0;
3024  for (unsigned int d=0; d<dim; ++d)
3025  if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
3026  return_value
3027  += fe_values.shape_gradients[shape_function_data[shape_function].row_index[d]][q_point][d];
3028 
3029  return return_value;
3030  }
3031  }
3032 
3033 
3034 
3035  template <int dim, int spacedim>
3036  inline
3038  Vector<dim,spacedim>::curl (const unsigned int shape_function, const unsigned int q_point) const
3039  {
3040  // this function works like in the case above
3041  typedef FEValuesBase<dim,spacedim> FVB;
3042 
3043  Assert (shape_function < fe_values.fe->dofs_per_cell,
3044  ExcIndexRange (shape_function, 0, fe_values.fe->dofs_per_cell));
3045  Assert (fe_values.update_flags & update_gradients,
3046  typename FVB::ExcAccessToUninitializedField("update_gradients"));
3047  // same as for the scalar case except that we have one more index
3048  const int snc = shape_function_data[shape_function].single_nonzero_component;
3049 
3050  if (snc == -2)
3051  return curl_type ();
3052 
3053  else
3054  switch (dim)
3055  {
3056  case 1:
3057  {
3058  Assert (false, ExcMessage("Computing the curl in 1d is not a useful operation"));
3059  return curl_type ();
3060  }
3061 
3062  case 2:
3063  {
3064  if (snc != -1)
3065  {
3066  curl_type return_value;
3067 
3068  // the single
3069  // nonzero component
3070  // can only be zero
3071  // or one in 2d
3072  if (shape_function_data[shape_function].single_nonzero_component_index == 0)
3073  return_value[0] = -1.0 * fe_values.shape_gradients[snc][q_point][1];
3074  else
3075  return_value[0] = fe_values.shape_gradients[snc][q_point][0];
3076 
3077  return return_value;
3078  }
3079 
3080  else
3081  {
3082  curl_type return_value;
3083 
3084  return_value[0] = 0.0;
3085 
3086  if (shape_function_data[shape_function].is_nonzero_shape_function_component[0])
3087  return_value[0]
3088  -= fe_values.shape_gradients[shape_function_data[shape_function].row_index[0]][q_point][1];
3089 
3090  if (shape_function_data[shape_function].is_nonzero_shape_function_component[1])
3091  return_value[0]
3092  += fe_values.shape_gradients[shape_function_data[shape_function].row_index[1]][q_point][0];
3093 
3094  return return_value;
3095  }
3096  }
3097 
3098  case 3:
3099  {
3100  if (snc != -1)
3101  {
3102  curl_type return_value;
3103 
3104  switch (shape_function_data[shape_function].single_nonzero_component_index)
3105  {
3106  case 0:
3107  {
3108  return_value[0] = 0;
3109  return_value[1] = fe_values.shape_gradients[snc][q_point][2];
3110  return_value[2] = -1.0 * fe_values.shape_gradients[snc][q_point][1];
3111  return return_value;
3112  }
3113 
3114  case 1:
3115  {
3116  return_value[0] = -1.0 * fe_values.shape_gradients[snc][q_point][2];
3117  return_value[1] = 0;
3118  return_value[2] = fe_values.shape_gradients[snc][q_point][0];
3119  return return_value;
3120  }
3121 
3122  default:
3123  {
3124  return_value[0] = fe_values.shape_gradients[snc][q_point][1];
3125  return_value[1] = -1.0 * fe_values.shape_gradients[snc][q_point][0];
3126  return_value[2] = 0;
3127  return return_value;
3128  }
3129  }
3130  }
3131 
3132  else
3133  {
3134  curl_type return_value;
3135 
3136  for (unsigned int i = 0; i < dim; ++i)
3137  return_value[i] = 0.0;
3138 
3139  if (shape_function_data[shape_function].is_nonzero_shape_function_component[0])
3140  {
3141  return_value[1]
3142  += fe_values.shape_gradients[shape_function_data[shape_function].row_index[0]][q_point][2];
3143  return_value[2]
3144  -= fe_values.shape_gradients[shape_function_data[shape_function].row_index[0]][q_point][1];
3145  }
3146 
3147  if (shape_function_data[shape_function].is_nonzero_shape_function_component[1])
3148  {
3149  return_value[0]
3150  -= fe_values.shape_gradients[shape_function_data[shape_function].row_index[1]][q_point][2];
3151  return_value[2]
3152  += fe_values.shape_gradients[shape_function_data[shape_function].row_index[1]][q_point][0];
3153  }
3154 
3155  if (shape_function_data[shape_function].is_nonzero_shape_function_component[2])
3156  {
3157  return_value[0]
3158  += fe_values.shape_gradients[shape_function_data[shape_function].row_index[2]][q_point][1];
3159  return_value[1]
3160  -= fe_values.shape_gradients[shape_function_data[shape_function].row_index[2]][q_point][0];
3161  }
3162 
3163  return return_value;
3164  }
3165  }
3166  }
3167  // should not end up here
3168  Assert (false, ExcInternalError());
3169  return curl_type();
3170  }
3171 
3172  template <int dim, int spacedim>
3173  inline
3175  Vector<dim,spacedim>::hessian (const unsigned int shape_function,
3176  const unsigned int q_point) const
3177  {
3178  // this function works like in
3179  // the case above
3180  typedef FEValuesBase<dim,spacedim> FVB;
3181  Assert (shape_function < fe_values.fe->dofs_per_cell,
3182  ExcIndexRange (shape_function, 0, fe_values.fe->dofs_per_cell));
3183  Assert (fe_values.update_flags & update_hessians,
3184  typename FVB::ExcAccessToUninitializedField("update_hessians"));
3185 
3186  // same as for the scalar case except
3187  // that we have one more index
3188  const int snc = shape_function_data[shape_function].single_nonzero_component;
3189  if (snc == -2)
3190  return hessian_type();
3191  else if (snc != -1)
3192  {
3193  hessian_type return_value;
3194  return_value[shape_function_data[shape_function].single_nonzero_component_index]
3195  = fe_values.shape_hessians[snc][q_point];
3196  return return_value;
3197  }
3198  else
3199  {
3200  hessian_type return_value;
3201  for (unsigned int d=0; d<dim; ++d)
3202  if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
3203  return_value[d]
3204  = fe_values.shape_hessians[shape_function_data[shape_function].row_index[d]][q_point];
3205 
3206  return return_value;
3207  }
3208  }
3209 
3210 
3211  namespace
3212  {
3219  inline
3220  ::SymmetricTensor<2,1>
3221  symmetrize_single_row (const unsigned int n,
3222  const ::Tensor<1,1> &t)
3223  {
3224  Assert (n < 1, ExcIndexRange (n, 0, 1));
3225  (void)n; // removes -Wunused-parameter warning in optimized mode
3226 
3227  const double array[1] = { t[0] };
3228  return ::SymmetricTensor<2,1>(array);
3229  }
3230 
3231 
3232  inline
3233  ::SymmetricTensor<2,2>
3234  symmetrize_single_row (const unsigned int n,
3235  const ::Tensor<1,2> &t)
3236  {
3237  switch (n)
3238  {
3239  case 0:
3240  {
3241  const double array[3] = { t[0], 0, t[1]/2 };
3242  return ::SymmetricTensor<2,2>(array);
3243  }
3244  case 1:
3245  {
3246  const double array[3] = { 0, t[1], t[0]/2 };
3247  return ::SymmetricTensor<2,2>(array);
3248  }
3249  default:
3250  {
3251  Assert (false, ExcIndexRange (n, 0, 2));
3252  return ::SymmetricTensor<2,2>();
3253  }
3254  }
3255  }
3256 
3257 
3258  inline
3259  ::SymmetricTensor<2,3>
3260  symmetrize_single_row (const unsigned int n,
3261  const ::Tensor<1,3> &t)
3262  {
3263  switch (n)
3264  {
3265  case 0:
3266  {
3267  const double array[6] = { t[0], 0, 0, t[1]/2, t[2]/2, 0 };
3268  return ::SymmetricTensor<2,3>(array);
3269  }
3270  case 1:
3271  {
3272  const double array[6] = { 0, t[1], 0, t[0]/2, 0, t[2]/2 };
3273  return ::SymmetricTensor<2,3>(array);
3274  }
3275  case 2:
3276  {
3277  const double array[6] = { 0, 0, t[2], 0, t[0]/2, t[1]/2 };
3278  return ::SymmetricTensor<2,3>(array);
3279  }
3280  default:
3281  {
3282  Assert (false, ExcIndexRange (n, 0, 3));
3283  return ::SymmetricTensor<2,3>();
3284  }
3285  }
3286  }
3287  }
3288 
3289 
3290  template <int dim, int spacedim>
3291  inline
3293  Vector<dim,spacedim>::symmetric_gradient (const unsigned int shape_function,
3294  const unsigned int q_point) const
3295  {
3296  typedef FEValuesBase<dim,spacedim> FVB;
3297  Assert (shape_function < fe_values.fe->dofs_per_cell,
3298  ExcIndexRange (shape_function, 0, fe_values.fe->dofs_per_cell));
3299  Assert (fe_values.update_flags & update_gradients,
3300  typename FVB::ExcAccessToUninitializedField("update_gradients"));
3301 
3302  // same as for the scalar case except
3303  // that we have one more index
3304  const int snc = shape_function_data[shape_function].single_nonzero_component;
3305  if (snc == -2)
3306  return symmetric_gradient_type();
3307  else if (snc != -1)
3308  return symmetrize_single_row (shape_function_data[shape_function].single_nonzero_component_index,
3309  fe_values.shape_gradients[snc][q_point]);
3310  else
3311  {
3312  gradient_type return_value;
3313  for (unsigned int d=0; d<dim; ++d)
3314  if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
3315  return_value[d]
3316  = fe_values.shape_gradients[shape_function_data[shape_function].row_index[d]][q_point];
3317 
3318  return symmetrize(return_value);
3319  }
3320  }
3321 
3322 
3323 
3324  template <int dim, int spacedim>
3325  inline
3327  SymmetricTensor<2, dim, spacedim>::value (const unsigned int shape_function,
3328  const unsigned int q_point) const
3329  {
3330  typedef FEValuesBase<dim,spacedim> FVB;
3331  Assert (shape_function < fe_values.fe->dofs_per_cell,
3332  ExcIndexRange (shape_function, 0, fe_values.fe->dofs_per_cell));
3333  Assert (fe_values.update_flags & update_values,
3334  typename FVB::ExcAccessToUninitializedField("update_values"));
3335 
3336  // similar to the vector case where we
3337  // have more then one index and we need
3338  // to convert between unrolled and
3339  // component indexing for tensors
3340  const int snc
3341  = shape_function_data[shape_function].single_nonzero_component;
3342 
3343  if (snc == -2)
3344  {
3345  // shape function is zero for the
3346  // selected components
3347  return value_type();
3348 
3349  }
3350  else if (snc != -1)
3351  {
3352  value_type return_value;
3353  const unsigned int comp =
3354  shape_function_data[shape_function].single_nonzero_component_index;
3355  return_value[value_type::unrolled_to_component_indices(comp)]
3356  = fe_values.shape_values(snc,q_point);
3357  return return_value;
3358  }
3359  else
3360  {
3361  value_type return_value;
3362  for (unsigned int d = 0; d < value_type::n_independent_components; ++d)
3363  if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
3364  return_value[value_type::unrolled_to_component_indices(d)]
3365  = fe_values.shape_values(shape_function_data[shape_function].row_index[d],q_point);
3366  return return_value;
3367  }
3368  }
3369 
3370 
3371  template <int dim, int spacedim>
3372  inline
3374  SymmetricTensor<2, dim, spacedim>::divergence(const unsigned int shape_function,
3375  const unsigned int q_point) const
3376  {
3377  typedef FEValuesBase<dim,spacedim> FVB;
3378  Assert (shape_function < fe_values.fe->dofs_per_cell,
3379  ExcIndexRange (shape_function, 0, fe_values.fe->dofs_per_cell));
3380  Assert (fe_values.update_flags & update_gradients,
3381  typename FVB::ExcAccessToUninitializedField("update_gradients"));
3382 
3383  const int snc = shape_function_data[shape_function].single_nonzero_component;
3384 
3385  if (snc == -2)
3386  {
3387  // shape function is zero for the
3388  // selected components
3389  return divergence_type();
3390  }
3391  else if (snc != -1)
3392  {
3393  // we have a single non-zero component
3394  // when the symmetric tensor is
3395  // represented in unrolled form.
3396  // this implies we potentially have
3397  // two non-zero components when
3398  // represented in component form! we
3399  // will only have one non-zero entry
3400  // if the non-zero component lies on
3401  // the diagonal of the tensor.
3402  //
3403  // the divergence of a second-order tensor
3404  // is a first order tensor.
3405  //
3406  // assume the second-order tensor is
3407  // A with components A_{ij}. then
3408  // A_{ij} = A_{ji} and there is only
3409  // one (if diagonal) or two non-zero
3410  // entries in the tensorial
3411  // representation. define the
3412  // divergence as:
3413  // b_i := \dfrac{\partial phi_{ij}}{\partial x_j}.
3414  // (which is incidentally also
3415  // b_j := \dfrac{\partial phi_{ij}}{\partial x_i}).
3416  // In both cases, a sum is implied.
3417  //
3418  // Now, we know the nonzero component
3419  // in unrolled form: it is indicated
3420  // by 'snc'. we can figure out which
3421  // tensor components belong to this:
3422  const unsigned int comp =
3423  shape_function_data[shape_function].single_nonzero_component_index;
3424  const unsigned int ii = value_type::unrolled_to_component_indices(comp)[0];
3425  const unsigned int jj = value_type::unrolled_to_component_indices(comp)[1];
3426 
3427  // given the form of the divergence
3428  // above, if ii=jj there is only a
3429  // single nonzero component of the
3430  // full tensor and the gradient
3431  // equals
3432  // b_ii := \dfrac{\partial phi_{ii,ii}}{\partial x_ii}.
3433  // all other entries of 'b' are zero
3434  //
3435  // on the other hand, if ii!=jj, then
3436  // there are two nonzero entries in
3437  // the full tensor and
3438  // b_ii := \dfrac{\partial phi_{ii,jj}}{\partial x_ii}.
3439  // b_jj := \dfrac{\partial phi_{ii,jj}}{\partial x_jj}.
3440  // again, all other entries of 'b' are
3441  // zero
3442  const ::Tensor<1, spacedim> phi_grad = fe_values.shape_gradients[snc][q_point];
3443 
3444  divergence_type return_value;
3445  return_value[ii] = phi_grad[jj];
3446 
3447  if (ii != jj)
3448  return_value[jj] = phi_grad[ii];
3449 
3450  return return_value;
3451 
3452  }
3453  else
3454  {
3455  Assert (false, ExcNotImplemented());
3456  divergence_type return_value;
3457  return return_value;
3458  }
3459  }
3460 
3461  template <int dim, int spacedim>
3462  inline
3464  Tensor<2, dim, spacedim>::value (const unsigned int shape_function,
3465  const unsigned int q_point) const
3466  {
3467  typedef FEValuesBase<dim,spacedim> FVB;
3468  Assert (shape_function < fe_values.fe->dofs_per_cell,
3469  ExcIndexRange (shape_function, 0, fe_values.fe->dofs_per_cell));
3470  Assert (fe_values.update_flags & update_values,
3471  typename FVB::ExcAccessToUninitializedField("update_values"));
3472 
3473  // similar to the vector case where we
3474  // have more then one index and we need
3475  // to convert between unrolled and
3476  // component indexing for tensors
3477  const int snc
3478  = shape_function_data[shape_function].single_nonzero_component;
3479 
3480  if (snc == -2)
3481  {
3482  // shape function is zero for the
3483  // selected components
3484  return value_type();
3485 
3486  }
3487  else if (snc != -1)
3488  {
3489  value_type return_value;
3490  const unsigned int comp =
3491  shape_function_data[shape_function].single_nonzero_component_index;
3493  return_value[indices] = fe_values.shape_values(snc,q_point);
3494  return return_value;
3495  }
3496  else
3497  {
3498  value_type return_value;
3499  for (unsigned int d = 0; d < dim*dim; ++d)
3500  if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
3501  {
3503  return_value[indices]
3504  = fe_values.shape_values(shape_function_data[shape_function].row_index[d],q_point);
3505  }
3506  return return_value;
3507  }
3508  }
3509 
3510 
3511  template <int dim, int spacedim>
3512  inline
3514  Tensor<2, dim, spacedim>::divergence(const unsigned int shape_function,
3515  const unsigned int q_point) const
3516  {
3517  typedef FEValuesBase<dim,spacedim> FVB;
3518  Assert (shape_function < fe_values.fe->dofs_per_cell,
3519  ExcIndexRange (shape_function, 0, fe_values.fe->dofs_per_cell));
3520  Assert (fe_values.update_flags & update_gradients,
3521  typename FVB::ExcAccessToUninitializedField("update_gradients"));
3522 
3523  const int snc = shape_function_data[shape_function].single_nonzero_component;
3524 
3525  if (snc == -2)
3526  {
3527  // shape function is zero for the
3528  // selected components
3529  return divergence_type();
3530  }
3531  else if (snc != -1)
3532  {
3533  // we have a single non-zero component
3534  // when the tensor is
3535  // represented in unrolled form.
3536  //
3537  // the divergence of a second-order tensor
3538  // is a first order tensor.
3539  //
3540  // assume the second-order tensor is
3541  // A with components A_{ij}.
3542  // divergence as:
3543  // b_j := \dfrac{\partial phi_{ij}}{\partial x_i}.
3544  //
3545  // Now, we know the nonzero component
3546  // in unrolled form: it is indicated
3547  // by 'snc'. we can figure out which
3548  // tensor components belong to this:
3549  const unsigned int comp =
3550  shape_function_data[shape_function].single_nonzero_component_index;
3552  const unsigned int ii = indices[0];
3553  const unsigned int jj = indices[1];
3554 
3555  const ::Tensor<1, spacedim> phi_grad = fe_values.shape_gradients[snc][q_point];
3556 
3557  divergence_type return_value;
3558  return_value[jj] = phi_grad[ii];
3559 
3560  return return_value;
3561 
3562  }
3563  else
3564  {
3565  Assert (false, ExcNotImplemented());
3566  divergence_type return_value;
3567  return return_value;
3568  }
3569  }
3570 }
3571 
3572 
3573 
3574 /*------------------------ Inline functions: FEValuesBase ------------------------*/
3575 
3576 
3577 
3578 template <int dim, int spacedim>
3579 inline
3582 operator[] (const FEValuesExtractors::Scalar &scalar) const
3583 {
3584  Assert (scalar.component < fe_values_views_cache.scalars.size(),
3585  ExcIndexRange (scalar.component,
3586  0, fe_values_views_cache.scalars.size()));
3587 
3588  return fe_values_views_cache.scalars[scalar.component];
3589 }
3590 
3591 
3592 
3593 template <int dim, int spacedim>
3594 inline
3597 operator[] (const FEValuesExtractors::Vector &vector) const
3598 {
3599  Assert (vector.first_vector_component <
3600  fe_values_views_cache.vectors.size(),
3602  0, fe_values_views_cache.vectors.size()));
3603 
3604  return fe_values_views_cache.vectors[vector.first_vector_component];
3605 }
3606 
3607 template <int dim, int spacedim>
3608 inline
3612 {
3613  Assert (tensor.first_tensor_component <
3614  fe_values_views_cache.symmetric_second_order_tensors.size(),
3616  0, fe_values_views_cache.symmetric_second_order_tensors.size()));
3617 
3618  return fe_values_views_cache.symmetric_second_order_tensors[tensor.first_tensor_component];
3619 }
3620 
3621 template <int dim, int spacedim>
3622 inline
3625 operator[] (const FEValuesExtractors::Tensor<2> &tensor) const
3626 {
3627  Assert (tensor.first_tensor_component <
3628  fe_values_views_cache.second_order_tensors.size(),
3630  0, fe_values_views_cache.second_order_tensors.size()));
3631 
3632  return fe_values_views_cache.second_order_tensors[tensor.first_tensor_component];
3633 }
3634 
3635 
3636 
3637 
3638 template <int dim, int spacedim>
3639 inline
3640 const double &
3641 FEValuesBase<dim,spacedim>::shape_value (const unsigned int i,
3642  const unsigned int j) const
3643 {
3644  Assert (i < fe->dofs_per_cell,
3645  ExcIndexRange (i, 0, fe->dofs_per_cell));
3647  ExcAccessToUninitializedField("update_values"));
3648  Assert (fe->is_primitive (i),
3649  ExcShapeFunctionNotPrimitive(i));
3650 
3651  // if the entire FE is primitive,
3652  // then we can take a short-cut:
3653  if (fe->is_primitive())
3654  return this->shape_values(i,j);
3655  else
3656  {
3657  // otherwise, use the mapping
3658  // between shape function
3659  // numbers and rows. note that
3660  // by the assertions above, we
3661  // know that this particular
3662  // shape function is primitive,
3663  // so we can call
3664  // system_to_component_index
3665  const unsigned int
3666  row = this->shape_function_to_row_table[i * fe->n_components() + fe->system_to_component_index(i).first];
3667  return this->shape_values(row, j);
3668  }
3669 }
3670 
3671 
3672 
3673 template <int dim, int spacedim>
3674 inline
3675 double
3677  const unsigned int j,
3678  const unsigned int component) const
3679 {
3680  Assert (i < fe->dofs_per_cell,
3681  ExcIndexRange (i, 0, fe->dofs_per_cell));
3682  Assert (this->update_flags & update_values,
3683  ExcAccessToUninitializedField("update_values"));
3684  Assert (component < fe->n_components(),
3685  ExcIndexRange(component, 0, fe->n_components()));
3686 
3687  // check whether the shape function
3688  // is non-zero at all within
3689  // this component:
3690  if (fe->get_nonzero_components(i)[component] == false)
3691  return 0;
3692 
3693  // look up the right row in the
3694  // table and take the data from
3695  // there
3696  const unsigned int
3697  row = this->shape_function_to_row_table[i * fe->n_components() + component];
3698  return this->shape_values(row, j);
3699 }
3700 
3701 
3702 
3703 template <int dim, int spacedim>
3704 inline
3705 const Tensor<1,spacedim> &
3706 FEValuesBase<dim,spacedim>::shape_grad (const unsigned int i,
3707  const unsigned int j) const
3708 {
3709  Assert (i < fe->dofs_per_cell,
3710  ExcIndexRange (i, 0, fe->dofs_per_cell));
3712  ExcAccessToUninitializedField("update_gradients"));
3713  Assert (fe->is_primitive (i),
3714  ExcShapeFunctionNotPrimitive(i));
3715  Assert (i<this->shape_gradients.size(),
3716  ExcIndexRange (i, 0, this->shape_gradients.size()));
3717  Assert (j<this->shape_gradients[0].size(),
3718  ExcIndexRange (j, 0, this->shape_gradients[0].size()));
3719 
3720  // if the entire FE is primitive,
3721  // then we can take a short-cut:
3722  if (fe->is_primitive())
3723  return this->shape_gradients[i][j];
3724  else
3725  {
3726  // otherwise, use the mapping
3727  // between shape function
3728  // numbers and rows. note that
3729  // by the assertions above, we
3730  // know that this particular
3731  // shape function is primitive,
3732  // so we can call
3733  // system_to_component_index
3734  const unsigned int
3735  row = this->shape_function_to_row_table[i * fe->n_components() + fe->system_to_component_index(i).first];
3736  return this->shape_gradients[row][j];
3737  }
3738 }
3739 
3740 
3741 
3742 template <int dim, int spacedim>
3743 inline
3746  const unsigned int j,
3747  const unsigned int component) const
3748 {
3749  Assert (i < fe->dofs_per_cell,
3750  ExcIndexRange (i, 0, fe->dofs_per_cell));
3751  Assert (this->update_flags & update_gradients,
3752  ExcAccessToUninitializedField("update_gradients"));
3753  Assert (component < fe->n_components(),
3754  ExcIndexRange(component, 0, fe->n_components()));
3755 
3756  // check whether the shape function
3757  // is non-zero at all within
3758  // this component:
3759  if (fe->get_nonzero_components(i)[component] == false)
3760  return Tensor<1,spacedim>();
3761 
3762  // look up the right row in the
3763  // table and take the data from
3764  // there
3765  const unsigned int
3766  row = this->shape_function_to_row_table[i * fe->n_components() + component];
3767  return this->shape_gradients[row][j];
3768 }
3769 
3770 
3771 
3772 template <int dim, int spacedim>
3773 inline
3774 const Tensor<2,spacedim> &
3775 FEValuesBase<dim,spacedim>::shape_hessian (const unsigned int i,
3776  const unsigned int j) const
3777 {
3778  Assert (i < fe->dofs_per_cell,
3779  ExcIndexRange (i, 0, fe->dofs_per_cell));
3781  ExcAccessToUninitializedField("update_hessians"));
3782  Assert (fe->is_primitive (i),
3783  ExcShapeFunctionNotPrimitive(i));
3784  Assert (i<this->shape_hessians.size(),
3785  ExcIndexRange (i, 0, this->shape_hessians.size()));
3786  Assert (j<this->shape_hessians[0].size(),
3787  ExcIndexRange (j, 0, this->shape_hessians[0].size()));
3788 
3789  // if the entire FE is primitive,
3790  // then we can take a short-cut:
3791  if (fe->is_primitive())
3792  return this->shape_hessians[i][j];
3793  else
3794  {
3795  // otherwise, use the mapping
3796  // between shape function
3797  // numbers and rows. note that
3798  // by the assertions above, we
3799  // know that this particular
3800  // shape function is primitive,
3801  // so we can call
3802  // system_to_component_index
3803  const unsigned int
3804  row = this->shape_function_to_row_table[i * fe->n_components() + fe->system_to_component_index(i).first];
3805  return this->shape_hessians[row][j];
3806  }
3807 }
3808 
3809 
3810 
3811 template <int dim, int spacedim>
3812 inline
3813 const Tensor<2,spacedim> &
3815  const unsigned int j) const
3816 {
3817  return shape_hessian(i,j);
3818 }
3819 
3820 
3821 
3822 template <int dim, int spacedim>
3823 inline
3826  const unsigned int j,
3827  const unsigned int component) const
3828 {
3829  Assert (i < fe->dofs_per_cell,
3830  ExcIndexRange (i, 0, fe->dofs_per_cell));
3831  Assert (this->update_flags & update_hessians,
3832  ExcAccessToUninitializedField("update_hessians"));
3833  Assert (component < fe->n_components(),
3834  ExcIndexRange(component, 0, fe->n_components()));
3835 
3836  // check whether the shape function
3837  // is non-zero at all within
3838  // this component:
3839  if (fe->get_nonzero_components(i)[component] == false)
3840  return Tensor<2,spacedim>();
3841 
3842  // look up the right row in the
3843  // table and take the data from
3844  // there
3845  const unsigned int
3846  row = this->shape_function_to_row_table[i * fe->n_components() + component];
3847  return this->shape_hessians[row][j];
3848 }
3849 
3850 
3851 
3852 template <int dim, int spacedim>
3853 inline
3856  const unsigned int j,
3857  const unsigned int component) const
3858 {
3859  return shape_hessian_component(i,j,component);
3860 }
3861 
3862 
3863 
3864 template <int dim, int spacedim>
3865 inline
3868 {
3869  return *fe;
3870 }
3871 
3872 
3873 template <int dim, int spacedim>
3874 inline
3875 const Mapping<dim,spacedim> &
3877 {
3878  return *mapping;
3879 }
3880 
3881 
3882 
3883 template <int dim, int spacedim>
3884 inline
3887 {
3888  return this->update_flags;
3889 }
3890 
3891 
3892 
3893 template <int dim, int spacedim>
3894 inline
3895 const std::vector<Point<spacedim> > &
3897 {
3899  ExcAccessToUninitializedField("update_quadrature_points"));
3900  return this->quadrature_points;
3901 }
3902 
3903 
3904 
3905 template <int dim, int spacedim>
3906 inline
3907 const std::vector<double> &
3909 {
3911  ExcAccessToUninitializedField("update_JxW_values"));
3912  return this->JxW_values;
3913 }
3914 
3915 
3916 
3917 template <int dim, int spacedim>
3918 inline
3919 const std::vector<DerivativeForm<1,dim,spacedim> > &
3921 {
3923  ExcAccessToUninitializedField("update_jacobians"));
3924  return this->jacobians;
3925 }
3926 
3927 
3928 
3929 template <int dim, int spacedim>
3930 inline
3931 const std::vector<DerivativeForm<2,dim,spacedim> > &
3933 {
3935  ExcAccessToUninitializedField("update_jacobians_grads"));
3936  return this->jacobian_grads;
3937 }
3938 
3939 
3940 
3941 template <int dim, int spacedim>
3942 inline
3943 const std::vector<DerivativeForm<1,spacedim,dim> > &
3945 {
3947  ExcAccessToUninitializedField("update_inverse_jacobians"));
3948  return this->inverse_jacobians;
3949 }
3950 
3951 
3952 
3953 template <int dim, int spacedim>
3954 inline
3955 const Point<spacedim> &
3956 FEValuesBase<dim,spacedim>::quadrature_point (const unsigned int i) const
3957 {
3959  ExcAccessToUninitializedField("update_quadrature_points"));
3960  Assert (i<this->quadrature_points.size(), ExcIndexRange(i, 0, this->quadrature_points.size()));
3961 
3962  return this->quadrature_points[i];
3963 }
3964 
3965 
3966 
3967 
3968 template <int dim, int spacedim>
3969 inline
3970 double
3971 FEValuesBase<dim,spacedim>::JxW (const unsigned int i) const
3972 {
3974  ExcAccessToUninitializedField("update_JxW_values"));
3975  Assert (i<this->JxW_values.size(), ExcIndexRange(i, 0, this->JxW_values.size()));
3976 
3977  return this->JxW_values[i];
3978 }
3979 
3980 
3981 
3982 template <int dim, int spacedim>
3983 inline
3985 FEValuesBase<dim,spacedim>::jacobian (const unsigned int i) const
3986 {
3988  ExcAccessToUninitializedField("update_jacobians"));
3989  Assert (i<this->jacobians.size(), ExcIndexRange(i, 0, this->jacobians.size()));
3990 
3991  return this->jacobians[i];
3992 }
3993 
3994 
3995 
3996 template <int dim, int spacedim>
3997 inline
3999 FEValuesBase<dim,spacedim>::jacobian_grad (const unsigned int i) const
4000 {
4002  ExcAccessToUninitializedField("update_jacobians_grads"));
4003  Assert (i<this->jacobian_grads.size(), ExcIndexRange(i, 0, this->jacobian_grads.size()));
4004 
4005  return this->jacobian_grads[i];
4006 }
4007 
4008 
4009 
4010 template <int dim, int spacedim>
4011 inline
4013 FEValuesBase<dim,spacedim>::inverse_jacobian (const unsigned int i) const
4014 {
4016  ExcAccessToUninitializedField("update_inverse_jacobians"));
4017  Assert (i<this->inverse_jacobians.size(), ExcIndexRange(i, 0, this->inverse_jacobians.size()));
4018 
4019  return this->inverse_jacobians[i];
4020 }
4021 
4022 
4023 template <int dim, int spacedim>
4024 template <class InputVector>
4025 inline
4026 void
4027 FEValuesBase<dim,spacedim>::get_function_grads (const InputVector &fe_function,
4028  std::vector<Tensor<1,spacedim> > &gradients) const
4029 {
4030  get_function_gradients(fe_function, gradients);
4031 }
4032 
4033 
4034 
4035 template <int dim, int spacedim>
4036 template <class InputVector>
4037 inline
4038 void
4040  const InputVector &fe_function,
4041  const VectorSlice<const std::vector<types::global_dof_index> > &indices,
4042  std::vector<Tensor<1,spacedim> > &values) const
4043 {
4044  get_function_gradients(fe_function, indices, values);
4045 }
4046 
4047 
4048 
4049 template <int dim, int spacedim>
4050 template <class InputVector>
4051 inline
4052 void
4054 get_function_grads (const InputVector &fe_function,
4055  std::vector<std::vector<Tensor<1,spacedim> > > &gradients) const
4056 {
4057  get_function_gradients(fe_function, gradients);
4058 }
4059 
4060 
4061 
4062 template <int dim, int spacedim>
4063 template <class InputVector>
4064 inline
4065 void
4067  const InputVector &fe_function,
4068  const VectorSlice<const std::vector<types::global_dof_index> > &indices,
4069  std::vector<std::vector<Tensor<1,spacedim> > > &values,
4070  bool q_points_fastest) const
4071 {
4072  get_function_gradients(fe_function, indices, values, q_points_fastest);
4073 }
4074 
4075 
4076 
4077 template <int dim, int spacedim>
4078 template <class InputVector>
4079 inline
4080 void
4082 get_function_2nd_derivatives (const InputVector &fe_function,
4083  std::vector<Tensor<2,spacedim> > &hessians) const
4084 {
4085  get_function_hessians(fe_function, hessians);
4086 }
4087 
4088 
4089 
4090 template <int dim, int spacedim>
4091 template <class InputVector>
4092 inline
4093 void
4095 get_function_2nd_derivatives (const InputVector &fe_function,
4096  std::vector<std::vector<Tensor<2,spacedim> > > &hessians,
4097  bool quadrature_points_fastest) const
4098 {
4099  get_function_hessians(fe_function, hessians, quadrature_points_fastest);
4100 }
4101 
4102 
4103 
4104 template <int dim, int spacedim>
4105 inline
4106 const Point<spacedim> &
4107 FEValuesBase<dim,spacedim>::normal_vector (const unsigned int i) const
4108 {
4109  typedef FEValuesBase<dim,spacedim> FVB;
4111  typename FVB::ExcAccessToUninitializedField("update_normal_vectors"));
4112  Assert (i<this->normal_vectors.size(),
4113  ExcIndexRange(i, 0, this->normal_vectors.size()));
4114 
4115  return this->normal_vectors[i];
4116 }
4117 
4118 
4119 
4120 template <int dim, int spacedim>
4121 inline
4122 const Point<spacedim> &
4123 FEValuesBase<dim,spacedim>::cell_normal_vector (const unsigned int i) const
4124 {
4125  return this->normal_vector(i);
4126 }
4127 
4128 
4129 
4130 
4131 /*------------------------ Inline functions: FEValues ----------------------------*/
4132 
4133 
4134 template <int dim, int spacedim>
4135 inline
4136 const Quadrature<dim> &
4138 {
4139  return quadrature;
4140 }
4141 
4142 
4143 
4144 template <int dim, int spacedim>
4145 inline
4146 const FEValues<dim,spacedim> &
4148 {
4149  return *this;
4150 }
4151 
4152 
4153 /*------------------------ Inline functions: FEFaceValuesBase --------------------*/
4154 
4155 
4156 template <int dim, int spacedim>
4157 inline
4158 unsigned int
4160 {
4161  return present_face_index;
4162 }
4163 
4164 
4165 /*------------------------ Inline functions: FE*FaceValues --------------------*/
4166 
4167 template <int dim, int spacedim>
4168 inline
4169 const Quadrature<dim-1> &
4171 {
4172  return quadrature;
4173 }
4174 
4175 
4176 
4177 template <int dim, int spacedim>
4178 inline
4181 {
4182  return *this;
4183 }
4184 
4185 
4186 
4187 template <int dim, int spacedim>
4188 inline
4191 {
4192  return *this;
4193 }
4194 
4195 
4196 
4197 template <int dim, int spacedim>
4198 inline
4199 const Tensor<1,spacedim> &
4200 FEFaceValuesBase<dim,spacedim>::boundary_form (const unsigned int i) const
4201 {
4202  typedef FEValuesBase<dim,spacedim> FVB;
4203  Assert (i<this->boundary_forms.size(),
4204  ExcIndexRange(i, 0, this->boundary_forms.size()));
4206  typename FVB::ExcAccessToUninitializedField("update_boundary_forms"));
4207 
4208  return this->boundary_forms[i];
4209 }
4210 
4211 #endif // DOXYGEN
4212 
4213 DEAL_II_NAMESPACE_CLOSE
4214 
4215 #endif
Transformed quadrature weights.
void get_function_symmetric_gradients(const InputVector &fe_function, std::vector< symmetric_gradient_type > &symmetric_gradients) const
const std::vector< DerivativeForm< 2, dim, spacedim > > & get_jacobian_grads() const
HessianVector shape_hessians
Definition: fe_values.h:1185
Shape function values.
double shape_value_component(const unsigned int function_no, const unsigned int point_no, const unsigned int component) const
const Tensor< 1, spacedim > & shape_grad(const unsigned int function_no, const unsigned int quadrature_point) const
hessian_type hessian(const unsigned int shape_function, const unsigned int q_point) const
Number value_type
Definition: vector.h:119
const FEValuesViews::Scalar< dim, spacedim > & operator[](const FEValuesExtractors::Scalar &scalar) const
CellSimilarity::Similarity cell_similarity
Definition: fe_values.h:2377
::Tensor< 1, spacedim > divergence_type
Definition: fe_values.h:919
Vector & operator=(const Vector< dim, spacedim > &)
Cache(const FEValuesBase< dim, spacedim > &fe_values)
void get_function_grads(const InputVector &fe_function, std::vector< Tensor< 1, spacedim > > &gradients) const DEAL_II_DEPRECATED
const FEValuesBase< dim, spacedim > & fe_values
Definition: fe_values.h:1056
Tensor< 1, spacedim > shape_grad_component(const unsigned int function_no, const unsigned int point_no, const unsigned int component) const
unsigned int present_face_index
Definition: fe_values.h:2598
const Quadrature< dim-1 > & get_quadrature() const
const Tensor< 2, spacedim > & shape_2nd_derivative(const unsigned int function_no, const unsigned int point_no) const DEAL_II_DEPRECATED
std::vector< std::vector< Tensor< 2, spacedim > > > HessianVector
Definition: fe_values.h:1165
const Tensor< 2, spacedim > & shape_hessian(const unsigned int function_no, const unsigned int point_no) const
std::vector< ShapeFunctionData > shape_function_data
Definition: fe_values.h:881
void get_function_hessians(const InputVector &fe_function, std::vector< hessian_type > &hessians) const
void get_function_gradients(const InputVector &fe_function, std::vector< gradient_type > &gradients) const
curl_type curl(const unsigned int shape_function, const unsigned int q_point) const
std::vector<::FEValuesViews::Scalar< dim, spacedim > > scalars
Definition: fe_values.h:1091
static TableIndices< rank_ > unrolled_to_component_indices(const unsigned int i)
MappingType
Definition: mapping.h:58
const unsigned int dofs_per_cell
Definition: fe_values.h:1418
GradientVector shape_gradients
Definition: fe_values.h:1178
const Point< spacedim > & normal_vector(const unsigned int i) const
const DerivativeForm< 1, dim, spacedim > & jacobian(const unsigned int quadrature_point) const
const unsigned int component
Definition: fe_values.h:318
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)
::Tensor< 3, spacedim > hessian_type
Definition: fe_values.h:408
const Quadrature< dim-1 > quadrature
Definition: fe_values.h:2603
Volume element.
UpdateFlags update_flags
Definition: fe_values.h:1277
::ExceptionBase & ExcMessage(std::string arg1)
unsigned int get_face_index() const
const FESubfaceValues< dim, spacedim > & get_present_fe_values() const
void get_function_laplacians(const InputVector &fe_function, std::vector< value_type > &laplacians) const
const std::vector< DerivativeForm< 1, spacedim, dim > > & get_inverse_jacobians() const
::internal::CurlType< spacedim >::type curl_type
Definition: fe_values.h:401
Outer normal vector, not normalized.
void get_function_curls(const InputVector &fe_function, std::vector< curl_type > &curls) const
const FEFaceValues< dim, spacedim > & get_present_fe_values() const
hessian_type hessian(const unsigned int shape_function, const unsigned int q_point) const
STL namespace.
Transformed quadrature points.
const Tensor< 1, spacedim > & boundary_form(const unsigned int i) const
value_type value(const unsigned int shape_function, const unsigned int q_point) const
::Tensor< 1, spacedim > gradient_type
Definition: fe_values.h:157
::SymmetricTensor< 2, spacedim > value_type
Definition: fe_values.h:724
unsigned int row_index[spacedim]
Definition: fe_values.h:435
void transform(std::vector< Tensor< 1, spacedim > > &transformed, const std::vector< Tensor< 1, dim > > &original, MappingType mapping) const
void check_cell_similarity(const typename Triangulation< dim, spacedim >::cell_iterator &cell)
static const unsigned int space_dimension
Definition: fe_values.h:1406
::internal::FEValuesViews::Cache< dim, spacedim > fe_values_views_cache
Definition: fe_values.h:2403
std::vector< Point< spacedim > > normal_vectors
Definition: fe_values.h:1229
std::vector< DerivativeForm< 1, dim, spacedim > > jacobians
Definition: fe_values.h:1205
const std::vector< double > & get_JxW_values() const
void get_function_divergences(const InputVector &fe_function, std::vector< divergence_type > &divergences) const
UpdateFlags get_update_flags() const
const Point< spacedim > & quadrature_point(const unsigned int i) const
void get_function_laplacians(const InputVector &fe_function, std::vector< value_type > &laplacians) const
void get_function_hessians(const InputVector &fe_function, std::vector< Tensor< 2, spacedim > > &hessians) const
void get_function_gradients(const InputVector &fe_function, std::vector< gradient_type > &gradients) const
std::vector< double > JxW_values
Definition: fe_values.h:1200
const Quadrature< dim > & get_quadrature() const
UpdateFlags compute_update_flags(const UpdateFlags update_flags) const
void get_function_values(const InputVector &fe_function, std::vector< number > &values) const
const SmartPointer< const FiniteElement< dim, spacedim >, FEValuesBase< dim, spacedim > > fe
Definition: fe_values.h:2349
const std::vector< DerivativeForm< 1, dim, spacedim > > & get_jacobians() const
const std::vector< Point< spacedim > > & get_quadrature_points() const
std::vector< Tensor< 1, spacedim > > boundary_forms
Definition: fe_values.h:1235
symmetric_gradient_type symmetric_gradient(const unsigned int shape_function, const unsigned int q_point) const
divergence_type divergence(const unsigned int shape_function, const unsigned int q_point) const
::Tensor< 1, spacedim > value_type
Definition: fe_values.h:364
const std::vector< Point< spacedim > > & get_cell_normal_vectors() const DEAL_II_DEPRECATED
boost::signals2::connection tria_listener
Definition: fe_values.h:2320
Definition: types.h:29
std::vector< ShapeFunctionData > shape_function_data
Definition: fe_values.h:1067
std::vector< DerivativeForm< 1, spacedim, dim > > inverse_jacobians
Definition: fe_values.h:1216
#define Assert(cond, exc)
Definition: exceptions.h:299
::Tensor< 2, spacedim > value_type
Definition: fe_values.h:914
UpdateFlags
void get_function_laplacians(const InputVector &fe_function, std::vector< number > &laplacians) const
gradient_type gradient(const unsigned int shape_function, const unsigned int q_point) const
std::vector< ShapeFunctionData > shape_function_data
Definition: fe_values.h:685
const Quadrature< dim > quadrature
Definition: fe_values.h:2507
const unsigned int first_vector_component
Definition: fe_values.h:680
std::size_t memory_consumption() const
ShapeVector shape_values
Definition: fe_values.h:1171
std::vector< unsigned int > shape_function_to_row_table
Definition: fe_values.h:1272
void get_function_gradients(const InputVector &fe_function, std::vector< Tensor< 1, spacedim > > &gradients) const
void invalidate_present_cell()
const SmartPointer< const Mapping< dim, spacedim >, FEValuesBase< dim, spacedim > > mapping
Definition: fe_values.h:2344
::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
BlockCompressedSparsityPattern CompressedBlockSparsityPattern DEAL_II_DEPRECATED
gradient_type gradient(const unsigned int shape_function, const unsigned int q_point) const
bool is_nonzero_shape_function_component[spacedim]
Definition: fe_values.h:424
const Triangulation< dim, spacedim >::cell_iterator get_cell() const
Tensor< 2, spacedim > shape_hessian_component(const unsigned int function_no, const unsigned int point_no, const unsigned int component) const
Second derivatives of shape functions.
Gradient of volume element.
const Mapping< dim, spacedim > & get_mapping() const
std::vector< DerivativeForm< 2, dim, spacedim > > jacobian_grads
Definition: fe_values.h:1211
std::vector< std::vector< Tensor< 1, spacedim > > > GradientVector
Definition: fe_values.h:1160
DeclException0(ExcCannotInitializeField)
const DerivativeForm< 1, spacedim, dim > & inverse_jacobian(const unsigned int quadrature_point) const
std::auto_ptr< const CellIteratorBase > present_cell
Definition: fe_values.h:2304
const unsigned int n_quadrature_points
Definition: fe_values.h:1411
const Point< spacedim > & cell_normal_vector(const unsigned int i) const DEAL_II_DEPRECATED
const double & shape_value(const unsigned int function_no, const unsigned int point_no) const
void get_function_values(const InputVector &fe_function, std::vector< value_type > &values) const
void initialize(const unsigned int n_quadrature_points, const FiniteElement< dim, spacedim > &fe, const UpdateFlags flags)
SmartPointer< typename Mapping< dim, spacedim >::InternalDataBase, FEValuesBase< dim, spacedim > > fe_data
Definition: fe_values.h:2360
static const unsigned int dimension
Definition: fe_values.h:1401
::Tensor< 2, spacedim > hessian_type
Definition: fe_values.h:164
const FEValuesBase< dim, spacedim > & fe_values
Definition: fe_values.h:312
Scalar & operator=(const Scalar< dim, spacedim > &)
SmartPointer< typename Mapping< dim, spacedim >::InternalDataBase, FEValuesBase< dim, spacedim > > mapping_data
Definition: fe_values.h:2355
::SymmetricTensor< 2, spacedim > symmetric_gradient_type
Definition: fe_values.h:386
std::vector< Point< spacedim > > quadrature_points
Definition: fe_values.h:1223
Shape function gradients.
Normal vectors.
Table< 2, double > ShapeVector
Definition: fe_values.h:1154
void maybe_invalidate_previous_present_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell)
::Tensor< 2, spacedim > gradient_type
Definition: fe_values.h:374
const FEValuesBase< dim, spacedim > & fe_values
Definition: fe_values.h:674
const DerivativeForm< 2, dim, spacedim > & jacobian_grad(const unsigned int quadrature_point) const
value_type value(const unsigned int shape_function, const unsigned int q_point) const
::ExceptionBase & ExcNotImplemented()
const FiniteElement< dim, spacedim > & get_fe() const
Tensor< 2, spacedim > shape_2nd_derivative_component(const unsigned int function_no, const unsigned int point_no, const unsigned int component) const DEAL_II_DEPRECATED
::ExceptionBase & ExcInternalError()
void get_function_2nd_derivatives(const InputVector &, std::vector< Tensor< 2, spacedim > > &) const DEAL_II_DEPRECATED
void get_function_values(const InputVector &fe_function, std::vector< value_type > &values) const
void get_function_hessians(const InputVector &fe_function, std::vector< hessian_type > &hessians) const
double JxW(const unsigned int quadrature_point) const
const std::vector< Point< spacedim > > & get_normal_vectors() const
DeclException1(ExcAccessToUninitializedField, char *,<< ("You are requesting information from an FEValues/FEFaceValues/FESubfaceValues ""object for which this kind of information has not been computed. What ""information these objects compute is determined by the update_* flags you ""pass to the constructor. Here, the operation you are attempting requires ""the <")<< arg1<< "> flag to be set, but it was apparently not specified upon construction.")
const FEValues< dim, spacedim > & get_present_fe_values() const
CellSimilarity::Similarity get_cell_similarity() const
std::vector< ShapeFunctionData > shape_function_data
Definition: fe_values.h:323
const FEValuesBase< dim, spacedim > & fe_values
Definition: fe_values.h:870