Reference documentation for deal.II version 8.5.0
fe_values.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1998 - 2017 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 at
12 // the top level of the deal.II distribution.
13 //
14 // ---------------------------------------------------------------------
15 
16 #ifndef dealii__fe_values_h
17 #define dealii__fe_values_h
18 
19 
20 #include <deal.II/base/config.h>
21 #include <deal.II/base/exceptions.h>
22 #include <deal.II/base/subscriptor.h>
23 #include <deal.II/base/point.h>
24 #include <deal.II/base/derivative_form.h>
25 #include <deal.II/base/symmetric_tensor.h>
26 #include <deal.II/base/vector_slice.h>
27 #include <deal.II/base/quadrature.h>
28 #include <deal.II/base/table.h>
29 #include <deal.II/base/std_cxx11/unique_ptr.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 
40 #include <algorithm>
41 
42 // dummy include in order to have the
43 // definition of PetscScalar available
44 // without including other PETSc stuff
45 #ifdef DEAL_II_WITH_PETSC
46 # include <petsc.h>
47 #endif
48 
49 DEAL_II_NAMESPACE_OPEN
50 
51 template <int dim, int spacedim=dim> class FEValuesBase;
52 
53 namespace internal
54 {
59  template <int dim>
60  struct CurlType;
61 
68  template <>
69  struct CurlType<1>
70  {
71  typedef Tensor<1,1> type;
72  };
73 
80  template <>
81  struct CurlType<2>
82  {
83  typedef Tensor<1,1> type;
84  };
85 
92  template <>
93  struct CurlType<3>
94  {
95  typedef Tensor<1,3> type;
96  };
97 }
98 
99 
100 
101 
123 namespace FEValuesViews
124 {
136  template <int dim, int spacedim=dim>
137  class Scalar
138  {
139  public:
145  typedef double value_type;
146 
152  typedef ::Tensor<1,spacedim> gradient_type;
153 
159  typedef ::Tensor<2,spacedim> hessian_type;
160 
166  typedef ::Tensor<3,spacedim> third_derivative_type;
167 
173  {
183 
192  unsigned int row_index;
193  };
194 
198  Scalar ();
199 
205  Scalar (const FEValuesBase<dim,spacedim> &fe_values_base,
206  const unsigned int component);
207 
213 
227  value_type
228  value (const unsigned int shape_function,
229  const unsigned int q_point) const;
230 
241  gradient_type
242  gradient (const unsigned int shape_function,
243  const unsigned int q_point) const;
244 
255  hessian_type
256  hessian (const unsigned int shape_function,
257  const unsigned int q_point) const;
258 
269  third_derivative_type
270  third_derivative (const unsigned int shape_function,
271  const unsigned int q_point) const;
272 
290  template <class InputVector>
291  void get_function_values (const InputVector &fe_function,
292  std::vector<typename ProductType<value_type,typename InputVector::value_type>::type> &values) const;
293 
311  template <class InputVector>
312  void get_function_gradients (const InputVector &fe_function,
313  std::vector<typename ProductType<gradient_type,typename InputVector::value_type>::type> &gradients) const;
314 
332  template <class InputVector>
333  void get_function_hessians (const InputVector &fe_function,
334  std::vector<typename ProductType<hessian_type,typename InputVector::value_type>::type> &hessians) const;
335 
354  template <class InputVector>
355  void get_function_laplacians (const InputVector &fe_function,
356  std::vector<typename ProductType<value_type,typename InputVector::value_type>::type> &laplacians) const;
357 
376  template <class InputVector>
377  void get_function_third_derivatives (const InputVector &fe_function,
378  std::vector<typename ProductType<third_derivative_type,
379  typename InputVector::value_type>::type> &third_derivatives) const;
380 
381  private:
386 
391  const unsigned int component;
392 
396  std::vector<ShapeFunctionData> shape_function_data;
397  };
398 
399 
400 
430  template <int dim, int spacedim=dim>
431  class Vector
432  {
433  public:
439  typedef ::Tensor<1,spacedim> value_type;
440 
449  typedef ::Tensor<2,spacedim> gradient_type;
450 
461  typedef ::SymmetricTensor<2,spacedim> symmetric_gradient_type;
462 
468  typedef double divergence_type;
469 
476  typedef typename ::internal::CurlType<spacedim>::type curl_type;
477 
483  typedef ::Tensor<3,spacedim> hessian_type;
484 
490  typedef ::Tensor<4,spacedim> third_derivative_type;
491 
497  {
507 
517  unsigned int row_index[spacedim];
518 
528  unsigned int single_nonzero_component_index;
529  };
530 
534  Vector ();
535 
544  Vector (const FEValuesBase<dim,spacedim> &fe_values_base,
545  const unsigned int first_vector_component);
546 
552 
569  value_type
570  value (const unsigned int shape_function,
571  const unsigned int q_point) const;
572 
586  gradient_type
587  gradient (const unsigned int shape_function,
588  const unsigned int q_point) const;
589 
605  symmetric_gradient_type
606  symmetric_gradient (const unsigned int shape_function,
607  const unsigned int q_point) const;
608 
619  divergence_type
620  divergence (const unsigned int shape_function,
621  const unsigned int q_point) const;
622 
643  curl_type
644  curl (const unsigned int shape_function,
645  const unsigned int q_point) const;
646 
657  hessian_type
658  hessian (const unsigned int shape_function,
659  const unsigned int q_point) const;
660 
671  third_derivative_type
672  third_derivative (const unsigned int shape_function,
673  const unsigned int q_point) const;
674 
692  template <class InputVector>
693  void get_function_values (const InputVector &fe_function,
694  std::vector<typename ProductType<value_type,typename InputVector::value_type>::type> &values) const;
695 
713  template <class InputVector>
714  void get_function_gradients (const InputVector &fe_function,
715  std::vector<typename ProductType<gradient_type,typename InputVector::value_type>::type> &gradients) const;
716 
740  template <class InputVector>
741  void
742  get_function_symmetric_gradients (const InputVector &fe_function,
743  std::vector<typename ProductType<symmetric_gradient_type,typename InputVector::value_type>::type> &symmetric_gradients) const;
744 
763  template <class InputVector>
764  void get_function_divergences (const InputVector &fe_function,
765  std::vector<typename ProductType<divergence_type,typename InputVector::value_type>::type> &divergences) const;
766 
785  template <class InputVector>
786  void get_function_curls (const InputVector &fe_function,
787  std::vector<typename ProductType<curl_type,typename InputVector::value_type>::type> &curls) const;
788 
806  template <class InputVector>
807  void get_function_hessians (const InputVector &fe_function,
808  std::vector<typename ProductType<hessian_type,typename InputVector::value_type>::type> &hessians) const;
809 
828  template <class InputVector>
829  void get_function_laplacians (const InputVector &fe_function,
830  std::vector<typename ProductType<value_type,typename InputVector::value_type>::type> &laplacians) const;
831 
850  template <class InputVector>
851  void get_function_third_derivatives (const InputVector &fe_function,
852  std::vector<typename ProductType<third_derivative_type,
853  typename InputVector::value_type>::type> &third_derivatives) const;
854 
855  private:
860 
865  const unsigned int first_vector_component;
866 
870  std::vector<ShapeFunctionData> shape_function_data;
871  };
872 
873 
874  template <int rank, int dim, int spacedim = dim>
875  class SymmetricTensor;
876 
901  template <int dim, int spacedim>
902  class SymmetricTensor<2,dim,spacedim>
903  {
904  public:
911  typedef ::SymmetricTensor<2, spacedim> value_type;
912 
922  typedef ::Tensor<1, spacedim> divergence_type;
923 
928  struct ShapeFunctionData
929  {
938  bool is_nonzero_shape_function_component[value_type::n_independent_components];
939 
949  unsigned int row_index[value_type::n_independent_components];
950 
960  unsigned int single_nonzero_component_index;
961  };
962 
966  SymmetricTensor();
967 
977  SymmetricTensor(const FEValuesBase<dim, spacedim> &fe_values_base,
978  const unsigned int first_tensor_component);
979 
984  SymmetricTensor &operator=(const SymmetricTensor<2, dim, spacedim> &);
985 
1003  value_type
1004  value (const unsigned int shape_function,
1005  const unsigned int q_point) const;
1006 
1007 
1021  divergence_type
1022  divergence (const unsigned int shape_function,
1023  const unsigned int q_point) const;
1024 
1042  template <class InputVector>
1043  void get_function_values (const InputVector &fe_function,
1044  std::vector<typename ProductType<value_type,typename InputVector::value_type>::type> &values) const;
1045 
1067  template <class InputVector>
1068  void get_function_divergences (const InputVector &fe_function,
1069  std::vector<typename ProductType<divergence_type,typename InputVector::value_type>::type> &divergences) const;
1070 
1071  private:
1076 
1081  const unsigned int first_tensor_component;
1082 
1086  std::vector<ShapeFunctionData> shape_function_data;
1087  };
1088 
1089 
1090  template <int rank, int dim, int spacedim = dim>
1091  class Tensor;
1092 
1112  template <int dim, int spacedim>
1113  class Tensor<2,dim,spacedim>
1114  {
1115  public:
1116 
1121  typedef ::Tensor<2, spacedim> value_type;
1122 
1126  typedef ::Tensor<1, spacedim> divergence_type;
1127 
1132  struct ShapeFunctionData
1133  {
1142  bool is_nonzero_shape_function_component[value_type::n_independent_components];
1143 
1153  unsigned int row_index[value_type::n_independent_components];
1154 
1164  unsigned int single_nonzero_component_index;
1165  };
1166 
1170  Tensor();
1171 
1172 
1182  Tensor(const FEValuesBase<dim, spacedim> &fe_values_base,
1183  const unsigned int first_tensor_component);
1184 
1185 
1190  Tensor &operator=(const Tensor<2, dim, spacedim> &);
1191 
1208  value_type
1209  value (const unsigned int shape_function,
1210  const unsigned int q_point) const;
1211 
1225  divergence_type
1226  divergence (const unsigned int shape_function,
1227  const unsigned int q_point) const;
1228 
1246  template <class InputVector>
1247  void get_function_values (const InputVector &fe_function,
1248  std::vector<typename ProductType<value_type,typename InputVector::value_type>::type> &values) const;
1249 
1250 
1272  template <class InputVector>
1273  void get_function_divergences (const InputVector &fe_function,
1274  std::vector<typename ProductType<divergence_type,typename InputVector::value_type>::type> &divergences) const;
1275 
1276  private:
1281 
1286  const unsigned int first_tensor_component;
1287 
1291  std::vector<ShapeFunctionData> shape_function_data;
1292  };
1293 
1294 }
1295 
1296 
1297 namespace internal
1298 {
1299  namespace FEValuesViews
1300  {
1308  template <int dim, int spacedim>
1309  struct Cache
1310  {
1315  std::vector<::FEValuesViews::Scalar<dim,spacedim> > scalars;
1316  std::vector<::FEValuesViews::Vector<dim,spacedim> > vectors;
1317  std::vector<::FEValuesViews::SymmetricTensor<2,dim,spacedim> >
1318  symmetric_second_order_tensors;
1319  std::vector<::FEValuesViews::Tensor<2,dim,spacedim> >
1320  second_order_tensors;
1321 
1325  Cache (const FEValuesBase<dim,spacedim> &fe_values);
1326  };
1327  }
1328 }
1329 
1330 
1331 
1434 template <int dim, int spacedim>
1435 class FEValuesBase :
1436  public Subscriptor
1437 {
1438 public:
1442  static const unsigned int dimension = dim;
1443 
1447  static const unsigned int space_dimension = spacedim;
1448 
1452  const unsigned int n_quadrature_points;
1453 
1459  const unsigned int dofs_per_cell;
1460 
1461 
1469  FEValuesBase (const unsigned int n_q_points,
1470  const unsigned int dofs_per_cell,
1471  const UpdateFlags update_flags,
1474 
1475 
1479  ~FEValuesBase ();
1480 
1481 
1483 
1484 
1505  const double &shape_value (const unsigned int function_no,
1506  const unsigned int point_no) const;
1507 
1528  double shape_value_component (const unsigned int function_no,
1529  const unsigned int point_no,
1530  const unsigned int component) const;
1531 
1557  const Tensor<1,spacedim> &
1558  shape_grad (const unsigned int function_no,
1559  const unsigned int quadrature_point) const;
1560 
1578  shape_grad_component (const unsigned int function_no,
1579  const unsigned int point_no,
1580  const unsigned int component) const;
1581 
1601  const Tensor<2,spacedim> &
1602  shape_hessian (const unsigned int function_no,
1603  const unsigned int point_no) const;
1604 
1622  shape_hessian_component (const unsigned int function_no,
1623  const unsigned int point_no,
1624  const unsigned int component) const;
1625 
1645  const Tensor<3,spacedim> &
1646  shape_3rd_derivative (const unsigned int function_no,
1647  const unsigned int point_no) const;
1648 
1666  shape_3rd_derivative_component (const unsigned int function_no,
1667  const unsigned int point_no,
1668  const unsigned int component) const;
1669 
1671 
1673 
1712  template <class InputVector>
1713  void get_function_values (const InputVector &fe_function,
1714  std::vector<typename InputVector::value_type> &values) const;
1715 
1729  template <class InputVector>
1730  void get_function_values (const InputVector &fe_function,
1731  std::vector<Vector<typename InputVector::value_type> > &values) const;
1732 
1751  template <class InputVector>
1752  void get_function_values (const InputVector &fe_function,
1753  const VectorSlice<const std::vector<types::global_dof_index> > &indices,
1754  std::vector<typename InputVector::value_type> &values) const;
1755 
1777  template <class InputVector>
1778  void get_function_values (const InputVector &fe_function,
1779  const VectorSlice<const std::vector<types::global_dof_index> > &indices,
1780  std::vector<Vector<typename InputVector::value_type> > &values) const;
1781 
1782 
1813  template <class InputVector>
1814  void get_function_values (const InputVector &fe_function,
1815  const VectorSlice<const std::vector<types::global_dof_index> > &indices,
1816  VectorSlice<std::vector<std::vector<typename InputVector::value_type> > > values,
1817  const bool quadrature_points_fastest) const;
1818 
1820 
1822 
1861  template <class InputVector>
1862  void get_function_gradients (const InputVector &fe_function,
1863  std::vector<Tensor<1,spacedim,typename InputVector::value_type> > &gradients) const;
1864 
1881  template <class InputVector>
1882  void get_function_gradients (const InputVector &fe_function,
1883  std::vector<std::vector<Tensor<1,spacedim,typename InputVector::value_type> > > &gradients) const;
1884 
1891  template <class InputVector>
1892  void get_function_gradients (const InputVector &fe_function,
1893  const VectorSlice<const std::vector<types::global_dof_index> > &indices,
1894  std::vector<Tensor<1,spacedim,typename InputVector::value_type> > &gradients) const;
1895 
1902  template <class InputVector>
1903  void get_function_gradients (const InputVector &fe_function,
1904  const VectorSlice<const std::vector<types::global_dof_index> > &indices,
1905  VectorSlice<std::vector<std::vector<Tensor<1,spacedim,typename InputVector::value_type> > > > gradients,
1906  bool quadrature_points_fastest = false) const;
1907 
1909 
1911 
1951  template <class InputVector>
1952  void
1953  get_function_hessians (const InputVector &fe_function,
1954  std::vector<Tensor<2,spacedim,typename InputVector::value_type> > &hessians) const;
1955 
1973  template <class InputVector>
1974  void
1975  get_function_hessians (const InputVector &fe_function,
1976  std::vector<std::vector<Tensor<2,spacedim,typename InputVector::value_type> > > &hessians,
1977  bool quadrature_points_fastest = false) const;
1978 
1983  template <class InputVector>
1984  void get_function_hessians (
1985  const InputVector &fe_function,
1986  const VectorSlice<const std::vector<types::global_dof_index> > &indices,
1987  std::vector<Tensor<2,spacedim,typename InputVector::value_type> > &hessians) const;
1988 
1995  template <class InputVector>
1996  void get_function_hessians (
1997  const InputVector &fe_function,
1998  const VectorSlice<const std::vector<types::global_dof_index> > &indices,
1999  VectorSlice<std::vector<std::vector<Tensor<2,spacedim,typename InputVector::value_type> > > > hessians,
2000  bool quadrature_points_fastest = false) const;
2001 
2044  template <class InputVector>
2045  void
2046  get_function_laplacians (const InputVector &fe_function,
2047  std::vector<typename InputVector::value_type> &laplacians) const;
2048 
2068  template <class InputVector>
2069  void
2070  get_function_laplacians (const InputVector &fe_function,
2071  std::vector<Vector<typename InputVector::value_type> > &laplacians) const;
2072 
2079  template <class InputVector>
2081  const InputVector &fe_function,
2082  const VectorSlice<const std::vector<types::global_dof_index> > &indices,
2083  std::vector<typename InputVector::value_type> &laplacians) const;
2084 
2091  template <class InputVector>
2093  const InputVector &fe_function,
2094  const VectorSlice<const std::vector<types::global_dof_index> > &indices,
2095  std::vector<Vector<typename InputVector::value_type> > &laplacians) const;
2096 
2103  template <class InputVector>
2105  const InputVector &fe_function,
2106  const VectorSlice<const std::vector<types::global_dof_index> > &indices,
2107  std::vector<std::vector<typename InputVector::value_type> > &laplacians,
2108  bool quadrature_points_fastest = false) const;
2109 
2111 
2113 
2154  template <class InputVector>
2155  void
2156  get_function_third_derivatives (const InputVector &fe_function,
2157  std::vector<Tensor<3,spacedim,typename InputVector::value_type> > &third_derivatives) const;
2158 
2177  template <class InputVector>
2178  void
2179  get_function_third_derivatives (const InputVector &fe_function,
2180  std::vector<std::vector<Tensor<3,spacedim,typename InputVector::value_type> > > &third_derivatives,
2181  bool quadrature_points_fastest = false) const;
2182 
2187  template <class InputVector>
2189  const InputVector &fe_function,
2190  const VectorSlice<const std::vector<types::global_dof_index> > &indices,
2191  std::vector<Tensor<3,spacedim,typename InputVector::value_type> > &third_derivatives) const;
2192 
2199  template <class InputVector>
2201  const InputVector &fe_function,
2202  const VectorSlice<const std::vector<types::global_dof_index> > &indices,
2203  VectorSlice<std::vector<std::vector<Tensor<3,spacedim,typename InputVector::value_type> > > > third_derivatives,
2204  bool quadrature_points_fastest = false) const;
2206 
2208 
2209 
2215  const Point<spacedim> &
2216  quadrature_point (const unsigned int q) const;
2217 
2223  const std::vector<Point<spacedim> > &get_quadrature_points () const;
2224 
2240  double JxW (const unsigned int quadrature_point) const;
2241 
2245  const std::vector<double> &get_JxW_values () const;
2246 
2253  const DerivativeForm<1,dim,spacedim> &jacobian (const unsigned int quadrature_point) const;
2254 
2261  const std::vector<DerivativeForm<1,dim,spacedim> > &get_jacobians () const;
2262 
2270  const DerivativeForm<2,dim,spacedim> &jacobian_grad (const unsigned int quadrature_point) const;
2271 
2278  const std::vector<DerivativeForm<2,dim,spacedim> > &get_jacobian_grads () const;
2279 
2288  const Tensor<3,spacedim> &jacobian_pushed_forward_grad (const unsigned int quadrature_point) const;
2289 
2296  const std::vector<Tensor<3,spacedim> > &get_jacobian_pushed_forward_grads () const;
2297 
2306 
2313  const std::vector<DerivativeForm<3,dim,spacedim> > &get_jacobian_2nd_derivatives () const;
2314 
2325 
2332  const std::vector<Tensor<4,spacedim> > &get_jacobian_pushed_forward_2nd_derivatives () const;
2333 
2343 
2350  const std::vector<DerivativeForm<4,dim,spacedim> > &get_jacobian_3rd_derivatives () const;
2351 
2362 
2369  const std::vector<Tensor<5,spacedim> > &get_jacobian_pushed_forward_3rd_derivatives () const;
2370 
2377  const DerivativeForm<1,spacedim,dim> &inverse_jacobian (const unsigned int quadrature_point) const;
2378 
2385  const std::vector<DerivativeForm<1,spacedim,dim> > &get_inverse_jacobians () const;
2386 
2400  const Tensor<1,spacedim> &normal_vector (const unsigned int i) const;
2401 
2417  const std::vector<Tensor<1,spacedim> > &get_all_normal_vectors () const;
2418 
2429  std::vector<Point<spacedim> > get_normal_vectors () const DEAL_II_DEPRECATED;
2430 
2437  void transform (std::vector<Tensor<1,spacedim> > &transformed,
2438  const std::vector<Tensor<1,dim> > &original,
2439  MappingType mapping) const DEAL_II_DEPRECATED;
2440 
2442 
2444 
2445 
2454  const FEValuesViews::Scalar<dim,spacedim> &
2455  operator[] (const FEValuesExtractors::Scalar &scalar) const;
2456 
2465  const FEValuesViews::Vector<dim,spacedim> &
2466  operator[] (const FEValuesExtractors::Vector &vector) const;
2467 
2477  const FEValuesViews::SymmetricTensor<2,dim,spacedim> &
2478  operator[] (const FEValuesExtractors::SymmetricTensor<2> &tensor) const;
2479 
2480 
2489  const FEValuesViews::Tensor<2,dim,spacedim> &
2490  operator[] (const FEValuesExtractors::Tensor<2> &tensor) const;
2491 
2493 
2495 
2496 
2500  const Mapping<dim,spacedim> &get_mapping () const;
2501 
2505  const FiniteElement<dim,spacedim> &get_fe () const;
2506 
2510  UpdateFlags get_update_flags () const;
2511 
2515  const typename Triangulation<dim,spacedim>::cell_iterator get_cell () const;
2516 
2523 
2528  std::size_t memory_consumption () const;
2530 
2531 
2539  char *,
2540  << "You are requesting information from an FEValues/FEFaceValues/FESubfaceValues "
2541  << "object for which this kind of information has not been computed. What "
2542  << "information these objects compute is determined by the update_* flags you "
2543  << "pass to the constructor. Here, the operation you are attempting requires "
2544  << "the <" << arg1 << "> flag to be set, but it was apparently not specified "
2545  << "upon construction.");
2570  int,
2571  << "The shape function with index " << arg1
2572  << " is not primitive, i.e. it is vector-valued and "
2573  << "has more than one non-zero vector component. This "
2574  << "function cannot be called for these shape functions. "
2575  << "Maybe you want to use the same function with the "
2576  << "_component suffix?");
2583 
2584 protected:
2613  class CellIteratorBase;
2614 
2619  template <typename CI> class CellIterator;
2621 
2627  std_cxx11::unique_ptr<const CellIteratorBase> present_cell;
2628 
2636  boost::signals2::connection tria_listener_refinement;
2637 
2645  boost::signals2::connection tria_listener_mesh_transform;
2646 
2652  void invalidate_present_cell ();
2653 
2663  void
2664  maybe_invalidate_previous_present_cell (const typename Triangulation<dim,spacedim>::cell_iterator &cell);
2665 
2669  const SmartPointer<const Mapping<dim,spacedim>,FEValuesBase<dim,spacedim> > mapping;
2670 
2676  std_cxx11::unique_ptr<typename Mapping<dim,spacedim>::InternalDataBase> mapping_data;
2677 
2682  ::internal::FEValues::MappingRelatedData<dim, spacedim> mapping_output;
2683 
2684 
2689  const SmartPointer<const FiniteElement<dim,spacedim>,FEValuesBase<dim,spacedim> > fe;
2690 
2696  std_cxx11::unique_ptr<typename FiniteElement<dim,spacedim>::InternalDataBase> fe_data;
2697 
2702  ::internal::FEValues::FiniteElementRelatedData<dim, spacedim> finite_element_output;
2703 
2704 
2709 
2718  UpdateFlags compute_update_flags (const UpdateFlags update_flags) const;
2719 
2726 
2732  void
2733  check_cell_similarity (const typename Triangulation<dim,spacedim>::cell_iterator &cell);
2734 
2735 private:
2740  FEValuesBase (const FEValuesBase &);
2741 
2746  FEValuesBase &operator= (const FEValuesBase &);
2747 
2752 
2757  template <int, int> friend class FEValuesViews::Scalar;
2758  template <int, int> friend class FEValuesViews::Vector;
2759  template <int, int, int> friend class FEValuesViews::SymmetricTensor;
2760  template <int, int, int> friend class FEValuesViews::Tensor;
2761 };
2762 
2763 
2764 
2775 template <int dim, int spacedim=dim>
2776 class FEValues : public FEValuesBase<dim,spacedim>
2777 {
2778 public:
2783  static const unsigned int integral_dimension = dim;
2784 
2789  FEValues (const Mapping<dim,spacedim> &mapping,
2790  const FiniteElement<dim,spacedim> &fe,
2791  const Quadrature<dim> &quadrature,
2792  const UpdateFlags update_flags);
2793 
2800  const Quadrature<dim> &quadrature,
2801  const UpdateFlags update_flags);
2802 
2809  template <template <int, int> class DoFHandlerType, bool level_dof_access>
2810  void reinit (const TriaIterator<DoFCellAccessor<DoFHandlerType<dim,spacedim>,level_dof_access> > &cell);
2811 
2825  void reinit (const typename Triangulation<dim,spacedim>::cell_iterator &cell);
2826 
2831  const Quadrature<dim> &get_quadrature () const;
2832 
2837  std::size_t memory_consumption () const;
2838 
2853  const FEValues<dim,spacedim> &get_present_fe_values () const;
2854 
2855 private:
2860 
2864  void initialize (const UpdateFlags update_flags);
2865 
2872  void do_reinit ();
2873 };
2874 
2875 
2886 template <int dim, int spacedim=dim>
2887 class FEFaceValuesBase : public FEValuesBase<dim,spacedim>
2888 {
2889 public:
2894  static const unsigned int integral_dimension = dim-1;
2895 
2907  FEFaceValuesBase (const unsigned int n_q_points,
2908  const unsigned int dofs_per_cell,
2909  const UpdateFlags update_flags,
2912  const Quadrature<dim-1>& quadrature);
2913 
2921  const Tensor<1,spacedim> &boundary_form (const unsigned int i) const;
2922 
2929  const std::vector<Tensor<1,spacedim> > &get_boundary_forms () const;
2930 
2935  unsigned int get_face_index() const;
2936 
2941  const Quadrature<dim-1> & get_quadrature () const;
2942 
2947  std::size_t memory_consumption () const;
2948 
2949 protected:
2950 
2955  unsigned int present_face_index;
2956 
2960  const Quadrature<dim-1> quadrature;
2961 };
2962 
2963 
2964 
2979 template <int dim, int spacedim=dim>
2980 class FEFaceValues : public FEFaceValuesBase<dim,spacedim>
2981 {
2982 public:
2987  static const unsigned int dimension = dim;
2988 
2989  static const unsigned int space_dimension = spacedim;
2990 
2995  static const unsigned int integral_dimension = dim-1;
2996 
3003  const Quadrature<dim-1> &quadrature,
3004  const UpdateFlags update_flags);
3005 
3012  const Quadrature<dim-1> &quadrature,
3013  const UpdateFlags update_flags);
3014 
3019  template <template <int, int> class DoFHandlerType, bool level_dof_access>
3020  void reinit (const TriaIterator<DoFCellAccessor<DoFHandlerType<dim,spacedim>,level_dof_access> > &cell,
3021  const unsigned int face_no);
3022 
3036  void reinit (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
3037  const unsigned int face_no);
3038 
3053  const FEFaceValues<dim,spacedim> &get_present_fe_values () const;
3054 private:
3055 
3059  void initialize (const UpdateFlags update_flags);
3060 
3067  void do_reinit (const unsigned int face_no);
3068 };
3069 
3070 
3088 template <int dim, int spacedim=dim>
3089 class FESubfaceValues : public FEFaceValuesBase<dim,spacedim>
3090 {
3091 public:
3095  static const unsigned int dimension = dim;
3096 
3100  static const unsigned int space_dimension = spacedim;
3101 
3106  static const unsigned int integral_dimension = dim-1;
3107 
3114  const Quadrature<dim-1> &face_quadrature,
3115  const UpdateFlags update_flags);
3116 
3123  const Quadrature<dim-1> &face_quadrature,
3124  const UpdateFlags update_flags);
3125 
3132  template <template <int, int> class DoFHandlerType, bool level_dof_access>
3133  void reinit (const TriaIterator<DoFCellAccessor<DoFHandlerType<dim,spacedim>,level_dof_access> > &cell,
3134  const unsigned int face_no,
3135  const unsigned int subface_no);
3136 
3150  void reinit (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
3151  const unsigned int face_no,
3152  const unsigned int subface_no);
3153 
3168  const FESubfaceValues<dim,spacedim> &get_present_fe_values () const;
3169 
3175  DeclException0 (ExcReinitCalledWithBoundaryFace);
3176 
3182  DeclException0 (ExcFaceHasNoSubfaces);
3183 
3184 private:
3185 
3189  void initialize (const UpdateFlags update_flags);
3190 
3197  void do_reinit (const unsigned int face_no,
3198  const unsigned int subface_no);
3199 };
3200 
3201 
3202 #ifndef DOXYGEN
3203 
3204 
3205 /*------------------------ Inline functions: namespace FEValuesViews --------*/
3206 
3207 namespace FEValuesViews
3208 {
3209  template <int dim, int spacedim>
3210  inline
3211  typename Scalar<dim,spacedim>::value_type
3212  Scalar<dim,spacedim>::value (const unsigned int shape_function,
3213  const unsigned int q_point) const
3214  {
3215  Assert (shape_function < fe_values->fe->dofs_per_cell,
3216  ExcIndexRange (shape_function, 0, fe_values->fe->dofs_per_cell));
3217  Assert (fe_values->update_flags & update_values,
3218  ((typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_values"))));
3219 
3220  // an adaptation of the FEValuesBase::shape_value_component function
3221  // except that here we know the component as fixed and we have
3222  // pre-computed and cached a bunch of information. See the comments there.
3223  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
3224  return fe_values->finite_element_output.shape_values(shape_function_data[shape_function]
3225  .row_index,
3226  q_point);
3227  else
3228  return 0;
3229  }
3230 
3231 
3232 
3233 
3234  template <int dim, int spacedim>
3235  inline
3236  typename Scalar<dim,spacedim>::gradient_type
3237  Scalar<dim,spacedim>::gradient (const unsigned int shape_function,
3238  const unsigned int q_point) const
3239  {
3240  Assert (shape_function < fe_values->fe->dofs_per_cell,
3241  ExcIndexRange (shape_function, 0, fe_values->fe->dofs_per_cell));
3242  Assert (fe_values->update_flags & update_gradients,
3243  (typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_gradients")));
3244 
3245  // an adaptation of the
3246  // FEValuesBase::shape_grad_component
3247  // function except that here we know the
3248  // component as fixed and we have
3249  // pre-computed and cached a bunch of
3250  // information. See the comments there.
3251  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
3252  return fe_values->finite_element_output.shape_gradients[shape_function_data[shape_function]
3253  .row_index][q_point];
3254  else
3255  return gradient_type();
3256  }
3257 
3258 
3259 
3260  template <int dim, int spacedim>
3261  inline
3262  typename Scalar<dim,spacedim>::hessian_type
3263  Scalar<dim,spacedim>::hessian (const unsigned int shape_function,
3264  const unsigned int q_point) const
3265  {
3266  Assert (shape_function < fe_values->fe->dofs_per_cell,
3267  ExcIndexRange (shape_function, 0, fe_values->fe->dofs_per_cell));
3268  Assert (fe_values->update_flags & update_hessians,
3269  (typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_hessians")));
3270 
3271  // an adaptation of the
3272  // FEValuesBase::shape_hessian_component
3273  // function except that here we know the
3274  // component as fixed and we have
3275  // pre-computed and cached a bunch of
3276  // information. See the comments there.
3277  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
3278  return fe_values->finite_element_output.shape_hessians[shape_function_data[shape_function].row_index][q_point];
3279  else
3280  return hessian_type();
3281  }
3282 
3283 
3284 
3285  template <int dim, int spacedim>
3286  inline
3287  typename Scalar<dim,spacedim>::third_derivative_type
3288  Scalar<dim,spacedim>::third_derivative (const unsigned int shape_function,
3289  const unsigned int q_point) const
3290  {
3291  Assert (shape_function < fe_values->fe->dofs_per_cell,
3292  ExcIndexRange (shape_function, 0, fe_values->fe->dofs_per_cell));
3293  Assert (fe_values->update_flags & update_3rd_derivatives,
3294  (typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_3rd_derivatives")));
3295 
3296  // an adaptation of the
3297  // FEValuesBase::shape_3rdderivative_component
3298  // function except that here we know the
3299  // component as fixed and we have
3300  // pre-computed and cached a bunch of
3301  // information. See the comments there.
3302  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
3303  return fe_values->finite_element_output.shape_3rd_derivatives[shape_function_data[shape_function].row_index][q_point];
3304  else
3305  return third_derivative_type();
3306  }
3307 
3308 
3309 
3310  template <int dim, int spacedim>
3311  inline
3313  Vector<dim,spacedim>::value (const unsigned int shape_function,
3314  const unsigned int q_point) const
3315  {
3316  Assert (shape_function < fe_values->fe->dofs_per_cell,
3317  ExcIndexRange (shape_function, 0, fe_values->fe->dofs_per_cell));
3318  Assert (fe_values->update_flags & update_values,
3319  (typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_values")));
3320 
3321  // same as for the scalar case except
3322  // that we have one more index
3323  const int snc = shape_function_data[shape_function].single_nonzero_component;
3324  if (snc == -2)
3325  return value_type();
3326  else if (snc != -1)
3327  {
3328  value_type return_value;
3329  return_value[shape_function_data[shape_function].single_nonzero_component_index]
3330  = fe_values->finite_element_output.shape_values(snc,q_point);
3331  return return_value;
3332  }
3333  else
3334  {
3335  value_type return_value;
3336  for (unsigned int d=0; d<dim; ++d)
3337  if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
3338  return_value[d]
3339  = fe_values->finite_element_output.shape_values(shape_function_data[shape_function].row_index[d],q_point);
3340 
3341  return return_value;
3342  }
3343  }
3344 
3345 
3346 
3347  template <int dim, int spacedim>
3348  inline
3350  Vector<dim,spacedim>::gradient (const unsigned int shape_function,
3351  const unsigned int q_point) const
3352  {
3353  Assert (shape_function < fe_values->fe->dofs_per_cell,
3354  ExcIndexRange (shape_function, 0, fe_values->fe->dofs_per_cell));
3355  Assert (fe_values->update_flags & update_gradients,
3356  (typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_gradients")));
3357 
3358  // same as for the scalar case except
3359  // that we have one more index
3360  const int snc = shape_function_data[shape_function].single_nonzero_component;
3361  if (snc == -2)
3362  return gradient_type();
3363  else if (snc != -1)
3364  {
3365  gradient_type return_value;
3366  return_value[shape_function_data[shape_function].single_nonzero_component_index]
3367  = fe_values->finite_element_output.shape_gradients[snc][q_point];
3368  return return_value;
3369  }
3370  else
3371  {
3372  gradient_type return_value;
3373  for (unsigned int d=0; d<dim; ++d)
3374  if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
3375  return_value[d]
3376  = fe_values->finite_element_output.shape_gradients[shape_function_data[shape_function].row_index[d]][q_point];
3377 
3378  return return_value;
3379  }
3380  }
3381 
3382 
3383 
3384  template <int dim, int spacedim>
3385  inline
3387  Vector<dim,spacedim>::divergence (const unsigned int shape_function,
3388  const unsigned int q_point) const
3389  {
3390  // this function works like in
3391  // the case above
3392  Assert (shape_function < fe_values->fe->dofs_per_cell,
3393  ExcIndexRange (shape_function, 0, fe_values->fe->dofs_per_cell));
3394  Assert (fe_values->update_flags & update_gradients,
3395  (typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_gradients")));
3396 
3397  // same as for the scalar case except
3398  // that we have one more index
3399  const int snc = shape_function_data[shape_function].single_nonzero_component;
3400  if (snc == -2)
3401  return divergence_type();
3402  else if (snc != -1)
3403  return
3404  fe_values->finite_element_output.shape_gradients[snc][q_point][shape_function_data[shape_function].single_nonzero_component_index];
3405  else
3406  {
3407  divergence_type return_value = 0;
3408  for (unsigned int d=0; d<dim; ++d)
3409  if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
3410  return_value
3411  += fe_values->finite_element_output.shape_gradients[shape_function_data[shape_function].row_index[d]][q_point][d];
3412 
3413  return return_value;
3414  }
3415  }
3416 
3417 
3418 
3419  template <int dim, int spacedim>
3420  inline
3422  Vector<dim,spacedim>::curl (const unsigned int shape_function, const unsigned int q_point) const
3423  {
3424  // this function works like in the case above
3425 
3426  Assert (shape_function < fe_values->fe->dofs_per_cell,
3427  ExcIndexRange (shape_function, 0, fe_values->fe->dofs_per_cell));
3428  Assert (fe_values->update_flags & update_gradients,
3429  (typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_gradients")));
3430  // same as for the scalar case except that we have one more index
3431  const int snc = shape_function_data[shape_function].single_nonzero_component;
3432 
3433  if (snc == -2)
3434  return curl_type ();
3435 
3436  else
3437  switch (dim)
3438  {
3439  case 1:
3440  {
3441  Assert (false, ExcMessage("Computing the curl in 1d is not a useful operation"));
3442  return curl_type ();
3443  }
3444 
3445  case 2:
3446  {
3447  if (snc != -1)
3448  {
3449  curl_type return_value;
3450 
3451  // the single
3452  // nonzero component
3453  // can only be zero
3454  // or one in 2d
3455  if (shape_function_data[shape_function].single_nonzero_component_index == 0)
3456  return_value[0] = -1.0 * fe_values->finite_element_output.shape_gradients[snc][q_point][1];
3457  else
3458  return_value[0] = fe_values->finite_element_output.shape_gradients[snc][q_point][0];
3459 
3460  return return_value;
3461  }
3462 
3463  else
3464  {
3465  curl_type return_value;
3466 
3467  return_value[0] = 0.0;
3468 
3469  if (shape_function_data[shape_function].is_nonzero_shape_function_component[0])
3470  return_value[0]
3471  -= fe_values->finite_element_output.shape_gradients[shape_function_data[shape_function].row_index[0]][q_point][1];
3472 
3473  if (shape_function_data[shape_function].is_nonzero_shape_function_component[1])
3474  return_value[0]
3475  += fe_values->finite_element_output.shape_gradients[shape_function_data[shape_function].row_index[1]][q_point][0];
3476 
3477  return return_value;
3478  }
3479  }
3480 
3481  case 3:
3482  {
3483  if (snc != -1)
3484  {
3485  curl_type return_value;
3486 
3487  switch (shape_function_data[shape_function].single_nonzero_component_index)
3488  {
3489  case 0:
3490  {
3491  return_value[0] = 0;
3492  return_value[1] = fe_values->finite_element_output.shape_gradients[snc][q_point][2];
3493  return_value[2] = -1.0 * fe_values->finite_element_output.shape_gradients[snc][q_point][1];
3494  return return_value;
3495  }
3496 
3497  case 1:
3498  {
3499  return_value[0] = -1.0 * fe_values->finite_element_output.shape_gradients[snc][q_point][2];
3500  return_value[1] = 0;
3501  return_value[2] = fe_values->finite_element_output.shape_gradients[snc][q_point][0];
3502  return return_value;
3503  }
3504 
3505  default:
3506  {
3507  return_value[0] = fe_values->finite_element_output.shape_gradients[snc][q_point][1];
3508  return_value[1] = -1.0 * fe_values->finite_element_output.shape_gradients[snc][q_point][0];
3509  return_value[2] = 0;
3510  return return_value;
3511  }
3512  }
3513  }
3514 
3515  else
3516  {
3517  curl_type return_value;
3518 
3519  for (unsigned int i = 0; i < dim; ++i)
3520  return_value[i] = 0.0;
3521 
3522  if (shape_function_data[shape_function].is_nonzero_shape_function_component[0])
3523  {
3524  return_value[1]
3525  += fe_values->finite_element_output.shape_gradients[shape_function_data[shape_function].row_index[0]][q_point][2];
3526  return_value[2]
3527  -= fe_values->finite_element_output.shape_gradients[shape_function_data[shape_function].row_index[0]][q_point][1];
3528  }
3529 
3530  if (shape_function_data[shape_function].is_nonzero_shape_function_component[1])
3531  {
3532  return_value[0]
3533  -= fe_values->finite_element_output.shape_gradients[shape_function_data[shape_function].row_index[1]][q_point][2];
3534  return_value[2]
3535  += fe_values->finite_element_output.shape_gradients[shape_function_data[shape_function].row_index[1]][q_point][0];
3536  }
3537 
3538  if (shape_function_data[shape_function].is_nonzero_shape_function_component[2])
3539  {
3540  return_value[0]
3541  += fe_values->finite_element_output.shape_gradients[shape_function_data[shape_function].row_index[2]][q_point][1];
3542  return_value[1]
3543  -= fe_values->finite_element_output.shape_gradients[shape_function_data[shape_function].row_index[2]][q_point][0];
3544  }
3545 
3546  return return_value;
3547  }
3548  }
3549  }
3550  // should not end up here
3551  Assert (false, ExcInternalError());
3552  return curl_type();
3553  }
3554 
3555  template <int dim, int spacedim>
3556  inline
3558  Vector<dim,spacedim>::hessian (const unsigned int shape_function,
3559  const unsigned int q_point) const
3560  {
3561  // this function works like in
3562  // the case above
3563  Assert (shape_function < fe_values->fe->dofs_per_cell,
3564  ExcIndexRange (shape_function, 0, fe_values->fe->dofs_per_cell));
3565  Assert (fe_values->update_flags & update_hessians,
3566  (typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_hessians")));
3567 
3568  // same as for the scalar case except
3569  // that we have one more index
3570  const int snc = shape_function_data[shape_function].single_nonzero_component;
3571  if (snc == -2)
3572  return hessian_type();
3573  else if (snc != -1)
3574  {
3575  hessian_type return_value;
3576  return_value[shape_function_data[shape_function].single_nonzero_component_index]
3577  = fe_values->finite_element_output.shape_hessians[snc][q_point];
3578  return return_value;
3579  }
3580  else
3581  {
3582  hessian_type return_value;
3583  for (unsigned int d=0; d<dim; ++d)
3584  if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
3585  return_value[d]
3586  = fe_values->finite_element_output.shape_hessians[shape_function_data[shape_function].row_index[d]][q_point];
3587 
3588  return return_value;
3589  }
3590  }
3591 
3592  template <int dim, int spacedim>
3593  inline
3595  Vector<dim,spacedim>::third_derivative (const unsigned int shape_function,
3596  const unsigned int q_point) const
3597  {
3598  // this function works like in
3599  // the case above
3600  Assert (shape_function < fe_values->fe->dofs_per_cell,
3601  ExcIndexRange (shape_function, 0, fe_values->fe->dofs_per_cell));
3602  Assert (fe_values->update_flags & update_3rd_derivatives,
3603  (typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_3rd_derivatives")));
3604 
3605  // same as for the scalar case except
3606  // that we have one more index
3607  const int snc = shape_function_data[shape_function].single_nonzero_component;
3608  if (snc == -2)
3609  return third_derivative_type();
3610  else if (snc != -1)
3611  {
3612  third_derivative_type return_value;
3613  return_value[shape_function_data[shape_function].single_nonzero_component_index]
3614  = fe_values->finite_element_output.shape_3rd_derivatives[snc][q_point];
3615  return return_value;
3616  }
3617  else
3618  {
3619  third_derivative_type return_value;
3620  for (unsigned int d=0; d<dim; ++d)
3621  if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
3622  return_value[d]
3623  = fe_values->finite_element_output.shape_3rd_derivatives[shape_function_data[shape_function].row_index[d]][q_point];
3624 
3625  return return_value;
3626  }
3627  }
3628 
3629 
3630  namespace
3631  {
3636  inline
3637  ::SymmetricTensor<2,1>
3638  symmetrize_single_row (const unsigned int n,
3639  const ::Tensor<1,1> &t)
3640  {
3641  Assert (n < 1, ExcIndexRange (n, 0, 1));
3642  (void)n; // removes -Wunused-parameter warning in optimized mode
3643 
3644  const double array[1] = { t[0] };
3645  return ::SymmetricTensor<2,1>(array);
3646  }
3647 
3648 
3649  inline
3650  ::SymmetricTensor<2,2>
3651  symmetrize_single_row (const unsigned int n,
3652  const ::Tensor<1,2> &t)
3653  {
3654  switch (n)
3655  {
3656  case 0:
3657  {
3658  const double array[3] = { t[0], 0, t[1]/2 };
3659  return ::SymmetricTensor<2,2>(array);
3660  }
3661  case 1:
3662  {
3663  const double array[3] = { 0, t[1], t[0]/2 };
3664  return ::SymmetricTensor<2,2>(array);
3665  }
3666  default:
3667  {
3668  Assert (false, ExcIndexRange (n, 0, 2));
3669  return ::SymmetricTensor<2,2>();
3670  }
3671  }
3672  }
3673 
3674 
3675  inline
3676  ::SymmetricTensor<2,3>
3677  symmetrize_single_row (const unsigned int n,
3678  const ::Tensor<1,3> &t)
3679  {
3680  switch (n)
3681  {
3682  case 0:
3683  {
3684  const double array[6] = { t[0], 0, 0, t[1]/2, t[2]/2, 0 };
3685  return ::SymmetricTensor<2,3>(array);
3686  }
3687  case 1:
3688  {
3689  const double array[6] = { 0, t[1], 0, t[0]/2, 0, t[2]/2 };
3690  return ::SymmetricTensor<2,3>(array);
3691  }
3692  case 2:
3693  {
3694  const double array[6] = { 0, 0, t[2], 0, t[0]/2, t[1]/2 };
3695  return ::SymmetricTensor<2,3>(array);
3696  }
3697  default:
3698  {
3699  Assert (false, ExcIndexRange (n, 0, 3));
3700  return ::SymmetricTensor<2,3>();
3701  }
3702  }
3703  }
3704  }
3705 
3706 
3707  template <int dim, int spacedim>
3708  inline
3710  Vector<dim,spacedim>::symmetric_gradient (const unsigned int shape_function,
3711  const unsigned int q_point) const
3712  {
3713  Assert (shape_function < fe_values->fe->dofs_per_cell,
3714  ExcIndexRange (shape_function, 0, fe_values->fe->dofs_per_cell));
3715  Assert (fe_values->update_flags & update_gradients,
3716  (typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_gradients")));
3717 
3718  // same as for the scalar case except
3719  // that we have one more index
3720  const int snc = shape_function_data[shape_function].single_nonzero_component;
3721  if (snc == -2)
3722  return symmetric_gradient_type();
3723  else if (snc != -1)
3724  return symmetrize_single_row (shape_function_data[shape_function].single_nonzero_component_index,
3725  fe_values->finite_element_output.shape_gradients[snc][q_point]);
3726  else
3727  {
3728  gradient_type return_value;
3729  for (unsigned int d=0; d<dim; ++d)
3730  if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
3731  return_value[d]
3732  = fe_values->finite_element_output.shape_gradients[shape_function_data[shape_function].row_index[d]][q_point];
3733 
3734  return symmetrize(return_value);
3735  }
3736  }
3737 
3738 
3739 
3740  template <int dim, int spacedim>
3741  inline
3743  SymmetricTensor<2, dim, spacedim>::value (const unsigned int shape_function,
3744  const unsigned int q_point) const
3745  {
3746  Assert (shape_function < fe_values->fe->dofs_per_cell,
3747  ExcIndexRange (shape_function, 0, fe_values->fe->dofs_per_cell));
3748  Assert (fe_values->update_flags & update_values,
3749  (typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_values")));
3750 
3751  // similar to the vector case where we
3752  // have more then one index and we need
3753  // to convert between unrolled and
3754  // component indexing for tensors
3755  const int snc
3756  = shape_function_data[shape_function].single_nonzero_component;
3757 
3758  if (snc == -2)
3759  {
3760  // shape function is zero for the
3761  // selected components
3762  return value_type();
3763 
3764  }
3765  else if (snc != -1)
3766  {
3767  value_type return_value;
3768  const unsigned int comp =
3769  shape_function_data[shape_function].single_nonzero_component_index;
3770  return_value[value_type::unrolled_to_component_indices(comp)]
3771  = fe_values->finite_element_output.shape_values(snc,q_point);
3772  return return_value;
3773  }
3774  else
3775  {
3776  value_type return_value;
3777  for (unsigned int d = 0; d < value_type::n_independent_components; ++d)
3778  if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
3779  return_value[value_type::unrolled_to_component_indices(d)]
3780  = fe_values->finite_element_output.shape_values(shape_function_data[shape_function].row_index[d],q_point);
3781  return return_value;
3782  }
3783  }
3784 
3785 
3786  template <int dim, int spacedim>
3787  inline
3789  SymmetricTensor<2, dim, spacedim>::divergence(const unsigned int shape_function,
3790  const unsigned int q_point) const
3791  {
3792  Assert (shape_function < fe_values->fe->dofs_per_cell,
3793  ExcIndexRange (shape_function, 0, fe_values->fe->dofs_per_cell));
3794  Assert (fe_values->update_flags & update_gradients,
3795  (typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_gradients")));
3796 
3797  const int snc = shape_function_data[shape_function].single_nonzero_component;
3798 
3799  if (snc == -2)
3800  {
3801  // shape function is zero for the
3802  // selected components
3803  return divergence_type();
3804  }
3805  else if (snc != -1)
3806  {
3807  // we have a single non-zero component
3808  // when the symmetric tensor is
3809  // represented in unrolled form.
3810  // this implies we potentially have
3811  // two non-zero components when
3812  // represented in component form! we
3813  // will only have one non-zero entry
3814  // if the non-zero component lies on
3815  // the diagonal of the tensor.
3816  //
3817  // the divergence of a second-order tensor
3818  // is a first order tensor.
3819  //
3820  // assume the second-order tensor is
3821  // A with components A_{ij}. then
3822  // A_{ij} = A_{ji} and there is only
3823  // one (if diagonal) or two non-zero
3824  // entries in the tensorial
3825  // representation. define the
3826  // divergence as:
3827  // b_i := \dfrac{\partial phi_{ij}}{\partial x_j}.
3828  // (which is incidentally also
3829  // b_j := \dfrac{\partial phi_{ij}}{\partial x_i}).
3830  // In both cases, a sum is implied.
3831  //
3832  // Now, we know the nonzero component
3833  // in unrolled form: it is indicated
3834  // by 'snc'. we can figure out which
3835  // tensor components belong to this:
3836  const unsigned int comp =
3837  shape_function_data[shape_function].single_nonzero_component_index;
3838  const unsigned int ii = value_type::unrolled_to_component_indices(comp)[0];
3839  const unsigned int jj = value_type::unrolled_to_component_indices(comp)[1];
3840 
3841  // given the form of the divergence
3842  // above, if ii=jj there is only a
3843  // single nonzero component of the
3844  // full tensor and the gradient
3845  // equals
3846  // b_ii := \dfrac{\partial phi_{ii,ii}}{\partial x_ii}.
3847  // all other entries of 'b' are zero
3848  //
3849  // on the other hand, if ii!=jj, then
3850  // there are two nonzero entries in
3851  // the full tensor and
3852  // b_ii := \dfrac{\partial phi_{ii,jj}}{\partial x_ii}.
3853  // b_jj := \dfrac{\partial phi_{ii,jj}}{\partial x_jj}.
3854  // again, all other entries of 'b' are
3855  // zero
3856  const ::Tensor<1, spacedim> phi_grad = fe_values->finite_element_output.shape_gradients[snc][q_point];
3857 
3858  divergence_type return_value;
3859  return_value[ii] = phi_grad[jj];
3860 
3861  if (ii != jj)
3862  return_value[jj] = phi_grad[ii];
3863 
3864  return return_value;
3865 
3866  }
3867  else
3868  {
3869  Assert (false, ExcNotImplemented());
3870  divergence_type return_value;
3871  return return_value;
3872  }
3873  }
3874 
3875  template <int dim, int spacedim>
3876  inline
3878  Tensor<2, dim, spacedim>::value (const unsigned int shape_function,
3879  const unsigned int q_point) const
3880  {
3881  Assert (shape_function < fe_values->fe->dofs_per_cell,
3882  ExcIndexRange (shape_function, 0, fe_values->fe->dofs_per_cell));
3883  Assert (fe_values->update_flags & update_values,
3884  (typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_values")));
3885 
3886  // similar to the vector case where we
3887  // have more then one index and we need
3888  // to convert between unrolled and
3889  // component indexing for tensors
3890  const int snc
3891  = shape_function_data[shape_function].single_nonzero_component;
3892 
3893  if (snc == -2)
3894  {
3895  // shape function is zero for the
3896  // selected components
3897  return value_type();
3898 
3899  }
3900  else if (snc != -1)
3901  {
3902  value_type return_value;
3903  const unsigned int comp =
3904  shape_function_data[shape_function].single_nonzero_component_index;
3906  return_value[indices] = fe_values->finite_element_output.shape_values(snc,q_point);
3907  return return_value;
3908  }
3909  else
3910  {
3911  value_type return_value;
3912  for (unsigned int d = 0; d < dim*dim; ++d)
3913  if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
3914  {
3916  return_value[indices]
3917  = fe_values->finite_element_output.shape_values(shape_function_data[shape_function].row_index[d],q_point);
3918  }
3919  return return_value;
3920  }
3921  }
3922 
3923 
3924  template <int dim, int spacedim>
3925  inline
3927  Tensor<2, dim, spacedim>::divergence(const unsigned int shape_function,
3928  const unsigned int q_point) const
3929  {
3930  Assert (shape_function < fe_values->fe->dofs_per_cell,
3931  ExcIndexRange (shape_function, 0, fe_values->fe->dofs_per_cell));
3932  Assert (fe_values->update_flags & update_gradients,
3933  (typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_gradients")));
3934 
3935  const int snc = shape_function_data[shape_function].single_nonzero_component;
3936 
3937  if (snc == -2)
3938  {
3939  // shape function is zero for the
3940  // selected components
3941  return divergence_type();
3942  }
3943  else if (snc != -1)
3944  {
3945  // we have a single non-zero component
3946  // when the tensor is
3947  // represented in unrolled form.
3948  //
3949  // the divergence of a second-order tensor
3950  // is a first order tensor.
3951  //
3952  // assume the second-order tensor is
3953  // A with components A_{ij}.
3954  // divergence as:
3955  // b_j := \dfrac{\partial phi_{ij}}{\partial x_i}.
3956  //
3957  // Now, we know the nonzero component
3958  // in unrolled form: it is indicated
3959  // by 'snc'. we can figure out which
3960  // tensor components belong to this:
3961  const unsigned int comp =
3962  shape_function_data[shape_function].single_nonzero_component_index;
3964  const unsigned int ii = indices[0];
3965  const unsigned int jj = indices[1];
3966 
3967  const ::Tensor<1, spacedim> phi_grad = fe_values->finite_element_output.shape_gradients[snc][q_point];
3968 
3969  divergence_type return_value;
3970  return_value[jj] = phi_grad[ii];
3971 
3972  return return_value;
3973 
3974  }
3975  else
3976  {
3977  Assert (false, ExcNotImplemented());
3978  divergence_type return_value;
3979  return return_value;
3980  }
3981  }
3982 }
3983 
3984 
3985 
3986 /*------------------------ Inline functions: FEValuesBase ------------------------*/
3987 
3988 
3989 
3990 template <int dim, int spacedim>
3991 inline
3994 operator[] (const FEValuesExtractors::Scalar &scalar) const
3995 {
3996  Assert (scalar.component < fe_values_views_cache.scalars.size(),
3997  ExcIndexRange (scalar.component,
3998  0, fe_values_views_cache.scalars.size()));
3999 
4000  return fe_values_views_cache.scalars[scalar.component];
4001 }
4002 
4003 
4004 
4005 template <int dim, int spacedim>
4006 inline
4009 operator[] (const FEValuesExtractors::Vector &vector) const
4010 {
4011  Assert (vector.first_vector_component <
4012  fe_values_views_cache.vectors.size(),
4014  0, fe_values_views_cache.vectors.size()));
4015 
4016  return fe_values_views_cache.vectors[vector.first_vector_component];
4017 }
4018 
4019 template <int dim, int spacedim>
4020 inline
4024 {
4025  Assert (tensor.first_tensor_component <
4026  fe_values_views_cache.symmetric_second_order_tensors.size(),
4028  0, fe_values_views_cache.symmetric_second_order_tensors.size()));
4029 
4030  return fe_values_views_cache.symmetric_second_order_tensors[tensor.first_tensor_component];
4031 }
4032 
4033 template <int dim, int spacedim>
4034 inline
4037 operator[] (const FEValuesExtractors::Tensor<2> &tensor) const
4038 {
4039  Assert (tensor.first_tensor_component <
4040  fe_values_views_cache.second_order_tensors.size(),
4042  0, fe_values_views_cache.second_order_tensors.size()));
4043 
4044  return fe_values_views_cache.second_order_tensors[tensor.first_tensor_component];
4045 }
4046 
4047 
4048 
4049 
4050 template <int dim, int spacedim>
4051 inline
4052 const double &
4053 FEValuesBase<dim,spacedim>::shape_value (const unsigned int i,
4054  const unsigned int j) const
4055 {
4056  Assert (i < fe->dofs_per_cell,
4057  ExcIndexRange (i, 0, fe->dofs_per_cell));
4059  ExcAccessToUninitializedField("update_values"));
4060  Assert (fe->is_primitive (i),
4062 
4063  // if the entire FE is primitive,
4064  // then we can take a short-cut:
4065  if (fe->is_primitive())
4066  return this->finite_element_output.shape_values(i,j);
4067  else
4068  {
4069  // otherwise, use the mapping
4070  // between shape function
4071  // numbers and rows. note that
4072  // by the assertions above, we
4073  // know that this particular
4074  // shape function is primitive,
4075  // so we can call
4076  // system_to_component_index
4077  const unsigned int
4078  row = this->finite_element_output.shape_function_to_row_table[i * fe->n_components() + fe->system_to_component_index(i).first];
4079  return this->finite_element_output.shape_values(row, j);
4080  }
4081 }
4082 
4083 
4084 
4085 template <int dim, int spacedim>
4086 inline
4087 double
4089  const unsigned int j,
4090  const unsigned int component) const
4091 {
4092  Assert (i < fe->dofs_per_cell,
4093  ExcIndexRange (i, 0, fe->dofs_per_cell));
4094  Assert (this->update_flags & update_values,
4095  ExcAccessToUninitializedField("update_values"));
4096  Assert (component < fe->n_components(),
4097  ExcIndexRange(component, 0, fe->n_components()));
4098 
4099  // check whether the shape function
4100  // is non-zero at all within
4101  // this component:
4102  if (fe->get_nonzero_components(i)[component] == false)
4103  return 0;
4104 
4105  // look up the right row in the
4106  // table and take the data from
4107  // there
4108  const unsigned int
4109  row = this->finite_element_output.shape_function_to_row_table[i * fe->n_components() + component];
4110  return this->finite_element_output.shape_values(row, j);
4111 }
4112 
4113 
4114 
4115 template <int dim, int spacedim>
4116 inline
4117 const Tensor<1,spacedim> &
4118 FEValuesBase<dim,spacedim>::shape_grad (const unsigned int i,
4119  const unsigned int j) const
4120 {
4121  Assert (i < fe->dofs_per_cell,
4122  ExcIndexRange (i, 0, fe->dofs_per_cell));
4124  ExcAccessToUninitializedField("update_gradients"));
4125  Assert (fe->is_primitive (i),
4127 
4128  // if the entire FE is primitive,
4129  // then we can take a short-cut:
4130  if (fe->is_primitive())
4131  return this->finite_element_output.shape_gradients[i][j];
4132  else
4133  {
4134  // otherwise, use the mapping
4135  // between shape function
4136  // numbers and rows. note that
4137  // by the assertions above, we
4138  // know that this particular
4139  // shape function is primitive,
4140  // so we can call
4141  // system_to_component_index
4142  const unsigned int
4143  row = this->finite_element_output.shape_function_to_row_table[i * fe->n_components() + fe->system_to_component_index(i).first];
4144  return this->finite_element_output.shape_gradients[row][j];
4145  }
4146 }
4147 
4148 
4149 
4150 template <int dim, int spacedim>
4151 inline
4154  const unsigned int j,
4155  const unsigned int component) const
4156 {
4157  Assert (i < fe->dofs_per_cell,
4158  ExcIndexRange (i, 0, fe->dofs_per_cell));
4159  Assert (this->update_flags & update_gradients,
4160  ExcAccessToUninitializedField("update_gradients"));
4161  Assert (component < fe->n_components(),
4162  ExcIndexRange(component, 0, fe->n_components()));
4163 
4164  // check whether the shape function
4165  // is non-zero at all within
4166  // this component:
4167  if (fe->get_nonzero_components(i)[component] == false)
4168  return Tensor<1,spacedim>();
4169 
4170  // look up the right row in the
4171  // table and take the data from
4172  // there
4173  const unsigned int
4174  row = this->finite_element_output.shape_function_to_row_table[i * fe->n_components() + component];
4175  return this->finite_element_output.shape_gradients[row][j];
4176 }
4177 
4178 
4179 
4180 template <int dim, int spacedim>
4181 inline
4182 const Tensor<2,spacedim> &
4183 FEValuesBase<dim,spacedim>::shape_hessian (const unsigned int i,
4184  const unsigned int j) const
4185 {
4186  Assert (i < fe->dofs_per_cell,
4187  ExcIndexRange (i, 0, fe->dofs_per_cell));
4188  Assert (this->update_flags & update_hessians,
4189  ExcAccessToUninitializedField("update_hessians"));
4190  Assert (fe->is_primitive (i),
4191  ExcShapeFunctionNotPrimitive(i));
4192 
4193  // if the entire FE is primitive,
4194  // then we can take a short-cut:
4195  if (fe->is_primitive())
4196  return this->finite_element_output.shape_hessians[i][j];
4197  else
4198  {
4199  // otherwise, use the mapping
4200  // between shape function
4201  // numbers and rows. note that
4202  // by the assertions above, we
4203  // know that this particular
4204  // shape function is primitive,
4205  // so we can call
4206  // system_to_component_index
4207  const unsigned int
4208  row = this->finite_element_output.shape_function_to_row_table[i * fe->n_components() + fe->system_to_component_index(i).first];
4209  return this->finite_element_output.shape_hessians[row][j];
4210  }
4211 }
4212 
4213 
4214 
4215 template <int dim, int spacedim>
4216 inline
4219  const unsigned int j,
4220  const unsigned int component) const
4221 {
4222  Assert (i < fe->dofs_per_cell,
4223  ExcIndexRange (i, 0, fe->dofs_per_cell));
4224  Assert (this->update_flags & update_hessians,
4225  ExcAccessToUninitializedField("update_hessians"));
4226  Assert (component < fe->n_components(),
4227  ExcIndexRange(component, 0, fe->n_components()));
4228 
4229  // check whether the shape function
4230  // is non-zero at all within
4231  // this component:
4232  if (fe->get_nonzero_components(i)[component] == false)
4233  return Tensor<2,spacedim>();
4234 
4235  // look up the right row in the
4236  // table and take the data from
4237  // there
4238  const unsigned int
4239  row = this->finite_element_output.shape_function_to_row_table[i * fe->n_components() + component];
4240  return this->finite_element_output.shape_hessians[row][j];
4241 }
4242 
4243 
4244 
4245 template <int dim, int spacedim>
4246 inline
4247 const Tensor<3,spacedim> &
4249  const unsigned int j) const
4250 {
4251  Assert (i < fe->dofs_per_cell,
4252  ExcIndexRange (i, 0, fe->dofs_per_cell));
4253  Assert (this->update_flags & update_hessians,
4254  ExcAccessToUninitializedField("update_3rd_derivatives"));
4255  Assert (fe->is_primitive (i),
4256  ExcShapeFunctionNotPrimitive(i));
4257 
4258  // if the entire FE is primitive,
4259  // then we can take a short-cut:
4260  if (fe->is_primitive())
4261  return this->finite_element_output.shape_3rd_derivatives[i][j];
4262  else
4263  {
4264  // otherwise, use the mapping
4265  // between shape function
4266  // numbers and rows. note that
4267  // by the assertions above, we
4268  // know that this particular
4269  // shape function is primitive,
4270  // so we can call
4271  // system_to_component_index
4272  const unsigned int
4273  row = this->finite_element_output.shape_function_to_row_table[i * fe->n_components() + fe->system_to_component_index(i).first];
4274  return this->finite_element_output.shape_3rd_derivatives[row][j];
4275  }
4276 }
4277 
4278 
4279 
4280 template <int dim, int spacedim>
4281 inline
4284  const unsigned int j,
4285  const unsigned int component) const
4286 {
4287  Assert (i < fe->dofs_per_cell,
4288  ExcIndexRange (i, 0, fe->dofs_per_cell));
4289  Assert (this->update_flags & update_hessians,
4290  ExcAccessToUninitializedField("update_3rd_derivatives"));
4291  Assert (component < fe->n_components(),
4292  ExcIndexRange(component, 0, fe->n_components()));
4293 
4294  // check whether the shape function
4295  // is non-zero at all within
4296  // this component:
4297  if (fe->get_nonzero_components(i)[component] == false)
4298  return Tensor<3,spacedim>();
4299 
4300  // look up the right row in the
4301  // table and take the data from
4302  // there
4303  const unsigned int
4304  row = this->finite_element_output.shape_function_to_row_table[i * fe->n_components() + component];
4305  return this->finite_element_output.shape_3rd_derivatives[row][j];
4306 }
4307 
4308 
4309 
4310 template <int dim, int spacedim>
4311 inline
4314 {
4315  return *fe;
4316 }
4317 
4318 
4319 template <int dim, int spacedim>
4320 inline
4321 const Mapping<dim,spacedim> &
4323 {
4324  return *mapping;
4325 }
4326 
4327 
4328 
4329 template <int dim, int spacedim>
4330 inline
4333 {
4334  return this->update_flags;
4335 }
4336 
4337 
4338 
4339 template <int dim, int spacedim>
4340 inline
4341 const std::vector<Point<spacedim> > &
4343 {
4344  Assert (this->update_flags & update_quadrature_points,
4345  ExcAccessToUninitializedField("update_quadrature_points"));
4346  return this->mapping_output.quadrature_points;
4347 }
4348 
4349 
4350 
4351 template <int dim, int spacedim>
4352 inline
4353 const std::vector<double> &
4355 {
4356  Assert (this->update_flags & update_JxW_values,
4357  ExcAccessToUninitializedField("update_JxW_values"));
4358  return this->mapping_output.JxW_values;
4359 }
4360 
4361 
4362 
4363 template <int dim, int spacedim>
4364 inline
4365 const std::vector<DerivativeForm<1,dim,spacedim> > &
4367 {
4368  Assert (this->update_flags & update_jacobians,
4369  ExcAccessToUninitializedField("update_jacobians"));
4370  return this->mapping_output.jacobians;
4371 }
4372 
4373 
4374 
4375 template <int dim, int spacedim>
4376 inline
4377 const std::vector<DerivativeForm<2,dim,spacedim> > &
4379 {
4380  Assert (this->update_flags & update_jacobian_grads,
4381  ExcAccessToUninitializedField("update_jacobians_grads"));
4382  return this->mapping_output.jacobian_grads;
4383 }
4384 
4385 
4386 
4387 template <int dim, int spacedim>
4388 inline
4389 const Tensor<3,spacedim> &
4391 {
4392  Assert (this->update_flags & update_jacobian_pushed_forward_grads,
4393  ExcAccessToUninitializedField("update_jacobian_pushed_forward_grads"));
4394  return this->mapping_output.jacobian_pushed_forward_grads[i];
4395 }
4396 
4397 
4398 
4399 template <int dim, int spacedim>
4400 inline
4401 const std::vector<Tensor<3,spacedim> > &
4403 {
4404  Assert (this->update_flags & update_jacobian_pushed_forward_grads,
4405  ExcAccessToUninitializedField("update_jacobian_pushed_forward_grads"));
4406  return this->mapping_output.jacobian_pushed_forward_grads;
4407 }
4408 
4409 
4410 
4411 template <int dim, int spacedim>
4412 inline
4414 FEValuesBase<dim,spacedim>::jacobian_2nd_derivative (const unsigned int i) const
4415 {
4416  Assert (this->update_flags & update_jacobian_2nd_derivatives,
4417  ExcAccessToUninitializedField("update_jacobian_2nd_derivatives"));
4418  return this->mapping_output.jacobian_2nd_derivatives[i];
4419 }
4420 
4421 
4422 
4423 template <int dim, int spacedim>
4424 inline
4425 const std::vector<DerivativeForm<3,dim,spacedim> > &
4427 {
4428  Assert (this->update_flags & update_jacobian_2nd_derivatives,
4429  ExcAccessToUninitializedField("update_jacobian_2nd_derivatives"));
4430  return this->mapping_output.jacobian_2nd_derivatives;
4431 }
4432 
4433 
4434 
4435 template <int dim, int spacedim>
4436 inline
4437 const Tensor<4,spacedim> &
4439 {
4441  ExcAccessToUninitializedField("update_jacobian_pushed_forward_2nd_derivatives"));
4442  return this->mapping_output.jacobian_pushed_forward_2nd_derivatives[i];
4443 }
4444 
4445 
4446 
4447 template <int dim, int spacedim>
4448 inline
4449 const std::vector<Tensor<4,spacedim> > &
4451 {
4453  ExcAccessToUninitializedField("update_jacobian_pushed_forward_2nd_derivatives"));
4454  return this->mapping_output.jacobian_pushed_forward_2nd_derivatives;
4455 }
4456 
4457 
4458 
4459 template <int dim, int spacedim>
4460 inline
4462 FEValuesBase<dim,spacedim>::jacobian_3rd_derivative (const unsigned int i) const
4463 {
4464  Assert (this->update_flags & update_jacobian_3rd_derivatives,
4465  ExcAccessToUninitializedField("update_jacobian_3rd_derivatives"));
4466  return this->mapping_output.jacobian_3rd_derivatives[i];
4467 }
4468 
4469 
4470 
4471 template <int dim, int spacedim>
4472 inline
4473 const std::vector<DerivativeForm<4,dim,spacedim> > &
4475 {
4476  Assert (this->update_flags & update_jacobian_3rd_derivatives,
4477  ExcAccessToUninitializedField("update_jacobian_3rd_derivatives"));
4478  return this->mapping_output.jacobian_3rd_derivatives;
4479 }
4480 
4481 
4482 
4483 template <int dim, int spacedim>
4484 inline
4485 const Tensor<5,spacedim> &
4487 {
4489  ExcAccessToUninitializedField("update_jacobian_pushed_forward_3rd_derivatives"));
4490  return this->mapping_output.jacobian_pushed_forward_3rd_derivatives[i];
4491 }
4492 
4493 
4494 
4495 template <int dim, int spacedim>
4496 inline
4497 const std::vector<Tensor<5,spacedim> > &
4499 {
4501  ExcAccessToUninitializedField("update_jacobian_pushed_forward_3rd_derivatives"));
4502  return this->mapping_output.jacobian_pushed_forward_3rd_derivatives;
4503 }
4504 
4505 
4506 template <int dim, int spacedim>
4507 inline
4508 const std::vector<DerivativeForm<1,spacedim,dim> > &
4510 {
4511  Assert (this->update_flags & update_inverse_jacobians,
4512  ExcAccessToUninitializedField("update_inverse_jacobians"));
4513  return this->mapping_output.inverse_jacobians;
4514 }
4515 
4516 
4517 
4518 template <int dim, int spacedim>
4519 inline
4520 const Point<spacedim> &
4521 FEValuesBase<dim,spacedim>::quadrature_point (const unsigned int i) const
4522 {
4523  Assert (this->update_flags & update_quadrature_points,
4524  ExcAccessToUninitializedField("update_quadrature_points"));
4525  Assert (i<this->mapping_output.quadrature_points.size(),
4526  ExcIndexRange(i, 0, this->mapping_output.quadrature_points.size()));
4527 
4528  return this->mapping_output.quadrature_points[i];
4529 }
4530 
4531 
4532 
4533 
4534 template <int dim, int spacedim>
4535 inline
4536 double
4537 FEValuesBase<dim,spacedim>::JxW (const unsigned int i) const
4538 {
4539  Assert (this->update_flags & update_JxW_values,
4540  ExcAccessToUninitializedField("update_JxW_values"));
4541  Assert (i<this->mapping_output.JxW_values.size(),
4542  ExcIndexRange(i, 0, this->mapping_output.JxW_values.size()));
4543 
4544  return this->mapping_output.JxW_values[i];
4545 }
4546 
4547 
4548 
4549 template <int dim, int spacedim>
4550 inline
4552 FEValuesBase<dim,spacedim>::jacobian (const unsigned int i) const
4553 {
4554  Assert (this->update_flags & update_jacobians,
4555  ExcAccessToUninitializedField("update_jacobians"));
4556  Assert (i<this->mapping_output.jacobians.size(),
4557  ExcIndexRange(i, 0, this->mapping_output.jacobians.size()));
4558 
4559  return this->mapping_output.jacobians[i];
4560 }
4561 
4562 
4563 
4564 template <int dim, int spacedim>
4565 inline
4567 FEValuesBase<dim,spacedim>::jacobian_grad (const unsigned int i) const
4568 {
4569  Assert (this->update_flags & update_jacobian_grads,
4570  ExcAccessToUninitializedField("update_jacobians_grads"));
4571  Assert (i<this->mapping_output.jacobian_grads.size(),
4572  ExcIndexRange(i, 0, this->mapping_output.jacobian_grads.size()));
4573 
4574  return this->mapping_output.jacobian_grads[i];
4575 }
4576 
4577 
4578 
4579 template <int dim, int spacedim>
4580 inline
4582 FEValuesBase<dim,spacedim>::inverse_jacobian (const unsigned int i) const
4583 {
4584  Assert (this->update_flags & update_inverse_jacobians,
4585  ExcAccessToUninitializedField("update_inverse_jacobians"));
4586  Assert (i<this->mapping_output.inverse_jacobians.size(),
4587  ExcIndexRange(i, 0, this->mapping_output.inverse_jacobians.size()));
4588 
4589  return this->mapping_output.inverse_jacobians[i];
4590 }
4591 
4592 
4593 template <int dim, int spacedim>
4594 inline
4595 const Tensor<1,spacedim> &
4596 FEValuesBase<dim,spacedim>::normal_vector (const unsigned int i) const
4597 {
4598  Assert (this->update_flags & update_normal_vectors,
4599  (typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_normal_vectors")));
4600  Assert (i<this->mapping_output.normal_vectors.size(),
4601  ExcIndexRange(i, 0, this->mapping_output.normal_vectors.size()));
4602 
4603  return this->mapping_output.normal_vectors[i];
4604 }
4605 
4606 
4607 
4608 /*------------------------ Inline functions: FEValues ----------------------------*/
4609 
4610 
4611 template <int dim, int spacedim>
4612 inline
4613 const Quadrature<dim> &
4615 {
4616  return quadrature;
4617 }
4618 
4619 
4620 
4621 template <int dim, int spacedim>
4622 inline
4623 const FEValues<dim,spacedim> &
4625 {
4626  return *this;
4627 }
4628 
4629 
4630 /*------------------------ Inline functions: FEFaceValuesBase --------------------*/
4631 
4632 
4633 template <int dim, int spacedim>
4634 inline
4635 unsigned int
4637 {
4638  return present_face_index;
4639 }
4640 
4641 
4642 /*------------------------ Inline functions: FE*FaceValues --------------------*/
4643 
4644 template <int dim, int spacedim>
4645 inline
4646 const Quadrature<dim-1> &
4648 {
4649  return quadrature;
4650 }
4651 
4652 
4653 
4654 template <int dim, int spacedim>
4655 inline
4658 {
4659  return *this;
4660 }
4661 
4662 
4663 
4664 template <int dim, int spacedim>
4665 inline
4668 {
4669  return *this;
4670 }
4671 
4672 
4673 
4674 template <int dim, int spacedim>
4675 inline
4676 const Tensor<1,spacedim> &
4677 FEFaceValuesBase<dim,spacedim>::boundary_form (const unsigned int i) const
4678 {
4679  Assert (i<this->mapping_output.boundary_forms.size(),
4680  ExcIndexRange(i, 0, this->mapping_output.boundary_forms.size()));
4681  Assert (this->update_flags & update_boundary_forms,
4682  (typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_boundary_forms")));
4683 
4684  return this->mapping_output.boundary_forms[i];
4685 }
4686 
4687 #endif // DOXYGEN
4688 
4689 DEAL_II_NAMESPACE_CLOSE
4690 
4691 #endif
Transformed quadrature weights.
const std::vector< DerivativeForm< 2, dim, spacedim > > & get_jacobian_grads() const
void get_function_values(const InputVector &fe_function, std::vector< typename InputVector::value_type > &values) const
Definition: fe_values.cc:2689
void get_function_curls(const InputVector &fe_function, std::vector< typename ProductType< curl_type, typename InputVector::value_type >::type > &curls) const
Definition: fe_values.cc:1477
Shape function values.
double shape_value_component(const unsigned int function_no, const unsigned int point_no, const unsigned int component) const
std_cxx11::unique_ptr< typename FiniteElement< dim, spacedim >::InternalDataBase > fe_data
Definition: fe_values.h:2696
const Tensor< 1, spacedim > & shape_grad(const unsigned int function_no, const unsigned int quadrature_point) const
void get_function_laplacians(const InputVector &fe_function, std::vector< typename ProductType< value_type, typename InputVector::value_type >::type > &laplacians) const
Definition: fe_values.cc:1338
const DerivativeForm< 4, dim, spacedim > & jacobian_3rd_derivative(const unsigned int quadrature_point) const
static::ExceptionBase & ExcCannotInitializeField()
void get_function_values(const InputVector &fe_function, std::vector< typename ProductType< value_type, typename InputVector::value_type >::type > &values) const
Definition: fe_values.cc:1269
hessian_type hessian(const unsigned int shape_function, const unsigned int q_point) const
Number value_type
Definition: vector.h:137
const FEValuesViews::Scalar< dim, spacedim > & operator[](const FEValuesExtractors::Scalar &scalar) const
CellSimilarity::Similarity cell_similarity
Definition: fe_values.h:2725
void get_function_third_derivatives(const InputVector &fe_function, std::vector< Tensor< 3, spacedim, typename InputVector::value_type > > &third_derivatives) const
Definition: fe_values.cc:3246
::Tensor< 1, spacedim > divergence_type
Definition: fe_values.h:1126
const Tensor< 3, spacedim > & jacobian_pushed_forward_grad(const unsigned int quadrature_point) const
Cache(const FEValuesBase< dim, spacedim > &fe_values)
Definition: fe_values.cc:1663
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:2955
const Quadrature< dim-1 > & get_quadrature() const
const Tensor< 2, spacedim > & shape_hessian(const unsigned int function_no, const unsigned int point_no) const
void get_function_hessians(const InputVector &fe_function, std::vector< typename ProductType< hessian_type, typename InputVector::value_type >::type > &hessians) const
Definition: fe_values.cc:1315
std::vector< ShapeFunctionData > shape_function_data
Definition: fe_values.h:1086
static::ExceptionBase & ExcAccessToUninitializedField()
const Tensor< 1, spacedim > & normal_vector(const unsigned int i) 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:1315
MappingType
Definition: mapping.h:50
const unsigned int dofs_per_cell
Definition: fe_values.h:1459
const DerivativeForm< 1, dim, spacedim > & jacobian(const unsigned int quadrature_point) const
const unsigned int component
Definition: fe_values.h:391
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)
Definition: fe_values.cc:2254
::Tensor< 3, spacedim > hessian_type
Definition: fe_values.h:483
const Quadrature< dim-1 > quadrature
Definition: fe_values.h:2960
Volume element.
unsigned int get_face_index() const
const FESubfaceValues< dim, spacedim > & get_present_fe_values() const
const Tensor< 4, spacedim > & jacobian_pushed_forward_2nd_derivative(const unsigned int quadrature_point) const
const std::vector< DerivativeForm< 1, spacedim, dim > > & get_inverse_jacobians() const
::internal::CurlType< spacedim >::type curl_type
Definition: fe_values.h:476
Outer normal vector, not normalized.
const FEFaceValues< dim, spacedim > & get_present_fe_values() const
third_derivative_type third_derivative(const unsigned int shape_function, const unsigned int q_point) const
hessian_type hessian(const unsigned int shape_function, const unsigned int q_point) const
static::ExceptionBase & ExcFEDontMatch()
Scalar & operator=(const Scalar< dim, spacedim > &)
Definition: fe_values.cc:148
STL namespace.
Transformed quadrature points.
const Tensor< 1, spacedim > & boundary_form(const unsigned int i) const
const Triangulation< dim, spacedim >::cell_iterator get_cell() const
Definition: fe_values.cc:3364
value_type value(const unsigned int shape_function, const unsigned int q_point) const
::Tensor< 1, spacedim > gradient_type
Definition: fe_values.h:152
::SymmetricTensor< 2, spacedim > value_type
Definition: fe_values.h:911
void transform(std::vector< Tensor< 1, spacedim > > &transformed, const std::vector< Tensor< 1, dim > > &original, MappingType mapping) const 1
Definition: fe_values.cc:3402
unsigned int row_index[spacedim]
Definition: fe_values.h:517
void get_function_hessians(const InputVector &fe_function, std::vector< Tensor< 2, spacedim, typename InputVector::value_type > > &hessians) const
Definition: fe_values.cc:2973
void check_cell_similarity(const typename Triangulation< dim, spacedim >::cell_iterator &cell)
Definition: fe_values.cc:3516
static const unsigned int space_dimension
Definition: fe_values.h:1447
::internal::FEValuesViews::Cache< dim, spacedim > fe_values_views_cache
Definition: fe_values.h:2751
static::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
void get_function_third_derivatives(const InputVector &fe_function, std::vector< typename ProductType< third_derivative_type, typename InputVector::value_type >::type > &third_derivatives) const
Definition: fe_values.cc:1361
static::ExceptionBase & ExcShapeFunctionNotPrimitive(int arg1)
void get_function_symmetric_gradients(const InputVector &fe_function, std::vector< typename ProductType< symmetric_gradient_type, typename InputVector::value_type >::type > &symmetric_gradients) const
Definition: fe_values.cc:1431
void get_function_gradients(const InputVector &fe_function, std::vector< typename ProductType< gradient_type, typename InputVector::value_type >::type > &gradients) const
Definition: fe_values.cc:1408
const std::vector< double > & get_JxW_values() const
UpdateFlags get_update_flags() const
static::ExceptionBase & ExcInvalidUpdateFlag()
const Quadrature< dim > & get_quadrature() const
UpdateFlags compute_update_flags(const UpdateFlags update_flags) const
Definition: fe_values.cc:3435
void get_function_laplacians(const InputVector &fe_function, std::vector< typename InputVector::value_type > &laplacians) const
Definition: fe_values.cc:3091
const SmartPointer< const FiniteElement< dim, spacedim >, FEValuesBase< dim, spacedim > > fe
Definition: fe_values.h:2689
static::ExceptionBase & ExcMessage(std::string arg1)
const std::vector< DerivativeForm< 1, dim, spacedim > > & get_jacobians() const
const std::vector< Point< spacedim > > & get_quadrature_points() const
void get_function_hessians(const InputVector &fe_function, std::vector< typename ProductType< hessian_type, typename InputVector::value_type >::type > &hessians) const
Definition: fe_values.cc:1499
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:439
#define DeclException1(Exception1, type1, outsequence)
Definition: exceptions.h:564
std_cxx11::unique_ptr< const CellIteratorBase > present_cell
Definition: fe_values.h:2620
third_derivative_type third_derivative(const unsigned int shape_function, const unsigned int q_point) const
Third derivatives of shape functions.
void get_function_gradients(const InputVector &fe_function, std::vector< typename ProductType< gradient_type, typename InputVector::value_type >::type > &gradients) const
Definition: fe_values.cc:1292
std::vector< ShapeFunctionData > shape_function_data
Definition: fe_values.h:1291
#define Assert(cond, exc)
Definition: exceptions.h:313
::Tensor< 2, spacedim > value_type
Definition: fe_values.h:1121
UpdateFlags
void get_function_third_derivatives(const InputVector &fe_function, std::vector< typename ProductType< third_derivative_type, typename InputVector::value_type >::type > &third_derivatives) const
Definition: fe_values.cc:1547
gradient_type gradient(const unsigned int shape_function, const unsigned int q_point) const
Abstract base class for mapping classes.
Definition: dof_tools.h:46
std::vector< ShapeFunctionData > shape_function_data
Definition: fe_values.h:870
const Quadrature< dim > quadrature
Definition: fe_values.h:2859
const unsigned int first_vector_component
Definition: fe_values.h:865
std::size_t memory_consumption() const
Definition: fe_values.cc:3415
#define DeclException0(Exception0)
Definition: exceptions.h:541
::internal::FEValues::FiniteElementRelatedData< dim, spacedim > finite_element_output
Definition: fe_values.h:2702
void invalidate_present_cell()
Definition: fe_values.cc:3453
const SmartPointer< const Mapping< dim, spacedim >, FEValuesBase< dim, spacedim > > mapping
Definition: fe_values.h:2669
const std::vector< Tensor< 3, spacedim > > & get_jacobian_pushed_forward_grads() const
Vector & operator=(const Vector< dim, spacedim > &)
Definition: fe_values.cc:243
Tensor()
Definition: tensor.h:745
const DerivativeForm< 3, dim, spacedim > & jacobian_2nd_derivative(const unsigned int quadrature_point) const
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:506
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.
static::ExceptionBase & ExcFENotPrimitive()
static::ExceptionBase & ExcAccessToUninitializedField(char *arg1)
Gradient of volume element.
const Mapping< dim, spacedim > & get_mapping() const
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
static TableIndices< rank_ > unrolled_to_component_indices(const unsigned int i)
Definition: tensor.h:1051
void get_function_values(const InputVector &fe_function, std::vector< typename ProductType< value_type, typename InputVector::value_type >::type > &values) const
Definition: fe_values.cc:1384
const DerivativeForm< 1, spacedim, dim > & inverse_jacobian(const unsigned int quadrature_point) const
const unsigned int n_quadrature_points
Definition: fe_values.h:1452
const double & shape_value(const unsigned int function_no, const unsigned int point_no) const
const SmartPointer< const FEValuesBase< dim, spacedim > > fe_values
Definition: fe_values.h:859
static const unsigned int dimension
Definition: fe_values.h:1442
::Tensor< 2, spacedim > hessian_type
Definition: fe_values.h:159
boost::signals2::connection tria_listener_mesh_transform
Definition: fe_values.h:2645
Tensor< 3, spacedim > shape_3rd_derivative_component(const unsigned int function_no, const unsigned int point_no, const unsigned int component) const
const SmartPointer< const FEValuesBase< dim, spacedim > > fe_values
Definition: fe_values.h:385
Definition: mpi.h:41
::SymmetricTensor< 2, spacedim > symmetric_gradient_type
Definition: fe_values.h:461
const std::vector< Tensor< 4, spacedim > > & get_jacobian_pushed_forward_2nd_derivatives() const
Shape function gradients.
Normal vectors.
void maybe_invalidate_previous_present_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell)
Definition: fe_values.cc:3471
::Tensor< 2, spacedim > gradient_type
Definition: fe_values.h:449
std_cxx11::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > mapping_data
Definition: fe_values.h:2676
void get_function_divergences(const InputVector &fe_function, std::vector< typename ProductType< divergence_type, typename InputVector::value_type >::type > &divergences) const
Definition: fe_values.cc:1455
const SmartPointer< const FEValuesBase< dim, spacedim > > fe_values
Definition: fe_values.h:1280
Definition: fe.h:30
const SmartPointer< const FEValuesBase< dim, spacedim > > fe_values
Definition: fe_values.h:1075
const DerivativeForm< 2, dim, spacedim > & jacobian_grad(const unsigned int quadrature_point) const
const std::vector< Tensor< 1, spacedim > > & get_all_normal_vectors() const
Definition: fe_values.cc:3373
value_type value(const unsigned int shape_function, const unsigned int q_point) const
void get_function_laplacians(const InputVector &fe_function, std::vector< typename ProductType< value_type, typename InputVector::value_type >::type > &laplacians) const
Definition: fe_values.cc:1522
static::ExceptionBase & ExcNotImplemented()
const FiniteElement< dim, spacedim > & get_fe() const
boost::signals2::connection tria_listener_refinement
Definition: fe_values.h:2636
::internal::FEValues::MappingRelatedData< dim, spacedim > mapping_output
Definition: fe_values.h:2682
const std::vector< DerivativeForm< 3, dim, spacedim > > & get_jacobian_2nd_derivatives() const
unsigned int n_components(const DoFHandler< dim, spacedim > &dh)
const Point< spacedim > & quadrature_point(const unsigned int q) const
std::vector< Point< spacedim > > get_normal_vectors() const 1
Definition: fe_values.cc:3384
const std::vector< DerivativeForm< 4, dim, spacedim > > & get_jacobian_3rd_derivatives() const
const Tensor< 3, spacedim > & shape_3rd_derivative(const unsigned int function_no, const unsigned int point_no) const
::Tensor< 3, spacedim > third_derivative_type
Definition: fe_values.h:166
const std::vector< Tensor< 5, spacedim > > & get_jacobian_pushed_forward_3rd_derivatives() const
const Tensor< 5, spacedim > & jacobian_pushed_forward_3rd_derivative(const unsigned int quadrature_point) const
double JxW(const unsigned int quadrature_point) const
const FEValues< dim, spacedim > & get_present_fe_values() const
CellSimilarity::Similarity get_cell_similarity() const
Definition: fe_values.cc:3570
std::vector< ShapeFunctionData > shape_function_data
Definition: fe_values.h:396
UpdateFlags update_flags
Definition: fe_values.h:2708
static::ExceptionBase & ExcInternalError()
void get_function_gradients(const InputVector &fe_function, std::vector< Tensor< 1, spacedim, typename InputVector::value_type > > &gradients) const
Definition: fe_values.cc:2850
::Tensor< 4, spacedim > third_derivative_type
Definition: fe_values.h:490