Reference documentation for deal.II version 8.4.1
fe_values.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1998 - 2016 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> class Quadrature;
52 template <int dim, int spacedim=dim> class FEValuesBase;
53 
54 template <typename Number> class Vector;
55 template <typename Number> class BlockVector;
56 
57 
58 namespace internal
59 {
64  template <int dim>
65  struct CurlType;
66 
73  template <>
74  struct CurlType<1>
75  {
76  typedef Tensor<1,1> type;
77  };
78 
85  template <>
86  struct CurlType<2>
87  {
88  typedef Tensor<1,1> type;
89  };
90 
97  template <>
98  struct CurlType<3>
99  {
100  typedef Tensor<1,3> type;
101  };
102 }
103 
104 
105 
106 
128 namespace FEValuesViews
129 {
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  typedef ::Tensor<3,spacedim> third_derivative_type;
172 
178  {
188 
197  unsigned int row_index;
198  };
199 
203  Scalar ();
204 
210  Scalar (const FEValuesBase<dim,spacedim> &fe_values_base,
211  const unsigned int component);
212 
218 
232  value_type
233  value (const unsigned int shape_function,
234  const unsigned int q_point) const;
235 
246  gradient_type
247  gradient (const unsigned int shape_function,
248  const unsigned int q_point) const;
249 
260  hessian_type
261  hessian (const unsigned int shape_function,
262  const unsigned int q_point) const;
263 
274  third_derivative_type
275  third_derivative (const unsigned int shape_function,
276  const unsigned int q_point) const;
277 
295  template <class InputVector>
296  void get_function_values (const InputVector &fe_function,
297  std::vector<typename ProductType<value_type,typename InputVector::value_type>::type> &values) const;
298 
316  template <class InputVector>
317  void get_function_gradients (const InputVector &fe_function,
318  std::vector<typename ProductType<gradient_type,typename InputVector::value_type>::type> &gradients) const;
319 
337  template <class InputVector>
338  void get_function_hessians (const InputVector &fe_function,
339  std::vector<typename ProductType<hessian_type,typename InputVector::value_type>::type> &hessians) const;
340 
359  template <class InputVector>
360  void get_function_laplacians (const InputVector &fe_function,
361  std::vector<typename ProductType<value_type,typename InputVector::value_type>::type> &laplacians) const;
362 
381  template <class InputVector>
382  void get_function_third_derivatives (const InputVector &fe_function,
383  std::vector<typename ProductType<third_derivative_type,
384  typename InputVector::value_type>::type> &third_derivatives) const;
385 
386  private:
391 
396  const unsigned int component;
397 
401  std::vector<ShapeFunctionData> shape_function_data;
402  };
403 
404 
405 
435  template <int dim, int spacedim=dim>
436  class Vector
437  {
438  public:
444  typedef ::Tensor<1,spacedim> value_type;
445 
454  typedef ::Tensor<2,spacedim> gradient_type;
455 
466  typedef ::SymmetricTensor<2,spacedim> symmetric_gradient_type;
467 
473  typedef double divergence_type;
474 
481  typedef typename ::internal::CurlType<spacedim>::type curl_type;
482 
488  typedef ::Tensor<3,spacedim> hessian_type;
489 
495  typedef ::Tensor<4,spacedim> third_derivative_type;
496 
502  {
512 
522  unsigned int row_index[spacedim];
523 
533  unsigned int single_nonzero_component_index;
534  };
535 
539  Vector ();
540 
549  Vector (const FEValuesBase<dim,spacedim> &fe_values_base,
550  const unsigned int first_vector_component);
551 
557 
574  value_type
575  value (const unsigned int shape_function,
576  const unsigned int q_point) const;
577 
591  gradient_type
592  gradient (const unsigned int shape_function,
593  const unsigned int q_point) const;
594 
610  symmetric_gradient_type
611  symmetric_gradient (const unsigned int shape_function,
612  const unsigned int q_point) const;
613 
624  divergence_type
625  divergence (const unsigned int shape_function,
626  const unsigned int q_point) const;
627 
648  curl_type
649  curl (const unsigned int shape_function,
650  const unsigned int q_point) const;
651 
662  hessian_type
663  hessian (const unsigned int shape_function,
664  const unsigned int q_point) const;
665 
676  third_derivative_type
677  third_derivative (const unsigned int shape_function,
678  const unsigned int q_point) const;
679 
697  template <class InputVector>
698  void get_function_values (const InputVector &fe_function,
699  std::vector<typename ProductType<value_type,typename InputVector::value_type>::type> &values) const;
700 
718  template <class InputVector>
719  void get_function_gradients (const InputVector &fe_function,
720  std::vector<typename ProductType<gradient_type,typename InputVector::value_type>::type> &gradients) const;
721 
745  template <class InputVector>
746  void
747  get_function_symmetric_gradients (const InputVector &fe_function,
748  std::vector<typename ProductType<symmetric_gradient_type,typename InputVector::value_type>::type> &symmetric_gradients) const;
749 
768  template <class InputVector>
769  void get_function_divergences (const InputVector &fe_function,
770  std::vector<typename ProductType<divergence_type,typename InputVector::value_type>::type> &divergences) const;
771 
790  template <class InputVector>
791  void get_function_curls (const InputVector &fe_function,
792  std::vector<typename ProductType<curl_type,typename InputVector::value_type>::type> &curls) const;
793 
811  template <class InputVector>
812  void get_function_hessians (const InputVector &fe_function,
813  std::vector<typename ProductType<hessian_type,typename InputVector::value_type>::type> &hessians) const;
814 
833  template <class InputVector>
834  void get_function_laplacians (const InputVector &fe_function,
835  std::vector<typename ProductType<value_type,typename InputVector::value_type>::type> &laplacians) const;
836 
855  template <class InputVector>
856  void get_function_third_derivatives (const InputVector &fe_function,
857  std::vector<typename ProductType<third_derivative_type,
858  typename InputVector::value_type>::type> &third_derivatives) const;
859 
860  private:
865 
870  const unsigned int first_vector_component;
871 
875  std::vector<ShapeFunctionData> shape_function_data;
876  };
877 
878 
879  template <int rank, int dim, int spacedim = dim>
880  class SymmetricTensor;
881 
906  template <int dim, int spacedim>
907  class SymmetricTensor<2,dim,spacedim>
908  {
909  public:
916  typedef ::SymmetricTensor<2, spacedim> value_type;
917 
927  typedef ::Tensor<1, spacedim> divergence_type;
928 
933  struct ShapeFunctionData
934  {
943  bool is_nonzero_shape_function_component[value_type::n_independent_components];
944 
954  unsigned int row_index[value_type::n_independent_components];
955 
965  unsigned int single_nonzero_component_index;
966  };
967 
971  SymmetricTensor();
972 
982  SymmetricTensor(const FEValuesBase<dim, spacedim> &fe_values_base,
983  const unsigned int first_tensor_component);
984 
989  SymmetricTensor &operator=(const SymmetricTensor<2, dim, spacedim> &);
990 
1008  value_type
1009  value (const unsigned int shape_function,
1010  const unsigned int q_point) const;
1011 
1012 
1026  divergence_type
1027  divergence (const unsigned int shape_function,
1028  const unsigned int q_point) const;
1029 
1047  template <class InputVector>
1048  void get_function_values (const InputVector &fe_function,
1049  std::vector<typename ProductType<value_type,typename InputVector::value_type>::type> &values) const;
1050 
1072  template <class InputVector>
1073  void get_function_divergences (const InputVector &fe_function,
1074  std::vector<typename ProductType<divergence_type,typename InputVector::value_type>::type> &divergences) const;
1075 
1076  private:
1081 
1086  const unsigned int first_tensor_component;
1087 
1091  std::vector<ShapeFunctionData> shape_function_data;
1092  };
1093 
1094 
1095  template <int rank, int dim, int spacedim = dim>
1096  class Tensor;
1097 
1117  template <int dim, int spacedim>
1118  class Tensor<2,dim,spacedim>
1119  {
1120  public:
1121 
1126  typedef ::Tensor<2, spacedim> value_type;
1127 
1131  typedef ::Tensor<1, spacedim> divergence_type;
1132 
1137  struct ShapeFunctionData
1138  {
1147  bool is_nonzero_shape_function_component[value_type::n_independent_components];
1148 
1158  unsigned int row_index[value_type::n_independent_components];
1159 
1169  unsigned int single_nonzero_component_index;
1170  };
1171 
1175  Tensor();
1176 
1177 
1187  Tensor(const FEValuesBase<dim, spacedim> &fe_values_base,
1188  const unsigned int first_tensor_component);
1189 
1190 
1195  Tensor &operator=(const Tensor<2, dim, spacedim> &);
1196 
1213  value_type
1214  value (const unsigned int shape_function,
1215  const unsigned int q_point) const;
1216 
1230  divergence_type
1231  divergence (const unsigned int shape_function,
1232  const unsigned int q_point) const;
1233 
1251  template <class InputVector>
1252  void get_function_values (const InputVector &fe_function,
1253  std::vector<typename ProductType<value_type,typename InputVector::value_type>::type> &values) const;
1254 
1255 
1277  template <class InputVector>
1278  void get_function_divergences (const InputVector &fe_function,
1279  std::vector<typename ProductType<divergence_type,typename InputVector::value_type>::type> &divergences) const;
1280 
1281  private:
1286 
1291  const unsigned int first_tensor_component;
1292 
1296  std::vector<ShapeFunctionData> shape_function_data;
1297  };
1298 
1299 }
1300 
1301 
1302 namespace internal
1303 {
1304  namespace FEValuesViews
1305  {
1313  template <int dim, int spacedim>
1314  struct Cache
1315  {
1320  std::vector<::FEValuesViews::Scalar<dim,spacedim> > scalars;
1321  std::vector<::FEValuesViews::Vector<dim,spacedim> > vectors;
1322  std::vector<::FEValuesViews::SymmetricTensor<2,dim,spacedim> >
1323  symmetric_second_order_tensors;
1324  std::vector<::FEValuesViews::Tensor<2,dim,spacedim> >
1325  second_order_tensors;
1326 
1330  Cache (const FEValuesBase<dim,spacedim> &fe_values);
1331  };
1332  }
1333 }
1334 
1335 
1336 
1439 template <int dim, int spacedim>
1440 class FEValuesBase :
1441  public Subscriptor
1442 {
1443 public:
1447  static const unsigned int dimension = dim;
1448 
1452  static const unsigned int space_dimension = spacedim;
1453 
1457  const unsigned int n_quadrature_points;
1458 
1464  const unsigned int dofs_per_cell;
1465 
1466 
1474  FEValuesBase (const unsigned int n_q_points,
1475  const unsigned int dofs_per_cell,
1476  const UpdateFlags update_flags,
1479 
1480 
1484  ~FEValuesBase ();
1485 
1486 
1488 
1489 
1510  const double &shape_value (const unsigned int function_no,
1511  const unsigned int point_no) const;
1512 
1533  double shape_value_component (const unsigned int function_no,
1534  const unsigned int point_no,
1535  const unsigned int component) const;
1536 
1562  const Tensor<1,spacedim> &
1563  shape_grad (const unsigned int function_no,
1564  const unsigned int quadrature_point) const;
1565 
1583  shape_grad_component (const unsigned int function_no,
1584  const unsigned int point_no,
1585  const unsigned int component) const;
1586 
1606  const Tensor<2,spacedim> &
1607  shape_hessian (const unsigned int function_no,
1608  const unsigned int point_no) const;
1609 
1627  shape_hessian_component (const unsigned int function_no,
1628  const unsigned int point_no,
1629  const unsigned int component) const;
1630 
1650  const Tensor<3,spacedim> &
1651  shape_3rd_derivative (const unsigned int function_no,
1652  const unsigned int point_no) const;
1653 
1671  shape_3rd_derivative_component (const unsigned int function_no,
1672  const unsigned int point_no,
1673  const unsigned int component) const;
1674 
1676 
1678 
1717  template <class InputVector>
1718  void get_function_values (const InputVector &fe_function,
1719  std::vector<typename InputVector::value_type> &values) const;
1720 
1734  template <class InputVector>
1735  void get_function_values (const InputVector &fe_function,
1736  std::vector<Vector<typename InputVector::value_type> > &values) const;
1737 
1756  template <class InputVector>
1757  void get_function_values (const InputVector &fe_function,
1758  const VectorSlice<const std::vector<types::global_dof_index> > &indices,
1759  std::vector<typename InputVector::value_type> &values) const;
1760 
1782  template <class InputVector>
1783  void get_function_values (const InputVector &fe_function,
1784  const VectorSlice<const std::vector<types::global_dof_index> > &indices,
1785  std::vector<Vector<typename InputVector::value_type> > &values) const;
1786 
1787 
1818  template <class InputVector>
1819  void get_function_values (const InputVector &fe_function,
1820  const VectorSlice<const std::vector<types::global_dof_index> > &indices,
1821  VectorSlice<std::vector<std::vector<typename InputVector::value_type> > > values,
1822  const bool quadrature_points_fastest) const;
1823 
1825 
1827 
1866  template <class InputVector>
1867  void get_function_gradients (const InputVector &fe_function,
1868  std::vector<Tensor<1,spacedim,typename InputVector::value_type> > &gradients) const;
1869 
1886  template <class InputVector>
1887  void get_function_gradients (const InputVector &fe_function,
1888  std::vector<std::vector<Tensor<1,spacedim,typename InputVector::value_type> > > &gradients) const;
1889 
1896  template <class InputVector>
1897  void get_function_gradients (const InputVector &fe_function,
1898  const VectorSlice<const std::vector<types::global_dof_index> > &indices,
1899  std::vector<Tensor<1,spacedim,typename InputVector::value_type> > &gradients) const;
1900 
1907  template <class InputVector>
1908  void get_function_gradients (const InputVector &fe_function,
1909  const VectorSlice<const std::vector<types::global_dof_index> > &indices,
1910  VectorSlice<std::vector<std::vector<Tensor<1,spacedim,typename InputVector::value_type> > > > gradients,
1911  bool quadrature_points_fastest = false) const;
1912 
1914 
1916 
1956  template <class InputVector>
1957  void
1958  get_function_hessians (const InputVector &fe_function,
1959  std::vector<Tensor<2,spacedim,typename InputVector::value_type> > &hessians) const;
1960 
1978  template <class InputVector>
1979  void
1980  get_function_hessians (const InputVector &fe_function,
1981  std::vector<std::vector<Tensor<2,spacedim,typename InputVector::value_type> > > &hessians,
1982  bool quadrature_points_fastest = false) const;
1983 
1988  template <class InputVector>
1989  void get_function_hessians (
1990  const InputVector &fe_function,
1991  const VectorSlice<const std::vector<types::global_dof_index> > &indices,
1992  std::vector<Tensor<2,spacedim,typename InputVector::value_type> > &hessians) const;
1993 
2000  template <class InputVector>
2001  void get_function_hessians (
2002  const InputVector &fe_function,
2003  const VectorSlice<const std::vector<types::global_dof_index> > &indices,
2004  VectorSlice<std::vector<std::vector<Tensor<2,spacedim,typename InputVector::value_type> > > > hessians,
2005  bool quadrature_points_fastest = false) const;
2006 
2051  template <class InputVector>
2052  void
2053  get_function_laplacians (const InputVector &fe_function,
2054  std::vector<typename InputVector::value_type> &laplacians) const;
2055 
2075  template <class InputVector>
2076  void
2077  get_function_laplacians (const InputVector &fe_function,
2078  std::vector<Vector<typename InputVector::value_type> > &laplacians) const;
2079 
2086  template <class InputVector>
2088  const InputVector &fe_function,
2089  const VectorSlice<const std::vector<types::global_dof_index> > &indices,
2090  std::vector<typename InputVector::value_type> &laplacians) const;
2091 
2098  template <class InputVector>
2100  const InputVector &fe_function,
2101  const VectorSlice<const std::vector<types::global_dof_index> > &indices,
2102  std::vector<Vector<typename InputVector::value_type> > &laplacians) const;
2103 
2110  template <class InputVector>
2112  const InputVector &fe_function,
2113  const VectorSlice<const std::vector<types::global_dof_index> > &indices,
2114  std::vector<std::vector<typename InputVector::value_type> > &laplacians,
2115  bool quadrature_points_fastest = false) const;
2116 
2118 
2120 
2161  template <class InputVector>
2162  void
2163  get_function_third_derivatives (const InputVector &fe_function,
2164  std::vector<Tensor<3,spacedim,typename InputVector::value_type> > &third_derivatives) const;
2165 
2184  template <class InputVector>
2185  void
2186  get_function_third_derivatives (const InputVector &fe_function,
2187  std::vector<std::vector<Tensor<3,spacedim,typename InputVector::value_type> > > &third_derivatives,
2188  bool quadrature_points_fastest = false) const;
2189 
2194  template <class InputVector>
2196  const InputVector &fe_function,
2197  const VectorSlice<const std::vector<types::global_dof_index> > &indices,
2198  std::vector<Tensor<3,spacedim,typename InputVector::value_type> > &third_derivatives) const;
2199 
2206  template <class InputVector>
2208  const InputVector &fe_function,
2209  const VectorSlice<const std::vector<types::global_dof_index> > &indices,
2210  VectorSlice<std::vector<std::vector<Tensor<3,spacedim,typename InputVector::value_type> > > > third_derivatives,
2211  bool quadrature_points_fastest = false) const;
2213 
2215 
2216 
2222  const Point<spacedim> &
2223  quadrature_point (const unsigned int q) const;
2224 
2230  const std::vector<Point<spacedim> > &get_quadrature_points () const;
2231 
2247  double JxW (const unsigned int quadrature_point) const;
2248 
2252  const std::vector<double> &get_JxW_values () const;
2253 
2260  const DerivativeForm<1,dim,spacedim> &jacobian (const unsigned int quadrature_point) const;
2261 
2268  const std::vector<DerivativeForm<1,dim,spacedim> > &get_jacobians () const;
2269 
2277  const DerivativeForm<2,dim,spacedim> &jacobian_grad (const unsigned int quadrature_point) const;
2278 
2285  const std::vector<DerivativeForm<2,dim,spacedim> > &get_jacobian_grads () const;
2286 
2295  const Tensor<3,spacedim> &jacobian_pushed_forward_grad (const unsigned int quadrature_point) const;
2296 
2303  const std::vector<Tensor<3,spacedim> > &get_jacobian_pushed_forward_grads () const;
2304 
2313 
2320  const std::vector<DerivativeForm<3,dim,spacedim> > &get_jacobian_2nd_derivatives () const;
2321 
2332 
2339  const std::vector<Tensor<4,spacedim> > &get_jacobian_pushed_forward_2nd_derivatives () const;
2340 
2350 
2357  const std::vector<DerivativeForm<4,dim,spacedim> > &get_jacobian_3rd_derivatives () const;
2358 
2369 
2376  const std::vector<Tensor<5,spacedim> > &get_jacobian_pushed_forward_3rd_derivatives () const;
2377 
2384  const DerivativeForm<1,spacedim,dim> &inverse_jacobian (const unsigned int quadrature_point) const;
2385 
2392  const std::vector<DerivativeForm<1,spacedim,dim> > &get_inverse_jacobians () const;
2393 
2407  const Tensor<1,spacedim> &normal_vector (const unsigned int i) const;
2408 
2424  const std::vector<Tensor<1,spacedim> > &get_all_normal_vectors () const;
2425 
2436  std::vector<Point<spacedim> > get_normal_vectors () const DEAL_II_DEPRECATED;
2437 
2444  void transform (std::vector<Tensor<1,spacedim> > &transformed,
2445  const std::vector<Tensor<1,dim> > &original,
2447 
2449 
2451 
2452 
2461  const FEValuesViews::Scalar<dim,spacedim> &
2462  operator[] (const FEValuesExtractors::Scalar &scalar) const;
2463 
2472  const FEValuesViews::Vector<dim,spacedim> &
2473  operator[] (const FEValuesExtractors::Vector &vector) const;
2474 
2484  const FEValuesViews::SymmetricTensor<2,dim,spacedim> &
2485  operator[] (const FEValuesExtractors::SymmetricTensor<2> &tensor) const;
2486 
2487 
2496  const FEValuesViews::Tensor<2,dim,spacedim> &
2497  operator[] (const FEValuesExtractors::Tensor<2> &tensor) const;
2498 
2500 
2502 
2503 
2507  const Mapping<dim,spacedim> &get_mapping () const;
2508 
2512  const FiniteElement<dim,spacedim> &get_fe () const;
2513 
2517  UpdateFlags get_update_flags () const;
2518 
2522  const typename Triangulation<dim,spacedim>::cell_iterator get_cell () const;
2523 
2529  CellSimilarity::Similarity get_cell_similarity () const;
2530 
2535  std::size_t memory_consumption () const;
2537 
2538 
2545  DeclException1 (ExcAccessToUninitializedField,
2546  char *,
2547  << ("You are requesting information from an FEValues/FEFaceValues/FESubfaceValues "
2548  "object for which this kind of information has not been computed. What "
2549  "information these objects compute is determined by the update_* flags you "
2550  "pass to the constructor. Here, the operation you are attempting requires "
2551  "the <")
2552  << arg1
2553  << "> flag to be set, but it was apparently not specified upon construction.");
2559  DeclException0 (ExcCannotInitializeField);
2565  DeclException0 (ExcInvalidUpdateFlag);
2571  DeclException0 (ExcFEDontMatch);
2577  DeclException1 (ExcShapeFunctionNotPrimitive,
2578  int,
2579  << "The shape function with index " << arg1
2580  << " is not primitive, i.e. it is vector-valued and "
2581  << "has more than one non-zero vector component. This "
2582  << "function cannot be called for these shape functions. "
2583  << "Maybe you want to use the same function with the "
2584  << "_component suffix?");
2590  DeclException0 (ExcFENotPrimitive);
2591 
2592 protected:
2621  class CellIteratorBase;
2622 
2627  template <typename CI> class CellIterator;
2629 
2635  std_cxx11::unique_ptr<const CellIteratorBase> present_cell;
2636 
2644  boost::signals2::connection tria_listener;
2645 
2651  void invalidate_present_cell ();
2652 
2662  void
2663  maybe_invalidate_previous_present_cell (const typename Triangulation<dim,spacedim>::cell_iterator &cell);
2664 
2668  const SmartPointer<const Mapping<dim,spacedim>,FEValuesBase<dim,spacedim> > mapping;
2669 
2675  std_cxx11::unique_ptr<typename Mapping<dim,spacedim>::InternalDataBase> mapping_data;
2676 
2681  ::internal::FEValues::MappingRelatedData<dim, spacedim> mapping_output;
2682 
2683 
2688  const SmartPointer<const FiniteElement<dim,spacedim>,FEValuesBase<dim,spacedim> > fe;
2689 
2695  std_cxx11::unique_ptr<typename FiniteElement<dim,spacedim>::InternalDataBase> fe_data;
2696 
2701  ::internal::FEValues::FiniteElementRelatedData<dim, spacedim> finite_element_output;
2702 
2703 
2708 
2717  UpdateFlags compute_update_flags (const UpdateFlags update_flags) const;
2718 
2725 
2731  void
2732  check_cell_similarity (const typename Triangulation<dim,spacedim>::cell_iterator &cell);
2733 
2734 private:
2739  FEValuesBase (const FEValuesBase &);
2740 
2745  FEValuesBase &operator= (const FEValuesBase &);
2746 
2751 
2756  template <int, int> friend class FEValuesViews::Scalar;
2757  template <int, int> friend class FEValuesViews::Vector;
2758  template <int, int, int> friend class FEValuesViews::SymmetricTensor;
2759  template <int, int, int> friend class FEValuesViews::Tensor;
2760 };
2761 
2762 
2763 
2774 template <int dim, int spacedim=dim>
2775 class FEValues : public FEValuesBase<dim,spacedim>
2776 {
2777 public:
2782  static const unsigned int integral_dimension = dim;
2783 
2788  FEValues (const Mapping<dim,spacedim> &mapping,
2789  const FiniteElement<dim,spacedim> &fe,
2790  const Quadrature<dim> &quadrature,
2791  const UpdateFlags update_flags);
2792 
2799  const Quadrature<dim> &quadrature,
2800  const UpdateFlags update_flags);
2801 
2808  template <template <int, int> class DoFHandlerType, bool level_dof_access>
2809  void reinit (const TriaIterator<DoFCellAccessor<DoFHandlerType<dim,spacedim>,level_dof_access> > &cell);
2810 
2824  void reinit (const typename Triangulation<dim,spacedim>::cell_iterator &cell);
2825 
2830  const Quadrature<dim> &get_quadrature () const;
2831 
2836  std::size_t memory_consumption () const;
2837 
2851  const FEValues<dim,spacedim> &get_present_fe_values () const;
2852 
2853 private:
2858 
2862  void initialize (const UpdateFlags update_flags);
2863 
2870  void do_reinit ();
2871 };
2872 
2873 
2884 template <int dim, int spacedim=dim>
2885 class FEFaceValuesBase : public FEValuesBase<dim,spacedim>
2886 {
2887 public:
2892  static const unsigned int integral_dimension = dim-1;
2893 
2905  FEFaceValuesBase (const unsigned int n_q_points,
2906  const unsigned int dofs_per_cell,
2907  const UpdateFlags update_flags,
2910  const Quadrature<dim-1>& quadrature);
2911 
2919  const Tensor<1,spacedim> &boundary_form (const unsigned int i) const;
2920 
2927  const std::vector<Tensor<1,spacedim> > &get_boundary_forms () const;
2928 
2933  unsigned int get_face_index() const;
2934 
2939  const Quadrature<dim-1> & get_quadrature () const;
2940 
2945  std::size_t memory_consumption () const;
2946 
2947 protected:
2948 
2953  unsigned int present_face_index;
2954 
2958  const Quadrature<dim-1> quadrature;
2959 };
2960 
2961 
2962 
2977 template <int dim, int spacedim=dim>
2978 class FEFaceValues : public FEFaceValuesBase<dim,spacedim>
2979 {
2980 public:
2985  static const unsigned int dimension = dim;
2986 
2987  static const unsigned int space_dimension = spacedim;
2988 
2993  static const unsigned int integral_dimension = dim-1;
2994 
3001  const Quadrature<dim-1> &quadrature,
3002  const UpdateFlags update_flags);
3003 
3010  const Quadrature<dim-1> &quadrature,
3011  const UpdateFlags update_flags);
3012 
3017  template <template <int, int> class DoFHandlerType, bool level_dof_access>
3018  void reinit (const TriaIterator<DoFCellAccessor<DoFHandlerType<dim,spacedim>,level_dof_access> > &cell,
3019  const unsigned int face_no);
3020 
3034  void reinit (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
3035  const unsigned int face_no);
3036 
3050  const FEFaceValues<dim,spacedim> &get_present_fe_values () const;
3051 private:
3052 
3056  void initialize (const UpdateFlags update_flags);
3057 
3064  void do_reinit (const unsigned int face_no);
3065 };
3066 
3067 
3085 template <int dim, int spacedim=dim>
3086 class FESubfaceValues : public FEFaceValuesBase<dim,spacedim>
3087 {
3088 public:
3092  static const unsigned int dimension = dim;
3093 
3097  static const unsigned int space_dimension = spacedim;
3098 
3103  static const unsigned int integral_dimension = dim-1;
3104 
3111  const Quadrature<dim-1> &face_quadrature,
3112  const UpdateFlags update_flags);
3113 
3120  const Quadrature<dim-1> &face_quadrature,
3121  const UpdateFlags update_flags);
3122 
3129  template <template <int, int> class DoFHandlerType, bool level_dof_access>
3130  void reinit (const TriaIterator<DoFCellAccessor<DoFHandlerType<dim,spacedim>,level_dof_access> > &cell,
3131  const unsigned int face_no,
3132  const unsigned int subface_no);
3133 
3147  void reinit (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
3148  const unsigned int face_no,
3149  const unsigned int subface_no);
3150 
3164  const FESubfaceValues<dim,spacedim> &get_present_fe_values () const;
3165 
3171  DeclException0 (ExcReinitCalledWithBoundaryFace);
3172 
3178  DeclException0 (ExcFaceHasNoSubfaces);
3179 
3180 private:
3181 
3185  void initialize (const UpdateFlags update_flags);
3186 
3193  void do_reinit (const unsigned int face_no,
3194  const unsigned int subface_no);
3195 };
3196 
3197 
3198 #ifndef DOXYGEN
3199 
3200 
3201 /*------------------------ Inline functions: namespace FEValuesViews --------*/
3202 
3203 namespace FEValuesViews
3204 {
3205  template <int dim, int spacedim>
3206  inline
3207  typename Scalar<dim,spacedim>::value_type
3208  Scalar<dim,spacedim>::value (const unsigned int shape_function,
3209  const unsigned int q_point) const
3210  {
3211  typedef FEValuesBase<dim,spacedim> FVB;
3212  Assert (shape_function < fe_values.fe->dofs_per_cell,
3213  ExcIndexRange (shape_function, 0, fe_values.fe->dofs_per_cell));
3214  Assert (fe_values.update_flags & update_values,
3215  typename FVB::ExcAccessToUninitializedField("update_values"));
3216 
3217  // an adaptation of the FEValuesBase::shape_value_component function
3218  // except that here we know the component as fixed and we have
3219  // pre-computed and cached a bunch of information. See the comments there.
3220  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
3221  return fe_values.finite_element_output.shape_values(shape_function_data[shape_function]
3222  .row_index,
3223  q_point);
3224  else
3225  return 0;
3226  }
3227 
3228 
3229 
3230 
3231  template <int dim, int spacedim>
3232  inline
3233  typename Scalar<dim,spacedim>::gradient_type
3234  Scalar<dim,spacedim>::gradient (const unsigned int shape_function,
3235  const unsigned int q_point) const
3236  {
3237  typedef FEValuesBase<dim,spacedim> FVB;
3238  Assert (shape_function < fe_values.fe->dofs_per_cell,
3239  ExcIndexRange (shape_function, 0, fe_values.fe->dofs_per_cell));
3240  Assert (fe_values.update_flags & update_gradients,
3241  typename FVB::ExcAccessToUninitializedField("update_gradients"));
3242 
3243  // an adaptation of the
3244  // FEValuesBase::shape_grad_component
3245  // function except that here we know the
3246  // component as fixed and we have
3247  // pre-computed and cached a bunch of
3248  // information. See the comments there.
3249  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
3250  return fe_values.finite_element_output.shape_gradients[shape_function_data[shape_function]
3251  .row_index][q_point];
3252  else
3253  return gradient_type();
3254  }
3255 
3256 
3257 
3258  template <int dim, int spacedim>
3259  inline
3260  typename Scalar<dim,spacedim>::hessian_type
3261  Scalar<dim,spacedim>::hessian (const unsigned int shape_function,
3262  const unsigned int q_point) const
3263  {
3264  typedef FEValuesBase<dim,spacedim> FVB;
3265  Assert (shape_function < fe_values.fe->dofs_per_cell,
3266  ExcIndexRange (shape_function, 0, fe_values.fe->dofs_per_cell));
3267  Assert (fe_values.update_flags & update_hessians,
3268  typename FVB::ExcAccessToUninitializedField("update_hessians"));
3269 
3270  // an adaptation of the
3271  // FEValuesBase::shape_hessian_component
3272  // function except that here we know the
3273  // component as fixed and we have
3274  // pre-computed and cached a bunch of
3275  // information. See the comments there.
3276  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
3277  return fe_values.finite_element_output.shape_hessians[shape_function_data[shape_function].row_index][q_point];
3278  else
3279  return hessian_type();
3280  }
3281 
3282 
3283 
3284  template <int dim, int spacedim>
3285  inline
3286  typename Scalar<dim,spacedim>::third_derivative_type
3287  Scalar<dim,spacedim>::third_derivative (const unsigned int shape_function,
3288  const unsigned int q_point) const
3289  {
3290  typedef FEValuesBase<dim,spacedim> FVB;
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 FVB::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  typedef FEValuesBase<dim,spacedim> FVB;
3317  Assert (shape_function < fe_values.fe->dofs_per_cell,
3318  ExcIndexRange (shape_function, 0, fe_values.fe->dofs_per_cell));
3319  Assert (fe_values.update_flags & update_values,
3320  typename FVB::ExcAccessToUninitializedField("update_values"));
3321 
3322  // same as for the scalar case except
3323  // that we have one more index
3324  const int snc = shape_function_data[shape_function].single_nonzero_component;
3325  if (snc == -2)
3326  return value_type();
3327  else if (snc != -1)
3328  {
3329  value_type return_value;
3330  return_value[shape_function_data[shape_function].single_nonzero_component_index]
3331  = fe_values.finite_element_output.shape_values(snc,q_point);
3332  return return_value;
3333  }
3334  else
3335  {
3336  value_type return_value;
3337  for (unsigned int d=0; d<dim; ++d)
3338  if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
3339  return_value[d]
3340  = fe_values.finite_element_output.shape_values(shape_function_data[shape_function].row_index[d],q_point);
3341 
3342  return return_value;
3343  }
3344  }
3345 
3346 
3347 
3348  template <int dim, int spacedim>
3349  inline
3351  Vector<dim,spacedim>::gradient (const unsigned int shape_function,
3352  const unsigned int q_point) const
3353  {
3354  typedef FEValuesBase<dim,spacedim> FVB;
3355  Assert (shape_function < fe_values.fe->dofs_per_cell,
3356  ExcIndexRange (shape_function, 0, fe_values.fe->dofs_per_cell));
3357  Assert (fe_values.update_flags & update_gradients,
3358  typename FVB::ExcAccessToUninitializedField("update_gradients"));
3359 
3360  // same as for the scalar case except
3361  // that we have one more index
3362  const int snc = shape_function_data[shape_function].single_nonzero_component;
3363  if (snc == -2)
3364  return gradient_type();
3365  else if (snc != -1)
3366  {
3367  gradient_type return_value;
3368  return_value[shape_function_data[shape_function].single_nonzero_component_index]
3369  = fe_values.finite_element_output.shape_gradients[snc][q_point];
3370  return return_value;
3371  }
3372  else
3373  {
3374  gradient_type return_value;
3375  for (unsigned int d=0; d<dim; ++d)
3376  if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
3377  return_value[d]
3378  = fe_values.finite_element_output.shape_gradients[shape_function_data[shape_function].row_index[d]][q_point];
3379 
3380  return return_value;
3381  }
3382  }
3383 
3384 
3385 
3386  template <int dim, int spacedim>
3387  inline
3389  Vector<dim,spacedim>::divergence (const unsigned int shape_function,
3390  const unsigned int q_point) const
3391  {
3392  // this function works like in
3393  // the case above
3394  typedef FEValuesBase<dim,spacedim> FVB;
3395  Assert (shape_function < fe_values.fe->dofs_per_cell,
3396  ExcIndexRange (shape_function, 0, fe_values.fe->dofs_per_cell));
3397  Assert (fe_values.update_flags & update_gradients,
3398  typename FVB::ExcAccessToUninitializedField("update_gradients"));
3399 
3400  // same as for the scalar case except
3401  // that we have one more index
3402  const int snc = shape_function_data[shape_function].single_nonzero_component;
3403  if (snc == -2)
3404  return divergence_type();
3405  else if (snc != -1)
3406  return
3407  fe_values.finite_element_output.shape_gradients[snc][q_point][shape_function_data[shape_function].single_nonzero_component_index];
3408  else
3409  {
3410  divergence_type return_value = 0;
3411  for (unsigned int d=0; d<dim; ++d)
3412  if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
3413  return_value
3414  += fe_values.finite_element_output.shape_gradients[shape_function_data[shape_function].row_index[d]][q_point][d];
3415 
3416  return return_value;
3417  }
3418  }
3419 
3420 
3421 
3422  template <int dim, int spacedim>
3423  inline
3425  Vector<dim,spacedim>::curl (const unsigned int shape_function, const unsigned int q_point) const
3426  {
3427  // this function works like in the case above
3428  typedef FEValuesBase<dim,spacedim> FVB;
3429 
3430  Assert (shape_function < fe_values.fe->dofs_per_cell,
3431  ExcIndexRange (shape_function, 0, fe_values.fe->dofs_per_cell));
3432  Assert (fe_values.update_flags & update_gradients,
3433  typename FVB::ExcAccessToUninitializedField("update_gradients"));
3434  // same as for the scalar case except that we have one more index
3435  const int snc = shape_function_data[shape_function].single_nonzero_component;
3436 
3437  if (snc == -2)
3438  return curl_type ();
3439 
3440  else
3441  switch (dim)
3442  {
3443  case 1:
3444  {
3445  Assert (false, ExcMessage("Computing the curl in 1d is not a useful operation"));
3446  return curl_type ();
3447  }
3448 
3449  case 2:
3450  {
3451  if (snc != -1)
3452  {
3453  curl_type return_value;
3454 
3455  // the single
3456  // nonzero component
3457  // can only be zero
3458  // or one in 2d
3459  if (shape_function_data[shape_function].single_nonzero_component_index == 0)
3460  return_value[0] = -1.0 * fe_values.finite_element_output.shape_gradients[snc][q_point][1];
3461  else
3462  return_value[0] = fe_values.finite_element_output.shape_gradients[snc][q_point][0];
3463 
3464  return return_value;
3465  }
3466 
3467  else
3468  {
3469  curl_type return_value;
3470 
3471  return_value[0] = 0.0;
3472 
3473  if (shape_function_data[shape_function].is_nonzero_shape_function_component[0])
3474  return_value[0]
3475  -= fe_values.finite_element_output.shape_gradients[shape_function_data[shape_function].row_index[0]][q_point][1];
3476 
3477  if (shape_function_data[shape_function].is_nonzero_shape_function_component[1])
3478  return_value[0]
3479  += fe_values.finite_element_output.shape_gradients[shape_function_data[shape_function].row_index[1]][q_point][0];
3480 
3481  return return_value;
3482  }
3483  }
3484 
3485  case 3:
3486  {
3487  if (snc != -1)
3488  {
3489  curl_type return_value;
3490 
3491  switch (shape_function_data[shape_function].single_nonzero_component_index)
3492  {
3493  case 0:
3494  {
3495  return_value[0] = 0;
3496  return_value[1] = fe_values.finite_element_output.shape_gradients[snc][q_point][2];
3497  return_value[2] = -1.0 * fe_values.finite_element_output.shape_gradients[snc][q_point][1];
3498  return return_value;
3499  }
3500 
3501  case 1:
3502  {
3503  return_value[0] = -1.0 * fe_values.finite_element_output.shape_gradients[snc][q_point][2];
3504  return_value[1] = 0;
3505  return_value[2] = fe_values.finite_element_output.shape_gradients[snc][q_point][0];
3506  return return_value;
3507  }
3508 
3509  default:
3510  {
3511  return_value[0] = fe_values.finite_element_output.shape_gradients[snc][q_point][1];
3512  return_value[1] = -1.0 * fe_values.finite_element_output.shape_gradients[snc][q_point][0];
3513  return_value[2] = 0;
3514  return return_value;
3515  }
3516  }
3517  }
3518 
3519  else
3520  {
3521  curl_type return_value;
3522 
3523  for (unsigned int i = 0; i < dim; ++i)
3524  return_value[i] = 0.0;
3525 
3526  if (shape_function_data[shape_function].is_nonzero_shape_function_component[0])
3527  {
3528  return_value[1]
3529  += fe_values.finite_element_output.shape_gradients[shape_function_data[shape_function].row_index[0]][q_point][2];
3530  return_value[2]
3531  -= fe_values.finite_element_output.shape_gradients[shape_function_data[shape_function].row_index[0]][q_point][1];
3532  }
3533 
3534  if (shape_function_data[shape_function].is_nonzero_shape_function_component[1])
3535  {
3536  return_value[0]
3537  -= fe_values.finite_element_output.shape_gradients[shape_function_data[shape_function].row_index[1]][q_point][2];
3538  return_value[2]
3539  += fe_values.finite_element_output.shape_gradients[shape_function_data[shape_function].row_index[1]][q_point][0];
3540  }
3541 
3542  if (shape_function_data[shape_function].is_nonzero_shape_function_component[2])
3543  {
3544  return_value[0]
3545  += fe_values.finite_element_output.shape_gradients[shape_function_data[shape_function].row_index[2]][q_point][1];
3546  return_value[1]
3547  -= fe_values.finite_element_output.shape_gradients[shape_function_data[shape_function].row_index[2]][q_point][0];
3548  }
3549 
3550  return return_value;
3551  }
3552  }
3553  }
3554  // should not end up here
3555  Assert (false, ExcInternalError());
3556  return curl_type();
3557  }
3558 
3559  template <int dim, int spacedim>
3560  inline
3562  Vector<dim,spacedim>::hessian (const unsigned int shape_function,
3563  const unsigned int q_point) const
3564  {
3565  // this function works like in
3566  // the case above
3567  typedef FEValuesBase<dim,spacedim> FVB;
3568  Assert (shape_function < fe_values.fe->dofs_per_cell,
3569  ExcIndexRange (shape_function, 0, fe_values.fe->dofs_per_cell));
3570  Assert (fe_values.update_flags & update_hessians,
3571  typename FVB::ExcAccessToUninitializedField("update_hessians"));
3572 
3573  // same as for the scalar case except
3574  // that we have one more index
3575  const int snc = shape_function_data[shape_function].single_nonzero_component;
3576  if (snc == -2)
3577  return hessian_type();
3578  else if (snc != -1)
3579  {
3580  hessian_type return_value;
3581  return_value[shape_function_data[shape_function].single_nonzero_component_index]
3582  = fe_values.finite_element_output.shape_hessians[snc][q_point];
3583  return return_value;
3584  }
3585  else
3586  {
3587  hessian_type return_value;
3588  for (unsigned int d=0; d<dim; ++d)
3589  if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
3590  return_value[d]
3591  = fe_values.finite_element_output.shape_hessians[shape_function_data[shape_function].row_index[d]][q_point];
3592 
3593  return return_value;
3594  }
3595  }
3596 
3597  template <int dim, int spacedim>
3598  inline
3600  Vector<dim,spacedim>::third_derivative (const unsigned int shape_function,
3601  const unsigned int q_point) const
3602  {
3603  // this function works like in
3604  // the case above
3605  typedef FEValuesBase<dim,spacedim> FVB;
3606  Assert (shape_function < fe_values.fe->dofs_per_cell,
3607  ExcIndexRange (shape_function, 0, fe_values.fe->dofs_per_cell));
3608  Assert (fe_values.update_flags & update_3rd_derivatives,
3609  typename FVB::ExcAccessToUninitializedField("update_3rd_derivatives"));
3610 
3611  // same as for the scalar case except
3612  // that we have one more index
3613  const int snc = shape_function_data[shape_function].single_nonzero_component;
3614  if (snc == -2)
3615  return third_derivative_type();
3616  else if (snc != -1)
3617  {
3618  third_derivative_type return_value;
3619  return_value[shape_function_data[shape_function].single_nonzero_component_index]
3620  = fe_values.finite_element_output.shape_3rd_derivatives[snc][q_point];
3621  return return_value;
3622  }
3623  else
3624  {
3625  third_derivative_type return_value;
3626  for (unsigned int d=0; d<dim; ++d)
3627  if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
3628  return_value[d]
3629  = fe_values.finite_element_output.shape_3rd_derivatives[shape_function_data[shape_function].row_index[d]][q_point];
3630 
3631  return return_value;
3632  }
3633  }
3634 
3635 
3636  namespace
3637  {
3642  inline
3643  ::SymmetricTensor<2,1>
3644  symmetrize_single_row (const unsigned int n,
3645  const ::Tensor<1,1> &t)
3646  {
3647  Assert (n < 1, ExcIndexRange (n, 0, 1));
3648  (void)n; // removes -Wunused-parameter warning in optimized mode
3649 
3650  const double array[1] = { t[0] };
3651  return ::SymmetricTensor<2,1>(array);
3652  }
3653 
3654 
3655  inline
3656  ::SymmetricTensor<2,2>
3657  symmetrize_single_row (const unsigned int n,
3658  const ::Tensor<1,2> &t)
3659  {
3660  switch (n)
3661  {
3662  case 0:
3663  {
3664  const double array[3] = { t[0], 0, t[1]/2 };
3665  return ::SymmetricTensor<2,2>(array);
3666  }
3667  case 1:
3668  {
3669  const double array[3] = { 0, t[1], t[0]/2 };
3670  return ::SymmetricTensor<2,2>(array);
3671  }
3672  default:
3673  {
3674  Assert (false, ExcIndexRange (n, 0, 2));
3675  return ::SymmetricTensor<2,2>();
3676  }
3677  }
3678  }
3679 
3680 
3681  inline
3682  ::SymmetricTensor<2,3>
3683  symmetrize_single_row (const unsigned int n,
3684  const ::Tensor<1,3> &t)
3685  {
3686  switch (n)
3687  {
3688  case 0:
3689  {
3690  const double array[6] = { t[0], 0, 0, t[1]/2, t[2]/2, 0 };
3691  return ::SymmetricTensor<2,3>(array);
3692  }
3693  case 1:
3694  {
3695  const double array[6] = { 0, t[1], 0, t[0]/2, 0, t[2]/2 };
3696  return ::SymmetricTensor<2,3>(array);
3697  }
3698  case 2:
3699  {
3700  const double array[6] = { 0, 0, t[2], 0, t[0]/2, t[1]/2 };
3701  return ::SymmetricTensor<2,3>(array);
3702  }
3703  default:
3704  {
3705  Assert (false, ExcIndexRange (n, 0, 3));
3706  return ::SymmetricTensor<2,3>();
3707  }
3708  }
3709  }
3710  }
3711 
3712 
3713  template <int dim, int spacedim>
3714  inline
3716  Vector<dim,spacedim>::symmetric_gradient (const unsigned int shape_function,
3717  const unsigned int q_point) const
3718  {
3719  typedef FEValuesBase<dim,spacedim> FVB;
3720  Assert (shape_function < fe_values.fe->dofs_per_cell,
3721  ExcIndexRange (shape_function, 0, fe_values.fe->dofs_per_cell));
3722  Assert (fe_values.update_flags & update_gradients,
3723  typename FVB::ExcAccessToUninitializedField("update_gradients"));
3724 
3725  // same as for the scalar case except
3726  // that we have one more index
3727  const int snc = shape_function_data[shape_function].single_nonzero_component;
3728  if (snc == -2)
3729  return symmetric_gradient_type();
3730  else if (snc != -1)
3731  return symmetrize_single_row (shape_function_data[shape_function].single_nonzero_component_index,
3732  fe_values.finite_element_output.shape_gradients[snc][q_point]);
3733  else
3734  {
3735  gradient_type return_value;
3736  for (unsigned int d=0; d<dim; ++d)
3737  if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
3738  return_value[d]
3739  = fe_values.finite_element_output.shape_gradients[shape_function_data[shape_function].row_index[d]][q_point];
3740 
3741  return symmetrize(return_value);
3742  }
3743  }
3744 
3745 
3746 
3747  template <int dim, int spacedim>
3748  inline
3750  SymmetricTensor<2, dim, spacedim>::value (const unsigned int shape_function,
3751  const unsigned int q_point) const
3752  {
3753  typedef FEValuesBase<dim,spacedim> FVB;
3754  Assert (shape_function < fe_values.fe->dofs_per_cell,
3755  ExcIndexRange (shape_function, 0, fe_values.fe->dofs_per_cell));
3756  Assert (fe_values.update_flags & update_values,
3757  typename FVB::ExcAccessToUninitializedField("update_values"));
3758 
3759  // similar to the vector case where we
3760  // have more then one index and we need
3761  // to convert between unrolled and
3762  // component indexing for tensors
3763  const int snc
3764  = shape_function_data[shape_function].single_nonzero_component;
3765 
3766  if (snc == -2)
3767  {
3768  // shape function is zero for the
3769  // selected components
3770  return value_type();
3771 
3772  }
3773  else if (snc != -1)
3774  {
3775  value_type return_value;
3776  const unsigned int comp =
3777  shape_function_data[shape_function].single_nonzero_component_index;
3778  return_value[value_type::unrolled_to_component_indices(comp)]
3779  = fe_values.finite_element_output.shape_values(snc,q_point);
3780  return return_value;
3781  }
3782  else
3783  {
3784  value_type return_value;
3785  for (unsigned int d = 0; d < value_type::n_independent_components; ++d)
3786  if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
3787  return_value[value_type::unrolled_to_component_indices(d)]
3788  = fe_values.finite_element_output.shape_values(shape_function_data[shape_function].row_index[d],q_point);
3789  return return_value;
3790  }
3791  }
3792 
3793 
3794  template <int dim, int spacedim>
3795  inline
3797  SymmetricTensor<2, dim, spacedim>::divergence(const unsigned int shape_function,
3798  const unsigned int q_point) const
3799  {
3800  typedef FEValuesBase<dim,spacedim> FVB;
3801  Assert (shape_function < fe_values.fe->dofs_per_cell,
3802  ExcIndexRange (shape_function, 0, fe_values.fe->dofs_per_cell));
3803  Assert (fe_values.update_flags & update_gradients,
3804  typename FVB::ExcAccessToUninitializedField("update_gradients"));
3805 
3806  const int snc = shape_function_data[shape_function].single_nonzero_component;
3807 
3808  if (snc == -2)
3809  {
3810  // shape function is zero for the
3811  // selected components
3812  return divergence_type();
3813  }
3814  else if (snc != -1)
3815  {
3816  // we have a single non-zero component
3817  // when the symmetric tensor is
3818  // represented in unrolled form.
3819  // this implies we potentially have
3820  // two non-zero components when
3821  // represented in component form! we
3822  // will only have one non-zero entry
3823  // if the non-zero component lies on
3824  // the diagonal of the tensor.
3825  //
3826  // the divergence of a second-order tensor
3827  // is a first order tensor.
3828  //
3829  // assume the second-order tensor is
3830  // A with components A_{ij}. then
3831  // A_{ij} = A_{ji} and there is only
3832  // one (if diagonal) or two non-zero
3833  // entries in the tensorial
3834  // representation. define the
3835  // divergence as:
3836  // b_i := \dfrac{\partial phi_{ij}}{\partial x_j}.
3837  // (which is incidentally also
3838  // b_j := \dfrac{\partial phi_{ij}}{\partial x_i}).
3839  // In both cases, a sum is implied.
3840  //
3841  // Now, we know the nonzero component
3842  // in unrolled form: it is indicated
3843  // by 'snc'. we can figure out which
3844  // tensor components belong to this:
3845  const unsigned int comp =
3846  shape_function_data[shape_function].single_nonzero_component_index;
3847  const unsigned int ii = value_type::unrolled_to_component_indices(comp)[0];
3848  const unsigned int jj = value_type::unrolled_to_component_indices(comp)[1];
3849 
3850  // given the form of the divergence
3851  // above, if ii=jj there is only a
3852  // single nonzero component of the
3853  // full tensor and the gradient
3854  // equals
3855  // b_ii := \dfrac{\partial phi_{ii,ii}}{\partial x_ii}.
3856  // all other entries of 'b' are zero
3857  //
3858  // on the other hand, if ii!=jj, then
3859  // there are two nonzero entries in
3860  // the full tensor and
3861  // b_ii := \dfrac{\partial phi_{ii,jj}}{\partial x_ii}.
3862  // b_jj := \dfrac{\partial phi_{ii,jj}}{\partial x_jj}.
3863  // again, all other entries of 'b' are
3864  // zero
3865  const ::Tensor<1, spacedim> phi_grad = fe_values.finite_element_output.shape_gradients[snc][q_point];
3866 
3867  divergence_type return_value;
3868  return_value[ii] = phi_grad[jj];
3869 
3870  if (ii != jj)
3871  return_value[jj] = phi_grad[ii];
3872 
3873  return return_value;
3874 
3875  }
3876  else
3877  {
3878  Assert (false, ExcNotImplemented());
3879  divergence_type return_value;
3880  return return_value;
3881  }
3882  }
3883 
3884  template <int dim, int spacedim>
3885  inline
3887  Tensor<2, dim, spacedim>::value (const unsigned int shape_function,
3888  const unsigned int q_point) const
3889  {
3890  typedef FEValuesBase<dim,spacedim> FVB;
3891  Assert (shape_function < fe_values.fe->dofs_per_cell,
3892  ExcIndexRange (shape_function, 0, fe_values.fe->dofs_per_cell));
3893  Assert (fe_values.update_flags & update_values,
3894  typename FVB::ExcAccessToUninitializedField("update_values"));
3895 
3896  // similar to the vector case where we
3897  // have more then one index and we need
3898  // to convert between unrolled and
3899  // component indexing for tensors
3900  const int snc
3901  = shape_function_data[shape_function].single_nonzero_component;
3902 
3903  if (snc == -2)
3904  {
3905  // shape function is zero for the
3906  // selected components
3907  return value_type();
3908 
3909  }
3910  else if (snc != -1)
3911  {
3912  value_type return_value;
3913  const unsigned int comp =
3914  shape_function_data[shape_function].single_nonzero_component_index;
3916  return_value[indices] = fe_values.finite_element_output.shape_values(snc,q_point);
3917  return return_value;
3918  }
3919  else
3920  {
3921  value_type return_value;
3922  for (unsigned int d = 0; d < dim*dim; ++d)
3923  if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
3924  {
3926  return_value[indices]
3927  = fe_values.finite_element_output.shape_values(shape_function_data[shape_function].row_index[d],q_point);
3928  }
3929  return return_value;
3930  }
3931  }
3932 
3933 
3934  template <int dim, int spacedim>
3935  inline
3937  Tensor<2, dim, spacedim>::divergence(const unsigned int shape_function,
3938  const unsigned int q_point) const
3939  {
3940  typedef FEValuesBase<dim,spacedim> FVB;
3941  Assert (shape_function < fe_values.fe->dofs_per_cell,
3942  ExcIndexRange (shape_function, 0, fe_values.fe->dofs_per_cell));
3943  Assert (fe_values.update_flags & update_gradients,
3944  typename FVB::ExcAccessToUninitializedField("update_gradients"));
3945 
3946  const int snc = shape_function_data[shape_function].single_nonzero_component;
3947 
3948  if (snc == -2)
3949  {
3950  // shape function is zero for the
3951  // selected components
3952  return divergence_type();
3953  }
3954  else if (snc != -1)
3955  {
3956  // we have a single non-zero component
3957  // when the tensor is
3958  // represented in unrolled form.
3959  //
3960  // the divergence of a second-order tensor
3961  // is a first order tensor.
3962  //
3963  // assume the second-order tensor is
3964  // A with components A_{ij}.
3965  // divergence as:
3966  // b_j := \dfrac{\partial phi_{ij}}{\partial x_i}.
3967  //
3968  // Now, we know the nonzero component
3969  // in unrolled form: it is indicated
3970  // by 'snc'. we can figure out which
3971  // tensor components belong to this:
3972  const unsigned int comp =
3973  shape_function_data[shape_function].single_nonzero_component_index;
3975  const unsigned int ii = indices[0];
3976  const unsigned int jj = indices[1];
3977 
3978  const ::Tensor<1, spacedim> phi_grad = fe_values.finite_element_output.shape_gradients[snc][q_point];
3979 
3980  divergence_type return_value;
3981  return_value[jj] = phi_grad[ii];
3982 
3983  return return_value;
3984 
3985  }
3986  else
3987  {
3988  Assert (false, ExcNotImplemented());
3989  divergence_type return_value;
3990  return return_value;
3991  }
3992  }
3993 }
3994 
3995 
3996 
3997 /*------------------------ Inline functions: FEValuesBase ------------------------*/
3998 
3999 
4000 
4001 template <int dim, int spacedim>
4002 inline
4005 operator[] (const FEValuesExtractors::Scalar &scalar) const
4006 {
4007  Assert (scalar.component < fe_values_views_cache.scalars.size(),
4008  ExcIndexRange (scalar.component,
4009  0, fe_values_views_cache.scalars.size()));
4010 
4011  return fe_values_views_cache.scalars[scalar.component];
4012 }
4013 
4014 
4015 
4016 template <int dim, int spacedim>
4017 inline
4020 operator[] (const FEValuesExtractors::Vector &vector) const
4021 {
4022  Assert (vector.first_vector_component <
4023  fe_values_views_cache.vectors.size(),
4024  ExcIndexRange (vector.first_vector_component,
4025  0, fe_values_views_cache.vectors.size()));
4026 
4027  return fe_values_views_cache.vectors[vector.first_vector_component];
4028 }
4029 
4030 template <int dim, int spacedim>
4031 inline
4035 {
4036  Assert (tensor.first_tensor_component <
4037  fe_values_views_cache.symmetric_second_order_tensors.size(),
4038  ExcIndexRange (tensor.first_tensor_component,
4039  0, fe_values_views_cache.symmetric_second_order_tensors.size()));
4040 
4041  return fe_values_views_cache.symmetric_second_order_tensors[tensor.first_tensor_component];
4042 }
4043 
4044 template <int dim, int spacedim>
4045 inline
4048 operator[] (const FEValuesExtractors::Tensor<2> &tensor) const
4049 {
4050  Assert (tensor.first_tensor_component <
4051  fe_values_views_cache.second_order_tensors.size(),
4052  ExcIndexRange (tensor.first_tensor_component,
4053  0, fe_values_views_cache.second_order_tensors.size()));
4054 
4055  return fe_values_views_cache.second_order_tensors[tensor.first_tensor_component];
4056 }
4057 
4058 
4059 
4060 
4061 template <int dim, int spacedim>
4062 inline
4063 const double &
4064 FEValuesBase<dim,spacedim>::shape_value (const unsigned int i,
4065  const unsigned int j) const
4066 {
4067  Assert (i < fe->dofs_per_cell,
4068  ExcIndexRange (i, 0, fe->dofs_per_cell));
4070  ExcAccessToUninitializedField("update_values"));
4071  Assert (fe->is_primitive (i),
4072  ExcShapeFunctionNotPrimitive(i));
4073 
4074  // if the entire FE is primitive,
4075  // then we can take a short-cut:
4076  if (fe->is_primitive())
4077  return this->finite_element_output.shape_values(i,j);
4078  else
4079  {
4080  // otherwise, use the mapping
4081  // between shape function
4082  // numbers and rows. note that
4083  // by the assertions above, we
4084  // know that this particular
4085  // shape function is primitive,
4086  // so we can call
4087  // system_to_component_index
4088  const unsigned int
4089  row = this->finite_element_output.shape_function_to_row_table[i * fe->n_components() + fe->system_to_component_index(i).first];
4090  return this->finite_element_output.shape_values(row, j);
4091  }
4092 }
4093 
4094 
4095 
4096 template <int dim, int spacedim>
4097 inline
4098 double
4100  const unsigned int j,
4101  const unsigned int component) const
4102 {
4103  Assert (i < fe->dofs_per_cell,
4104  ExcIndexRange (i, 0, fe->dofs_per_cell));
4105  Assert (this->update_flags & update_values,
4106  ExcAccessToUninitializedField("update_values"));
4107  Assert (component < fe->n_components(),
4108  ExcIndexRange(component, 0, fe->n_components()));
4109 
4110  // check whether the shape function
4111  // is non-zero at all within
4112  // this component:
4113  if (fe->get_nonzero_components(i)[component] == false)
4114  return 0;
4115 
4116  // look up the right row in the
4117  // table and take the data from
4118  // there
4119  const unsigned int
4120  row = this->finite_element_output.shape_function_to_row_table[i * fe->n_components() + component];
4121  return this->finite_element_output.shape_values(row, j);
4122 }
4123 
4124 
4125 
4126 template <int dim, int spacedim>
4127 inline
4128 const Tensor<1,spacedim> &
4129 FEValuesBase<dim,spacedim>::shape_grad (const unsigned int i,
4130  const unsigned int j) const
4131 {
4132  Assert (i < fe->dofs_per_cell,
4133  ExcIndexRange (i, 0, fe->dofs_per_cell));
4135  ExcAccessToUninitializedField("update_gradients"));
4136  Assert (fe->is_primitive (i),
4137  ExcShapeFunctionNotPrimitive(i));
4138 
4139  // if the entire FE is primitive,
4140  // then we can take a short-cut:
4141  if (fe->is_primitive())
4142  return this->finite_element_output.shape_gradients[i][j];
4143  else
4144  {
4145  // otherwise, use the mapping
4146  // between shape function
4147  // numbers and rows. note that
4148  // by the assertions above, we
4149  // know that this particular
4150  // shape function is primitive,
4151  // so we can call
4152  // system_to_component_index
4153  const unsigned int
4154  row = this->finite_element_output.shape_function_to_row_table[i * fe->n_components() + fe->system_to_component_index(i).first];
4155  return this->finite_element_output.shape_gradients[row][j];
4156  }
4157 }
4158 
4159 
4160 
4161 template <int dim, int spacedim>
4162 inline
4165  const unsigned int j,
4166  const unsigned int component) const
4167 {
4168  Assert (i < fe->dofs_per_cell,
4169  ExcIndexRange (i, 0, fe->dofs_per_cell));
4170  Assert (this->update_flags & update_gradients,
4171  ExcAccessToUninitializedField("update_gradients"));
4172  Assert (component < fe->n_components(),
4173  ExcIndexRange(component, 0, fe->n_components()));
4174 
4175  // check whether the shape function
4176  // is non-zero at all within
4177  // this component:
4178  if (fe->get_nonzero_components(i)[component] == false)
4179  return Tensor<1,spacedim>();
4180 
4181  // look up the right row in the
4182  // table and take the data from
4183  // there
4184  const unsigned int
4185  row = this->finite_element_output.shape_function_to_row_table[i * fe->n_components() + component];
4186  return this->finite_element_output.shape_gradients[row][j];
4187 }
4188 
4189 
4190 
4191 template <int dim, int spacedim>
4192 inline
4193 const Tensor<2,spacedim> &
4194 FEValuesBase<dim,spacedim>::shape_hessian (const unsigned int i,
4195  const unsigned int j) const
4196 {
4197  Assert (i < fe->dofs_per_cell,
4198  ExcIndexRange (i, 0, fe->dofs_per_cell));
4199  Assert (this->update_flags & update_hessians,
4200  ExcAccessToUninitializedField("update_hessians"));
4201  Assert (fe->is_primitive (i),
4202  ExcShapeFunctionNotPrimitive(i));
4203 
4204  // if the entire FE is primitive,
4205  // then we can take a short-cut:
4206  if (fe->is_primitive())
4207  return this->finite_element_output.shape_hessians[i][j];
4208  else
4209  {
4210  // otherwise, use the mapping
4211  // between shape function
4212  // numbers and rows. note that
4213  // by the assertions above, we
4214  // know that this particular
4215  // shape function is primitive,
4216  // so we can call
4217  // system_to_component_index
4218  const unsigned int
4219  row = this->finite_element_output.shape_function_to_row_table[i * fe->n_components() + fe->system_to_component_index(i).first];
4220  return this->finite_element_output.shape_hessians[row][j];
4221  }
4222 }
4223 
4224 
4225 
4226 template <int dim, int spacedim>
4227 inline
4230  const unsigned int j,
4231  const unsigned int component) const
4232 {
4233  Assert (i < fe->dofs_per_cell,
4234  ExcIndexRange (i, 0, fe->dofs_per_cell));
4235  Assert (this->update_flags & update_hessians,
4236  ExcAccessToUninitializedField("update_hessians"));
4237  Assert (component < fe->n_components(),
4238  ExcIndexRange(component, 0, fe->n_components()));
4239 
4240  // check whether the shape function
4241  // is non-zero at all within
4242  // this component:
4243  if (fe->get_nonzero_components(i)[component] == false)
4244  return Tensor<2,spacedim>();
4245 
4246  // look up the right row in the
4247  // table and take the data from
4248  // there
4249  const unsigned int
4250  row = this->finite_element_output.shape_function_to_row_table[i * fe->n_components() + component];
4251  return this->finite_element_output.shape_hessians[row][j];
4252 }
4253 
4254 
4255 
4256 template <int dim, int spacedim>
4257 inline
4258 const Tensor<3,spacedim> &
4260  const unsigned int j) const
4261 {
4262  Assert (i < fe->dofs_per_cell,
4263  ExcIndexRange (i, 0, fe->dofs_per_cell));
4264  Assert (this->update_flags & update_hessians,
4265  ExcAccessToUninitializedField("update_3rd_derivatives"));
4266  Assert (fe->is_primitive (i),
4267  ExcShapeFunctionNotPrimitive(i));
4268 
4269  // if the entire FE is primitive,
4270  // then we can take a short-cut:
4271  if (fe->is_primitive())
4272  return this->finite_element_output.shape_3rd_derivatives[i][j];
4273  else
4274  {
4275  // otherwise, use the mapping
4276  // between shape function
4277  // numbers and rows. note that
4278  // by the assertions above, we
4279  // know that this particular
4280  // shape function is primitive,
4281  // so we can call
4282  // system_to_component_index
4283  const unsigned int
4284  row = this->finite_element_output.shape_function_to_row_table[i * fe->n_components() + fe->system_to_component_index(i).first];
4285  return this->finite_element_output.shape_3rd_derivatives[row][j];
4286  }
4287 }
4288 
4289 
4290 
4291 template <int dim, int spacedim>
4292 inline
4295  const unsigned int j,
4296  const unsigned int component) const
4297 {
4298  Assert (i < fe->dofs_per_cell,
4299  ExcIndexRange (i, 0, fe->dofs_per_cell));
4300  Assert (this->update_flags & update_hessians,
4301  ExcAccessToUninitializedField("update_3rd_derivatives"));
4302  Assert (component < fe->n_components(),
4303  ExcIndexRange(component, 0, fe->n_components()));
4304 
4305  // check whether the shape function
4306  // is non-zero at all within
4307  // this component:
4308  if (fe->get_nonzero_components(i)[component] == false)
4309  return Tensor<3,spacedim>();
4310 
4311  // look up the right row in the
4312  // table and take the data from
4313  // there
4314  const unsigned int
4315  row = this->finite_element_output.shape_function_to_row_table[i * fe->n_components() + component];
4316  return this->finite_element_output.shape_3rd_derivatives[row][j];
4317 }
4318 
4319 
4320 
4321 template <int dim, int spacedim>
4322 inline
4325 {
4326  return *fe;
4327 }
4328 
4329 
4330 template <int dim, int spacedim>
4331 inline
4332 const Mapping<dim,spacedim> &
4334 {
4335  return *mapping;
4336 }
4337 
4338 
4339 
4340 template <int dim, int spacedim>
4341 inline
4344 {
4345  return this->update_flags;
4346 }
4347 
4348 
4349 
4350 template <int dim, int spacedim>
4351 inline
4352 const std::vector<Point<spacedim> > &
4354 {
4355  Assert (this->update_flags & update_quadrature_points,
4356  ExcAccessToUninitializedField("update_quadrature_points"));
4357  return this->mapping_output.quadrature_points;
4358 }
4359 
4360 
4361 
4362 template <int dim, int spacedim>
4363 inline
4364 const std::vector<double> &
4366 {
4367  Assert (this->update_flags & update_JxW_values,
4368  ExcAccessToUninitializedField("update_JxW_values"));
4369  return this->mapping_output.JxW_values;
4370 }
4371 
4372 
4373 
4374 template <int dim, int spacedim>
4375 inline
4376 const std::vector<DerivativeForm<1,dim,spacedim> > &
4378 {
4379  Assert (this->update_flags & update_jacobians,
4380  ExcAccessToUninitializedField("update_jacobians"));
4381  return this->mapping_output.jacobians;
4382 }
4383 
4384 
4385 
4386 template <int dim, int spacedim>
4387 inline
4388 const std::vector<DerivativeForm<2,dim,spacedim> > &
4390 {
4391  Assert (this->update_flags & update_jacobian_grads,
4392  ExcAccessToUninitializedField("update_jacobians_grads"));
4393  return this->mapping_output.jacobian_grads;
4394 }
4395 
4396 
4397 
4398 template <int dim, int spacedim>
4399 inline
4400 const Tensor<3,spacedim> &
4402 {
4403  Assert (this->update_flags & update_jacobian_pushed_forward_grads,
4404  ExcAccessToUninitializedField("update_jacobian_pushed_forward_grads"));
4405  return this->mapping_output.jacobian_pushed_forward_grads[i];
4406 }
4407 
4408 
4409 
4410 template <int dim, int spacedim>
4411 inline
4412 const std::vector<Tensor<3,spacedim> > &
4414 {
4415  Assert (this->update_flags & update_jacobian_pushed_forward_grads,
4416  ExcAccessToUninitializedField("update_jacobian_pushed_forward_grads"));
4417  return this->mapping_output.jacobian_pushed_forward_grads;
4418 }
4419 
4420 
4421 
4422 template <int dim, int spacedim>
4423 inline
4425 FEValuesBase<dim,spacedim>::jacobian_2nd_derivative (const unsigned int i) const
4426 {
4427  Assert (this->update_flags & update_jacobian_2nd_derivatives,
4428  ExcAccessToUninitializedField("update_jacobian_2nd_derivatives"));
4429  return this->mapping_output.jacobian_2nd_derivatives[i];
4430 }
4431 
4432 
4433 
4434 template <int dim, int spacedim>
4435 inline
4436 const std::vector<DerivativeForm<3,dim,spacedim> > &
4438 {
4439  Assert (this->update_flags & update_jacobian_2nd_derivatives,
4440  ExcAccessToUninitializedField("update_jacobian_2nd_derivatives"));
4441  return this->mapping_output.jacobian_2nd_derivatives;
4442 }
4443 
4444 
4445 
4446 template <int dim, int spacedim>
4447 inline
4448 const Tensor<4,spacedim> &
4450 {
4452  ExcAccessToUninitializedField("update_jacobian_pushed_forward_2nd_derivatives"));
4453  return this->mapping_output.jacobian_pushed_forward_2nd_derivatives[i];
4454 }
4455 
4456 
4457 
4458 template <int dim, int spacedim>
4459 inline
4460 const std::vector<Tensor<4,spacedim> > &
4462 {
4464  ExcAccessToUninitializedField("update_jacobian_pushed_forward_2nd_derivatives"));
4465  return this->mapping_output.jacobian_pushed_forward_2nd_derivatives;
4466 }
4467 
4468 
4469 
4470 template <int dim, int spacedim>
4471 inline
4473 FEValuesBase<dim,spacedim>::jacobian_3rd_derivative (const unsigned int i) const
4474 {
4475  Assert (this->update_flags & update_jacobian_3rd_derivatives,
4476  ExcAccessToUninitializedField("update_jacobian_3rd_derivatives"));
4477  return this->mapping_output.jacobian_3rd_derivatives[i];
4478 }
4479 
4480 
4481 
4482 template <int dim, int spacedim>
4483 inline
4484 const std::vector<DerivativeForm<4,dim,spacedim> > &
4486 {
4487  Assert (this->update_flags & update_jacobian_3rd_derivatives,
4488  ExcAccessToUninitializedField("update_jacobian_3rd_derivatives"));
4489  return this->mapping_output.jacobian_3rd_derivatives;
4490 }
4491 
4492 
4493 
4494 template <int dim, int spacedim>
4495 inline
4496 const Tensor<5,spacedim> &
4498 {
4500  ExcAccessToUninitializedField("update_jacobian_pushed_forward_3rd_derivatives"));
4501  return this->mapping_output.jacobian_pushed_forward_3rd_derivatives[i];
4502 }
4503 
4504 
4505 
4506 template <int dim, int spacedim>
4507 inline
4508 const std::vector<Tensor<5,spacedim> > &
4510 {
4512  ExcAccessToUninitializedField("update_jacobian_pushed_forward_3rd_derivatives"));
4513  return this->mapping_output.jacobian_pushed_forward_3rd_derivatives;
4514 }
4515 
4516 
4517 template <int dim, int spacedim>
4518 inline
4519 const std::vector<DerivativeForm<1,spacedim,dim> > &
4521 {
4522  Assert (this->update_flags & update_inverse_jacobians,
4523  ExcAccessToUninitializedField("update_inverse_jacobians"));
4524  return this->mapping_output.inverse_jacobians;
4525 }
4526 
4527 
4528 
4529 template <int dim, int spacedim>
4530 inline
4531 const Point<spacedim> &
4532 FEValuesBase<dim,spacedim>::quadrature_point (const unsigned int i) const
4533 {
4534  Assert (this->update_flags & update_quadrature_points,
4535  ExcAccessToUninitializedField("update_quadrature_points"));
4536  Assert (i<this->mapping_output.quadrature_points.size(),
4537  ExcIndexRange(i, 0, this->mapping_output.quadrature_points.size()));
4538 
4539  return this->mapping_output.quadrature_points[i];
4540 }
4541 
4542 
4543 
4544 
4545 template <int dim, int spacedim>
4546 inline
4547 double
4548 FEValuesBase<dim,spacedim>::JxW (const unsigned int i) const
4549 {
4550  Assert (this->update_flags & update_JxW_values,
4551  ExcAccessToUninitializedField("update_JxW_values"));
4552  Assert (i<this->mapping_output.JxW_values.size(),
4553  ExcIndexRange(i, 0, this->mapping_output.JxW_values.size()));
4554 
4555  return this->mapping_output.JxW_values[i];
4556 }
4557 
4558 
4559 
4560 template <int dim, int spacedim>
4561 inline
4563 FEValuesBase<dim,spacedim>::jacobian (const unsigned int i) const
4564 {
4565  Assert (this->update_flags & update_jacobians,
4566  ExcAccessToUninitializedField("update_jacobians"));
4567  Assert (i<this->mapping_output.jacobians.size(),
4568  ExcIndexRange(i, 0, this->mapping_output.jacobians.size()));
4569 
4570  return this->mapping_output.jacobians[i];
4571 }
4572 
4573 
4574 
4575 template <int dim, int spacedim>
4576 inline
4578 FEValuesBase<dim,spacedim>::jacobian_grad (const unsigned int i) const
4579 {
4580  Assert (this->update_flags & update_jacobian_grads,
4581  ExcAccessToUninitializedField("update_jacobians_grads"));
4582  Assert (i<this->mapping_output.jacobian_grads.size(),
4583  ExcIndexRange(i, 0, this->mapping_output.jacobian_grads.size()));
4584 
4585  return this->mapping_output.jacobian_grads[i];
4586 }
4587 
4588 
4589 
4590 template <int dim, int spacedim>
4591 inline
4593 FEValuesBase<dim,spacedim>::inverse_jacobian (const unsigned int i) const
4594 {
4595  Assert (this->update_flags & update_inverse_jacobians,
4596  ExcAccessToUninitializedField("update_inverse_jacobians"));
4597  Assert (i<this->mapping_output.inverse_jacobians.size(),
4598  ExcIndexRange(i, 0, this->mapping_output.inverse_jacobians.size()));
4599 
4600  return this->mapping_output.inverse_jacobians[i];
4601 }
4602 
4603 
4604 template <int dim, int spacedim>
4605 inline
4606 const Tensor<1,spacedim> &
4607 FEValuesBase<dim,spacedim>::normal_vector (const unsigned int i) const
4608 {
4609  typedef FEValuesBase<dim,spacedim> FVB;
4610  Assert (this->update_flags & update_normal_vectors,
4611  typename FVB::ExcAccessToUninitializedField("update_normal_vectors"));
4612  Assert (i<this->mapping_output.normal_vectors.size(),
4613  ExcIndexRange(i, 0, this->mapping_output.normal_vectors.size()));
4614 
4615  return this->mapping_output.normal_vectors[i];
4616 }
4617 
4618 
4619 
4620 /*------------------------ Inline functions: FEValues ----------------------------*/
4621 
4622 
4623 template <int dim, int spacedim>
4624 inline
4625 const Quadrature<dim> &
4627 {
4628  return quadrature;
4629 }
4630 
4631 
4632 
4633 template <int dim, int spacedim>
4634 inline
4635 const FEValues<dim,spacedim> &
4637 {
4638  return *this;
4639 }
4640 
4641 
4642 /*------------------------ Inline functions: FEFaceValuesBase --------------------*/
4643 
4644 
4645 template <int dim, int spacedim>
4646 inline
4647 unsigned int
4649 {
4650  return present_face_index;
4651 }
4652 
4653 
4654 /*------------------------ Inline functions: FE*FaceValues --------------------*/
4655 
4656 template <int dim, int spacedim>
4657 inline
4658 const Quadrature<dim-1> &
4660 {
4661  return quadrature;
4662 }
4663 
4664 
4665 
4666 template <int dim, int spacedim>
4667 inline
4670 {
4671  return *this;
4672 }
4673 
4674 
4675 
4676 template <int dim, int spacedim>
4677 inline
4680 {
4681  return *this;
4682 }
4683 
4684 
4685 
4686 template <int dim, int spacedim>
4687 inline
4688 const Tensor<1,spacedim> &
4689 FEFaceValuesBase<dim,spacedim>::boundary_form (const unsigned int i) const
4690 {
4691  typedef FEValuesBase<dim,spacedim> FVB;
4692  Assert (i<this->mapping_output.boundary_forms.size(),
4693  ExcIndexRange(i, 0, this->mapping_output.boundary_forms.size()));
4694  Assert (this->update_flags & update_boundary_forms,
4695  typename FVB::ExcAccessToUninitializedField("update_boundary_forms"));
4696 
4697  return this->mapping_output.boundary_forms[i];
4698 }
4699 
4700 #endif // DOXYGEN
4701 
4702 DEAL_II_NAMESPACE_CLOSE
4703 
4704 #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:2710
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:1485
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:2695
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:1340
const DerivativeForm< 4, dim, spacedim > & jacobian_3rd_derivative(const unsigned int quadrature_point) const
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:1268
hessian_type hessian(const unsigned int shape_function, const unsigned int q_point) const
Number value_type
Definition: vector.h:116
const FEValuesViews::Scalar< dim, spacedim > & operator[](const FEValuesExtractors::Scalar &scalar) const
CellSimilarity::Similarity cell_similarity
Definition: fe_values.h:2724
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:3267
::Tensor< 1, spacedim > divergence_type
Definition: fe_values.h:1131
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:1680
const FEValuesBase< dim, spacedim > & fe_values
Definition: fe_values.h:1285
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:2953
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:1316
std::vector< ShapeFunctionData > shape_function_data
Definition: fe_values.h:1091
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:1320
MappingType
Definition: mapping.h:50
const unsigned int dofs_per_cell
Definition: fe_values.h:1464
const DerivativeForm< 1, dim, spacedim > & jacobian(const unsigned int quadrature_point) const
const unsigned int component
Definition: fe_values.h:396
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:2276
::Tensor< 3, spacedim > hessian_type
Definition: fe_values.h:488
const Quadrature< dim-1 > quadrature
Definition: fe_values.h:2958
Volume element.
::ExceptionBase & ExcMessage(std::string arg1)
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:481
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
Scalar & operator=(const Scalar< dim, spacedim > &)
Definition: fe_values.cc:147
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:3385
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:916
unsigned int row_index[spacedim]
Definition: fe_values.h:522
void get_function_hessians(const InputVector &fe_function, std::vector< Tensor< 2, spacedim, typename InputVector::value_type > > &hessians) const
Definition: fe_values.cc:2994
void check_cell_similarity(const typename Triangulation< dim, spacedim >::cell_iterator &cell)
Definition: fe_values.cc:3531
static const unsigned int space_dimension
Definition: fe_values.h:1452
::internal::FEValuesViews::Cache< dim, spacedim > fe_values_views_cache
Definition: fe_values.h:2750
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:1364
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:1437
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:1413
const std::vector< double > & get_JxW_values() const
UpdateFlags get_update_flags() const
const Quadrature< dim > & get_quadrature() const
UpdateFlags compute_update_flags(const UpdateFlags update_flags) const
Definition: fe_values.cc:3458
BlockDynamicSparsityPattern BlockCompressedSparsityPattern DEAL_II_DEPRECATED
void get_function_laplacians(const InputVector &fe_function, std::vector< typename InputVector::value_type > &laplacians) const
Definition: fe_values.cc:3112
const SmartPointer< const FiniteElement< dim, spacedim >, FEValuesBase< dim, spacedim > > fe
Definition: fe_values.h:2688
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:1509
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:444
std_cxx11::unique_ptr< const CellIteratorBase > present_cell
Definition: fe_values.h:2628
third_derivative_type third_derivative(const unsigned int shape_function, const unsigned int q_point) const
boost::signals2::connection tria_listener
Definition: fe_values.h:2644
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:1296
#define Assert(cond, exc)
Definition: exceptions.h:294
::Tensor< 2, spacedim > value_type
Definition: fe_values.h:1126
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:1559
gradient_type gradient(const unsigned int shape_function, const unsigned int q_point) const
Abstract base class for mapping classes.
Definition: dof_tools.h:52
std::vector< ShapeFunctionData > shape_function_data
Definition: fe_values.h:875
const Quadrature< dim > quadrature
Definition: fe_values.h:2857
const unsigned int first_vector_component
Definition: fe_values.h:870
std::size_t memory_consumption() const
Definition: fe_values.cc:3438
::internal::FEValues::FiniteElementRelatedData< dim, spacedim > finite_element_output
Definition: fe_values.h:2701
void invalidate_present_cell()
Definition: fe_values.cc:3476
const SmartPointer< const Mapping< dim, spacedim >, FEValuesBase< dim, spacedim > > mapping
Definition: fe_values.h:2668
const std::vector< Tensor< 3, spacedim > > & get_jacobian_pushed_forward_grads() const
Vector & operator=(const Vector< dim, spacedim > &)
Definition: fe_values.cc:242
Tensor()
Definition: tensor.h:781
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:511
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
static TableIndices< rank_ > unrolled_to_component_indices(const unsigned int i)
Definition: tensor.h:1107
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:1388
void transform(std::vector< Tensor< 1, spacedim > > &transformed, const std::vector< Tensor< 1, dim > > &original, MappingType mapping) const DEAL_II_DEPRECATED
Definition: fe_values.cc:3425
DeclException0(ExcCannotInitializeField)
const DerivativeForm< 1, spacedim, dim > & inverse_jacobian(const unsigned int quadrature_point) const
const unsigned int n_quadrature_points
Definition: fe_values.h:1457
const double & shape_value(const unsigned int function_no, const unsigned int point_no) const
static const unsigned int dimension
Definition: fe_values.h:1447
::Tensor< 2, spacedim > hessian_type
Definition: fe_values.h:164
const FEValuesBase< dim, spacedim > & fe_values
Definition: fe_values.h:390
Tensor< 3, spacedim > shape_3rd_derivative_component(const unsigned int function_no, const unsigned int point_no, const unsigned int component) const
Definition: mpi.h:48
::SymmetricTensor< 2, spacedim > symmetric_gradient_type
Definition: fe_values.h:466
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:3493
::Tensor< 2, spacedim > gradient_type
Definition: fe_values.h:454
std_cxx11::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > mapping_data
Definition: fe_values.h:2675
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:1462
Definition: fe.h:31
const FEValuesBase< dim, spacedim > & fe_values
Definition: fe_values.h:864
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:3394
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:1533
const FiniteElement< dim, spacedim > & get_fe() const
::internal::FEValues::MappingRelatedData< dim, spacedim > mapping_output
Definition: fe_values.h:2681
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
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:171
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
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
Definition: fe_values.cc:3585
std::vector< ShapeFunctionData > shape_function_data
Definition: fe_values.h:401
UpdateFlags update_flags
Definition: fe_values.h:2707
const FEValuesBase< dim, spacedim > & fe_values
Definition: fe_values.h:1080
void get_function_gradients(const InputVector &fe_function, std::vector< Tensor< 1, spacedim, typename InputVector::value_type > > &gradients) const
Definition: fe_values.cc:2871
::Tensor< 4, spacedim > third_derivative_type
Definition: fe_values.h:495
std::vector< Point< spacedim > > get_normal_vectors() const DEAL_II_DEPRECATED
Definition: fe_values.cc:3406