Reference documentation for deal.II version Git 04ded8f 2017-09-19 10:11:38 +0200
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>
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/grid/tria.h>
30 #include <deal.II/grid/tria_iterator.h>
31 #include <deal.II/dofs/dof_handler.h>
32 #include <deal.II/dofs/dof_accessor.h>
33 #include <deal.II/hp/dof_handler.h>
34 #include <deal.II/fe/fe.h>
35 #include <deal.II/fe/fe_update_flags.h>
36 #include <deal.II/fe/fe_values_extractors.h>
37 #include <deal.II/fe/mapping.h>
38 
39 #include <algorithm>
40 #include <memory>
41 #include <type_traits>
42 
43 
44 // dummy include in order to have the
45 // definition of PetscScalar available
46 // without including other PETSc stuff
47 #ifdef DEAL_II_WITH_PETSC
48 # include <petsc.h>
49 #endif
50 
51 DEAL_II_NAMESPACE_OPEN
52 
53 template <int dim, int spacedim=dim> class FEValuesBase;
54 
55 namespace internal
56 {
61  template <int dim, class NumberType=double>
62  struct CurlType;
63 
70  template <class NumberType>
71  struct CurlType<1, NumberType>
72  {
74  };
75 
82  template <class NumberType>
83  struct CurlType<2,NumberType>
84  {
86  };
87 
94  template <class NumberType>
95  struct CurlType<3,NumberType>
96  {
98  };
99 }
100 
101 
102 
103 
125 namespace FEValuesViews
126 {
138  template <int dim, int spacedim=dim>
139  class Scalar
140  {
141  public:
147  typedef double value_type;
148 
154  typedef ::Tensor<1,spacedim> gradient_type;
155 
161  typedef ::Tensor<2,spacedim> hessian_type;
162 
168  typedef ::Tensor<3,spacedim> third_derivative_type;
169 
174  template <typename Number>
175  struct OutputType
176  {
182 
188 
194 
200 
206  };
207 
213  {
223 
232  unsigned int row_index;
233  };
234 
238  Scalar ();
239 
245  Scalar (const FEValuesBase<dim,spacedim> &fe_values_base,
246  const unsigned int component);
247 
253 
267  value_type
268  value (const unsigned int shape_function,
269  const unsigned int q_point) const;
270 
281  gradient_type
282  gradient (const unsigned int shape_function,
283  const unsigned int q_point) const;
284 
295  hessian_type
296  hessian (const unsigned int shape_function,
297  const unsigned int q_point) const;
298 
309  third_derivative_type
310  third_derivative (const unsigned int shape_function,
311  const unsigned int q_point) const;
312 
330  template <class InputVector>
331  void get_function_values (const InputVector &fe_function,
332  std::vector<typename ProductType<value_type,typename InputVector::value_type>::type> &values) const;
333 
354  template <class InputVector>
355  void get_function_values_from_local_dof_values (const InputVector &dof_values,
356  std::vector<typename OutputType<typename InputVector::value_type>::value_type> &values) const;
357 
375  template <class InputVector>
376  void get_function_gradients (const InputVector &fe_function,
377  std::vector<typename ProductType<gradient_type,typename InputVector::value_type>::type> &gradients) const;
378 
382  template <class InputVector>
383  void get_function_gradients_from_local_dof_values (const InputVector &dof_values,
384  std::vector<typename OutputType<typename InputVector::value_type>::gradient_type> &gradients) const;
385 
403  template <class InputVector>
404  void get_function_hessians (const InputVector &fe_function,
405  std::vector<typename ProductType<hessian_type,typename InputVector::value_type>::type> &hessians) const;
406 
410  template <class InputVector>
411  void get_function_hessians_from_local_dof_values (const InputVector &dof_values,
412  std::vector<typename OutputType<typename InputVector::value_type>::hessian_type> &hessians) const;
413 
414 
433  template <class InputVector>
434  void get_function_laplacians (const InputVector &fe_function,
435  std::vector<typename ProductType<value_type,typename InputVector::value_type>::type> &laplacians) const;
436 
440  template <class InputVector>
441  void get_function_laplacians_from_local_dof_values (const InputVector &dof_values,
442  std::vector<typename OutputType<typename InputVector::value_type>::laplacian_type> &laplacians) const;
443 
444 
463  template <class InputVector>
464  void get_function_third_derivatives (const InputVector &fe_function,
465  std::vector<typename ProductType<third_derivative_type,
466  typename InputVector::value_type>::type> &third_derivatives) const;
467 
471  template <class InputVector>
472  void get_function_third_derivatives_from_local_dof_values (const InputVector &dof_values,
473  std::vector<typename OutputType<typename InputVector::value_type>::third_derivative_type> &third_derivatives) const;
474 
475 
476  private:
481 
486  const unsigned int component;
487 
491  std::vector<ShapeFunctionData> shape_function_data;
492  };
493 
494 
495 
525  template <int dim, int spacedim=dim>
526  class Vector
527  {
528  public:
534  typedef ::Tensor<1,spacedim> value_type;
535 
544  typedef ::Tensor<2,spacedim> gradient_type;
545 
556  typedef ::SymmetricTensor<2,spacedim> symmetric_gradient_type;
557 
563  typedef double divergence_type;
564 
571  typedef typename ::internal::CurlType<spacedim>::type curl_type;
572 
578  typedef ::Tensor<3,spacedim> hessian_type;
579 
585  typedef ::Tensor<4,spacedim> third_derivative_type;
586 
591  template <typename Number>
592  struct OutputType
593  {
599 
605 
611 
617 
623 
629 
635 
641  };
642 
648  {
658 
668  unsigned int row_index[spacedim];
669 
679  unsigned int single_nonzero_component_index;
680  };
681 
685  Vector ();
686 
695  Vector (const FEValuesBase<dim,spacedim> &fe_values_base,
696  const unsigned int first_vector_component);
697 
703 
720  value_type
721  value (const unsigned int shape_function,
722  const unsigned int q_point) const;
723 
737  gradient_type
738  gradient (const unsigned int shape_function,
739  const unsigned int q_point) const;
740 
756  symmetric_gradient_type
757  symmetric_gradient (const unsigned int shape_function,
758  const unsigned int q_point) const;
759 
770  divergence_type
771  divergence (const unsigned int shape_function,
772  const unsigned int q_point) const;
773 
794  curl_type
795  curl (const unsigned int shape_function,
796  const unsigned int q_point) const;
797 
808  hessian_type
809  hessian (const unsigned int shape_function,
810  const unsigned int q_point) const;
811 
822  third_derivative_type
823  third_derivative (const unsigned int shape_function,
824  const unsigned int q_point) const;
825 
843  template <class InputVector>
844  void get_function_values (const InputVector &fe_function,
845  std::vector<typename ProductType<value_type,typename InputVector::value_type>::type> &values) const;
846 
867  template <class InputVector>
868  void get_function_values_from_local_dof_values (const InputVector &dof_values,
869  std::vector<typename OutputType<typename InputVector::value_type>::value_type> &values) const;
870 
888  template <class InputVector>
889  void get_function_gradients (const InputVector &fe_function,
890  std::vector<typename ProductType<gradient_type,typename InputVector::value_type>::type> &gradients) const;
891 
895  template <class InputVector>
896  void get_function_gradients_from_local_dof_values (const InputVector &dof_values,
897  std::vector<typename OutputType<typename InputVector::value_type>::gradient_type> &gradients) const;
898 
922  template <class InputVector>
923  void
924  get_function_symmetric_gradients (const InputVector &fe_function,
925  std::vector<typename ProductType<symmetric_gradient_type,typename InputVector::value_type>::type> &symmetric_gradients) const;
926 
930  template <class InputVector>
931  void
932  get_function_symmetric_gradients_from_local_dof_values (const InputVector &dof_values,
933  std::vector<typename OutputType<typename InputVector::value_type>::symmetric_gradient_type> &symmetric_gradients) const;
934 
953  template <class InputVector>
954  void get_function_divergences (const InputVector &fe_function,
955  std::vector<typename ProductType<divergence_type,typename InputVector::value_type>::type> &divergences) const;
956 
960  template <class InputVector>
961  void get_function_divergences_from_local_dof_values (const InputVector &dof_values,
962  std::vector<typename OutputType<typename InputVector::value_type>::divergence_type> &divergences) const;
963 
982  template <class InputVector>
983  void get_function_curls (const InputVector &fe_function,
984  std::vector<typename ProductType<curl_type,typename InputVector::value_type>::type> &curls) const;
985 
989  template <class InputVector>
990  void get_function_curls_from_local_dof_values (const InputVector &dof_values,
991  std::vector<typename OutputType<typename InputVector::value_type>::curl_type> &curls) const;
992 
1010  template <class InputVector>
1011  void get_function_hessians (const InputVector &fe_function,
1012  std::vector<typename ProductType<hessian_type,typename InputVector::value_type>::type> &hessians) const;
1013 
1017  template <class InputVector>
1018  void get_function_hessians_from_local_dof_values (const InputVector &dof_values,
1019  std::vector<typename OutputType<typename InputVector::value_type>::hessian_type> &hessians) const;
1020 
1039  template <class InputVector>
1040  void get_function_laplacians (const InputVector &fe_function,
1041  std::vector<typename ProductType<value_type,typename InputVector::value_type>::type> &laplacians) const;
1042 
1046  template <class InputVector>
1047  void get_function_laplacians_from_local_dof_values (const InputVector &dof_values,
1048  std::vector<typename OutputType<typename InputVector::value_type>::laplacian_type> &laplacians) const;
1049 
1068  template <class InputVector>
1069  void get_function_third_derivatives (const InputVector &fe_function,
1070  std::vector<typename ProductType<third_derivative_type,
1071  typename InputVector::value_type>::type> &third_derivatives) const;
1072 
1076  template <class InputVector>
1077  void get_function_third_derivatives_from_local_dof_values (const InputVector &dof_values,
1078  std::vector<typename OutputType<typename InputVector::value_type>::third_derivative_type> &third_derivatives) const;
1079 
1080  private:
1085 
1090  const unsigned int first_vector_component;
1091 
1095  std::vector<ShapeFunctionData> shape_function_data;
1096  };
1097 
1098 
1099  template <int rank, int dim, int spacedim = dim>
1100  class SymmetricTensor;
1101 
1126  template <int dim, int spacedim>
1127  class SymmetricTensor<2,dim,spacedim>
1128  {
1129  public:
1136  typedef ::SymmetricTensor<2, spacedim> value_type;
1137 
1147  typedef ::Tensor<1, spacedim> divergence_type;
1148 
1153  template <typename Number>
1154  struct OutputType
1155  {
1161 
1167  };
1168 
1173  struct ShapeFunctionData
1174  {
1183  bool is_nonzero_shape_function_component[value_type::n_independent_components];
1184 
1194  unsigned int row_index[value_type::n_independent_components];
1195 
1205  unsigned int single_nonzero_component_index;
1206  };
1207 
1211  SymmetricTensor();
1212 
1222  SymmetricTensor(const FEValuesBase<dim, spacedim> &fe_values_base,
1223  const unsigned int first_tensor_component);
1224 
1229  SymmetricTensor &operator=(const SymmetricTensor<2, dim, spacedim> &);
1230 
1248  value_type
1249  value (const unsigned int shape_function,
1250  const unsigned int q_point) const;
1251 
1252 
1266  divergence_type
1267  divergence (const unsigned int shape_function,
1268  const unsigned int q_point) const;
1269 
1287  template <class InputVector>
1288  void get_function_values (const InputVector &fe_function,
1289  std::vector<typename ProductType<value_type,typename InputVector::value_type>::type> &values) const;
1290 
1311  template <class InputVector>
1312  void get_function_values_from_local_dof_values (const InputVector &dof_values,
1313  std::vector<typename OutputType<typename InputVector::value_type>::value_type> &values) const;
1314 
1315 
1337  template <class InputVector>
1338  void get_function_divergences (const InputVector &fe_function,
1339  std::vector<typename ProductType<divergence_type,typename InputVector::value_type>::type> &divergences) const;
1340 
1344  template <class InputVector>
1345  void get_function_divergences_from_local_dof_values (const InputVector &dof_values,
1346  std::vector<typename OutputType<typename InputVector::value_type>::divergence_type> &divergences) const;
1347 
1348  private:
1353 
1358  const unsigned int first_tensor_component;
1359 
1363  std::vector<ShapeFunctionData> shape_function_data;
1364  };
1365 
1366 
1367  template <int rank, int dim, int spacedim = dim>
1368  class Tensor;
1369 
1389  template <int dim, int spacedim>
1390  class Tensor<2,dim,spacedim>
1391  {
1392  public:
1393 
1398  typedef ::Tensor<2, spacedim> value_type;
1399 
1403  typedef ::Tensor<1, spacedim> divergence_type;
1404 
1409  template <typename Number>
1410  struct OutputType
1411  {
1417 
1423  };
1424 
1429  struct ShapeFunctionData
1430  {
1439  bool is_nonzero_shape_function_component[value_type::n_independent_components];
1440 
1450  unsigned int row_index[value_type::n_independent_components];
1451 
1461  unsigned int single_nonzero_component_index;
1462  };
1463 
1467  Tensor();
1468 
1469 
1479  Tensor(const FEValuesBase<dim, spacedim> &fe_values_base,
1480  const unsigned int first_tensor_component);
1481 
1482 
1487  Tensor &operator=(const Tensor<2, dim, spacedim> &);
1488 
1505  value_type
1506  value (const unsigned int shape_function,
1507  const unsigned int q_point) const;
1508 
1522  divergence_type
1523  divergence (const unsigned int shape_function,
1524  const unsigned int q_point) const;
1525 
1543  template <class InputVector>
1544  void get_function_values (const InputVector &fe_function,
1545  std::vector<typename ProductType<value_type,typename InputVector::value_type>::type> &values) const;
1546 
1567  template <class InputVector>
1568  void get_function_values_from_local_dof_values (const InputVector &dof_values,
1569  std::vector<typename OutputType<typename InputVector::value_type>::value_type> &values) const;
1570 
1571 
1593  template <class InputVector>
1594  void get_function_divergences (const InputVector &fe_function,
1595  std::vector<typename ProductType<divergence_type,typename InputVector::value_type>::type> &divergences) const;
1596 
1600  template <class InputVector>
1601  void get_function_divergences_from_local_dof_values (const InputVector &dof_values,
1602  std::vector<typename OutputType<typename InputVector::value_type>::divergence_type> &values) const;
1603 
1604 
1605  private:
1610 
1615  const unsigned int first_tensor_component;
1616 
1620  std::vector<ShapeFunctionData> shape_function_data;
1621  };
1622 
1623 }
1624 
1625 
1626 namespace internal
1627 {
1628  namespace FEValuesViews
1629  {
1637  template <int dim, int spacedim>
1638  struct Cache
1639  {
1644  std::vector<::FEValuesViews::Scalar<dim,spacedim> > scalars;
1645  std::vector<::FEValuesViews::Vector<dim,spacedim> > vectors;
1646  std::vector<::FEValuesViews::SymmetricTensor<2,dim,spacedim> >
1647  symmetric_second_order_tensors;
1648  std::vector<::FEValuesViews::Tensor<2,dim,spacedim> >
1649  second_order_tensors;
1650 
1654  Cache (const FEValuesBase<dim,spacedim> &fe_values);
1655  };
1656  }
1657 }
1658 
1659 
1660 
1763 template <int dim, int spacedim>
1764 class FEValuesBase :
1765  public Subscriptor
1766 {
1767 public:
1771  static const unsigned int dimension = dim;
1772 
1776  static const unsigned int space_dimension = spacedim;
1777 
1781  const unsigned int n_quadrature_points;
1782 
1788  const unsigned int dofs_per_cell;
1789 
1790 
1798  FEValuesBase (const unsigned int n_q_points,
1799  const unsigned int dofs_per_cell,
1800  const UpdateFlags update_flags,
1803 
1804 
1808  ~FEValuesBase ();
1809 
1810 
1812 
1813 
1834  const double &shape_value (const unsigned int function_no,
1835  const unsigned int point_no) const;
1836 
1857  double shape_value_component (const unsigned int function_no,
1858  const unsigned int point_no,
1859  const unsigned int component) const;
1860 
1886  const Tensor<1,spacedim> &
1887  shape_grad (const unsigned int function_no,
1888  const unsigned int quadrature_point) const;
1889 
1907  shape_grad_component (const unsigned int function_no,
1908  const unsigned int point_no,
1909  const unsigned int component) const;
1910 
1930  const Tensor<2,spacedim> &
1931  shape_hessian (const unsigned int function_no,
1932  const unsigned int point_no) const;
1933 
1951  shape_hessian_component (const unsigned int function_no,
1952  const unsigned int point_no,
1953  const unsigned int component) const;
1954 
1974  const Tensor<3,spacedim> &
1975  shape_3rd_derivative (const unsigned int function_no,
1976  const unsigned int point_no) const;
1977 
1995  shape_3rd_derivative_component (const unsigned int function_no,
1996  const unsigned int point_no,
1997  const unsigned int component) const;
1998 
2000 
2002 
2041  template <class InputVector>
2042  void get_function_values (const InputVector &fe_function,
2043  std::vector<typename InputVector::value_type> &values) const;
2044 
2058  template <class InputVector>
2059  void get_function_values (const InputVector &fe_function,
2060  std::vector<Vector<typename InputVector::value_type> > &values) const;
2061 
2080  template <class InputVector>
2081  void get_function_values (const InputVector &fe_function,
2082  const VectorSlice<const std::vector<types::global_dof_index> > &indices,
2083  std::vector<typename InputVector::value_type> &values) const;
2084 
2106  template <class InputVector>
2107  void get_function_values (const InputVector &fe_function,
2108  const VectorSlice<const std::vector<types::global_dof_index> > &indices,
2109  std::vector<Vector<typename InputVector::value_type> > &values) const;
2110 
2111 
2142  template <class InputVector>
2143  void get_function_values (const InputVector &fe_function,
2144  const VectorSlice<const std::vector<types::global_dof_index> > &indices,
2145  VectorSlice<std::vector<std::vector<typename InputVector::value_type> > > values,
2146  const bool quadrature_points_fastest) const;
2147 
2149 
2151 
2190  template <class InputVector>
2191  void get_function_gradients (const InputVector &fe_function,
2192  std::vector<Tensor<1,spacedim,typename InputVector::value_type> > &gradients) const;
2193 
2210  template <class InputVector>
2211  void get_function_gradients (const InputVector &fe_function,
2212  std::vector<std::vector<Tensor<1,spacedim,typename InputVector::value_type> > > &gradients) const;
2213 
2220  template <class InputVector>
2221  void get_function_gradients (const InputVector &fe_function,
2222  const VectorSlice<const std::vector<types::global_dof_index> > &indices,
2223  std::vector<Tensor<1,spacedim,typename InputVector::value_type> > &gradients) const;
2224 
2231  template <class InputVector>
2232  void get_function_gradients (const InputVector &fe_function,
2233  const VectorSlice<const std::vector<types::global_dof_index> > &indices,
2234  VectorSlice<std::vector<std::vector<Tensor<1,spacedim,typename InputVector::value_type> > > > gradients,
2235  bool quadrature_points_fastest = false) const;
2236 
2238 
2240 
2280  template <class InputVector>
2281  void
2282  get_function_hessians (const InputVector &fe_function,
2283  std::vector<Tensor<2,spacedim,typename InputVector::value_type> > &hessians) const;
2284 
2302  template <class InputVector>
2303  void
2304  get_function_hessians (const InputVector &fe_function,
2305  std::vector<std::vector<Tensor<2,spacedim,typename InputVector::value_type> > > &hessians,
2306  bool quadrature_points_fastest = false) const;
2307 
2312  template <class InputVector>
2313  void get_function_hessians (
2314  const InputVector &fe_function,
2315  const VectorSlice<const std::vector<types::global_dof_index> > &indices,
2316  std::vector<Tensor<2,spacedim,typename InputVector::value_type> > &hessians) const;
2317 
2324  template <class InputVector>
2325  void get_function_hessians (
2326  const InputVector &fe_function,
2327  const VectorSlice<const std::vector<types::global_dof_index> > &indices,
2328  VectorSlice<std::vector<std::vector<Tensor<2,spacedim,typename InputVector::value_type> > > > hessians,
2329  bool quadrature_points_fastest = false) const;
2330 
2373  template <class InputVector>
2374  void
2375  get_function_laplacians (const InputVector &fe_function,
2376  std::vector<typename InputVector::value_type> &laplacians) const;
2377 
2397  template <class InputVector>
2398  void
2399  get_function_laplacians (const InputVector &fe_function,
2400  std::vector<Vector<typename InputVector::value_type> > &laplacians) const;
2401 
2408  template <class InputVector>
2410  const InputVector &fe_function,
2411  const VectorSlice<const std::vector<types::global_dof_index> > &indices,
2412  std::vector<typename InputVector::value_type> &laplacians) const;
2413 
2420  template <class InputVector>
2422  const InputVector &fe_function,
2423  const VectorSlice<const std::vector<types::global_dof_index> > &indices,
2424  std::vector<Vector<typename InputVector::value_type> > &laplacians) const;
2425 
2432  template <class InputVector>
2434  const InputVector &fe_function,
2435  const VectorSlice<const std::vector<types::global_dof_index> > &indices,
2436  std::vector<std::vector<typename InputVector::value_type> > &laplacians,
2437  bool quadrature_points_fastest = false) const;
2438 
2440 
2442 
2483  template <class InputVector>
2484  void
2485  get_function_third_derivatives (const InputVector &fe_function,
2486  std::vector<Tensor<3,spacedim,typename InputVector::value_type> > &third_derivatives) const;
2487 
2506  template <class InputVector>
2507  void
2508  get_function_third_derivatives (const InputVector &fe_function,
2509  std::vector<std::vector<Tensor<3,spacedim,typename InputVector::value_type> > > &third_derivatives,
2510  bool quadrature_points_fastest = false) const;
2511 
2516  template <class InputVector>
2518  const InputVector &fe_function,
2519  const VectorSlice<const std::vector<types::global_dof_index> > &indices,
2520  std::vector<Tensor<3,spacedim,typename InputVector::value_type> > &third_derivatives) const;
2521 
2528  template <class InputVector>
2530  const InputVector &fe_function,
2531  const VectorSlice<const std::vector<types::global_dof_index> > &indices,
2532  VectorSlice<std::vector<std::vector<Tensor<3,spacedim,typename InputVector::value_type> > > > third_derivatives,
2533  bool quadrature_points_fastest = false) const;
2535 
2537 
2538 
2544  const Point<spacedim> &
2545  quadrature_point (const unsigned int q) const;
2546 
2552  const std::vector<Point<spacedim> > &get_quadrature_points () const;
2553 
2569  double JxW (const unsigned int quadrature_point) const;
2570 
2574  const std::vector<double> &get_JxW_values () const;
2575 
2582  const DerivativeForm<1,dim,spacedim> &jacobian (const unsigned int quadrature_point) const;
2583 
2590  const std::vector<DerivativeForm<1,dim,spacedim> > &get_jacobians () const;
2591 
2599  const DerivativeForm<2,dim,spacedim> &jacobian_grad (const unsigned int quadrature_point) const;
2600 
2607  const std::vector<DerivativeForm<2,dim,spacedim> > &get_jacobian_grads () const;
2608 
2617  const Tensor<3,spacedim> &jacobian_pushed_forward_grad (const unsigned int quadrature_point) const;
2618 
2625  const std::vector<Tensor<3,spacedim> > &get_jacobian_pushed_forward_grads () const;
2626 
2635 
2642  const std::vector<DerivativeForm<3,dim,spacedim> > &get_jacobian_2nd_derivatives () const;
2643 
2654 
2661  const std::vector<Tensor<4,spacedim> > &get_jacobian_pushed_forward_2nd_derivatives () const;
2662 
2672 
2679  const std::vector<DerivativeForm<4,dim,spacedim> > &get_jacobian_3rd_derivatives () const;
2680 
2691 
2698  const std::vector<Tensor<5,spacedim> > &get_jacobian_pushed_forward_3rd_derivatives () const;
2699 
2706  const DerivativeForm<1,spacedim,dim> &inverse_jacobian (const unsigned int quadrature_point) const;
2707 
2714  const std::vector<DerivativeForm<1,spacedim,dim> > &get_inverse_jacobians () const;
2715 
2729  const Tensor<1,spacedim> &normal_vector (const unsigned int i) const;
2730 
2741  const std::vector<Tensor<1,spacedim> > &get_all_normal_vectors () const DEAL_II_DEPRECATED;
2742 
2750  const std::vector<Tensor<1,spacedim> > &get_normal_vectors () const;
2751 
2753 
2755 
2756 
2765  const FEValuesViews::Scalar<dim,spacedim> &
2766  operator[] (const FEValuesExtractors::Scalar &scalar) const;
2767 
2776  const FEValuesViews::Vector<dim,spacedim> &
2777  operator[] (const FEValuesExtractors::Vector &vector) const;
2778 
2788  const FEValuesViews::SymmetricTensor<2,dim,spacedim> &
2789  operator[] (const FEValuesExtractors::SymmetricTensor<2> &tensor) const;
2790 
2791 
2800  const FEValuesViews::Tensor<2,dim,spacedim> &
2801  operator[] (const FEValuesExtractors::Tensor<2> &tensor) const;
2802 
2804 
2806 
2807 
2811  const Mapping<dim,spacedim> &get_mapping () const;
2812 
2816  const FiniteElement<dim,spacedim> &get_fe () const;
2817 
2821  UpdateFlags get_update_flags () const;
2822 
2826  const typename Triangulation<dim,spacedim>::cell_iterator get_cell () const;
2827 
2834 
2839  std::size_t memory_consumption () const;
2841 
2842 
2850  char *,
2851  << "You are requesting information from an FEValues/FEFaceValues/FESubfaceValues "
2852  << "object for which this kind of information has not been computed. What "
2853  << "information these objects compute is determined by the update_* flags you "
2854  << "pass to the constructor. Here, the operation you are attempting requires "
2855  << "the <" << arg1 << "> flag to be set, but it was apparently not specified "
2856  << "upon construction.");
2881  int,
2882  << "The shape function with index " << arg1
2883  << " is not primitive, i.e. it is vector-valued and "
2884  << "has more than one non-zero vector component. This "
2885  << "function cannot be called for these shape functions. "
2886  << "Maybe you want to use the same function with the "
2887  << "_component suffix?");
2894 
2895 protected:
2924  class CellIteratorBase;
2925 
2930  template <typename CI> class CellIterator;
2932 
2938  std::unique_ptr<const CellIteratorBase> present_cell;
2939 
2947  boost::signals2::connection tria_listener_refinement;
2948 
2956  boost::signals2::connection tria_listener_mesh_transform;
2957 
2963  void invalidate_present_cell ();
2964 
2974  void
2975  maybe_invalidate_previous_present_cell (const typename Triangulation<dim,spacedim>::cell_iterator &cell);
2976 
2980  const SmartPointer<const Mapping<dim,spacedim>,FEValuesBase<dim,spacedim> > mapping;
2981 
2987  std::unique_ptr<typename Mapping<dim,spacedim>::InternalDataBase> mapping_data;
2988 
2993  ::internal::FEValues::MappingRelatedData<dim, spacedim> mapping_output;
2994 
2995 
3000  const SmartPointer<const FiniteElement<dim,spacedim>,FEValuesBase<dim,spacedim> > fe;
3001 
3007  std::unique_ptr<typename FiniteElement<dim,spacedim>::InternalDataBase> fe_data;
3008 
3013  ::internal::FEValues::FiniteElementRelatedData<dim, spacedim> finite_element_output;
3014 
3015 
3020 
3029  UpdateFlags compute_update_flags (const UpdateFlags update_flags) const;
3030 
3037 
3043  void
3044  check_cell_similarity (const typename Triangulation<dim,spacedim>::cell_iterator &cell);
3045 
3046 private:
3051  FEValuesBase (const FEValuesBase &);
3052 
3057  FEValuesBase &operator= (const FEValuesBase &);
3058 
3063 
3068  template <int, int> friend class FEValuesViews::Scalar;
3069  template <int, int> friend class FEValuesViews::Vector;
3070  template <int, int, int> friend class FEValuesViews::SymmetricTensor;
3071  template <int, int, int> friend class FEValuesViews::Tensor;
3072 };
3073 
3074 
3075 
3086 template <int dim, int spacedim=dim>
3087 class FEValues : public FEValuesBase<dim,spacedim>
3088 {
3089 public:
3094  static const unsigned int integral_dimension = dim;
3095 
3100  FEValues (const Mapping<dim,spacedim> &mapping,
3101  const FiniteElement<dim,spacedim> &fe,
3102  const Quadrature<dim> &quadrature,
3103  const UpdateFlags update_flags);
3104 
3111  const Quadrature<dim> &quadrature,
3112  const UpdateFlags update_flags);
3113 
3120  template <template <int, int> class DoFHandlerType, bool level_dof_access>
3121  void reinit (const TriaIterator<DoFCellAccessor<DoFHandlerType<dim,spacedim>,level_dof_access> > &cell);
3122 
3136  void reinit (const typename Triangulation<dim,spacedim>::cell_iterator &cell);
3137 
3142  const Quadrature<dim> &get_quadrature () const;
3143 
3148  std::size_t memory_consumption () const;
3149 
3164  const FEValues<dim,spacedim> &get_present_fe_values () const;
3165 
3166 private:
3171 
3175  void initialize (const UpdateFlags update_flags);
3176 
3183  void do_reinit ();
3184 };
3185 
3186 
3197 template <int dim, int spacedim=dim>
3198 class FEFaceValuesBase : public FEValuesBase<dim,spacedim>
3199 {
3200 public:
3205  static const unsigned int integral_dimension = dim-1;
3206 
3218  FEFaceValuesBase (const unsigned int n_q_points,
3219  const unsigned int dofs_per_cell,
3220  const UpdateFlags update_flags,
3223  const Quadrature<dim-1>& quadrature);
3224 
3232  const Tensor<1,spacedim> &boundary_form (const unsigned int i) const;
3233 
3240  const std::vector<Tensor<1,spacedim> > &get_boundary_forms () const;
3241 
3246  unsigned int get_face_index() const;
3247 
3252  const Quadrature<dim-1> & get_quadrature () const;
3253 
3258  std::size_t memory_consumption () const;
3259 
3260 protected:
3261 
3266  unsigned int present_face_index;
3267 
3271  const Quadrature<dim-1> quadrature;
3272 };
3273 
3274 
3275 
3290 template <int dim, int spacedim=dim>
3291 class FEFaceValues : public FEFaceValuesBase<dim,spacedim>
3292 {
3293 public:
3298  static const unsigned int dimension = dim;
3299 
3300  static const unsigned int space_dimension = spacedim;
3301 
3306  static const unsigned int integral_dimension = dim-1;
3307 
3314  const Quadrature<dim-1> &quadrature,
3315  const UpdateFlags update_flags);
3316 
3323  const Quadrature<dim-1> &quadrature,
3324  const UpdateFlags update_flags);
3325 
3330  template <template <int, int> class DoFHandlerType, bool level_dof_access>
3331  void reinit (const TriaIterator<DoFCellAccessor<DoFHandlerType<dim,spacedim>,level_dof_access> > &cell,
3332  const unsigned int face_no);
3333 
3347  void reinit (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
3348  const unsigned int face_no);
3349 
3364  const FEFaceValues<dim,spacedim> &get_present_fe_values () const;
3365 private:
3366 
3370  void initialize (const UpdateFlags update_flags);
3371 
3378  void do_reinit (const unsigned int face_no);
3379 };
3380 
3381 
3399 template <int dim, int spacedim=dim>
3400 class FESubfaceValues : public FEFaceValuesBase<dim,spacedim>
3401 {
3402 public:
3406  static const unsigned int dimension = dim;
3407 
3411  static const unsigned int space_dimension = spacedim;
3412 
3417  static const unsigned int integral_dimension = dim-1;
3418 
3425  const Quadrature<dim-1> &face_quadrature,
3426  const UpdateFlags update_flags);
3427 
3434  const Quadrature<dim-1> &face_quadrature,
3435  const UpdateFlags update_flags);
3436 
3443  template <template <int, int> class DoFHandlerType, bool level_dof_access>
3444  void reinit (const TriaIterator<DoFCellAccessor<DoFHandlerType<dim,spacedim>,level_dof_access> > &cell,
3445  const unsigned int face_no,
3446  const unsigned int subface_no);
3447 
3461  void reinit (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
3462  const unsigned int face_no,
3463  const unsigned int subface_no);
3464 
3479  const FESubfaceValues<dim,spacedim> &get_present_fe_values () const;
3480 
3486  DeclException0 (ExcReinitCalledWithBoundaryFace);
3487 
3493  DeclException0 (ExcFaceHasNoSubfaces);
3494 
3495 private:
3496 
3500  void initialize (const UpdateFlags update_flags);
3501 
3508  void do_reinit (const unsigned int face_no,
3509  const unsigned int subface_no);
3510 };
3511 
3512 
3513 #ifndef DOXYGEN
3514 
3515 
3516 /*------------------------ Inline functions: namespace FEValuesViews --------*/
3517 
3518 namespace FEValuesViews
3519 {
3520  template <int dim, int spacedim>
3521  inline
3522  typename Scalar<dim,spacedim>::value_type
3523  Scalar<dim,spacedim>::value (const unsigned int shape_function,
3524  const unsigned int q_point) const
3525  {
3526  Assert (shape_function < fe_values->fe->dofs_per_cell,
3527  ExcIndexRange (shape_function, 0, fe_values->fe->dofs_per_cell));
3528  Assert (fe_values->update_flags & update_values,
3529  ((typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_values"))));
3530 
3531  // an adaptation of the FEValuesBase::shape_value_component function
3532  // except that here we know the component as fixed and we have
3533  // pre-computed and cached a bunch of information. See the comments there.
3534  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
3535  return fe_values->finite_element_output.shape_values(shape_function_data[shape_function]
3536  .row_index,
3537  q_point);
3538  else
3539  return 0;
3540  }
3541 
3542 
3543 
3544 
3545  template <int dim, int spacedim>
3546  inline
3547  typename Scalar<dim,spacedim>::gradient_type
3548  Scalar<dim,spacedim>::gradient (const unsigned int shape_function,
3549  const unsigned int q_point) const
3550  {
3551  Assert (shape_function < fe_values->fe->dofs_per_cell,
3552  ExcIndexRange (shape_function, 0, fe_values->fe->dofs_per_cell));
3553  Assert (fe_values->update_flags & update_gradients,
3554  (typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_gradients")));
3555 
3556  // an adaptation of the
3557  // FEValuesBase::shape_grad_component
3558  // function except that here we know the
3559  // component as fixed and we have
3560  // pre-computed and cached a bunch of
3561  // information. See the comments there.
3562  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
3563  return fe_values->finite_element_output.shape_gradients[shape_function_data[shape_function]
3564  .row_index][q_point];
3565  else
3566  return gradient_type();
3567  }
3568 
3569 
3570 
3571  template <int dim, int spacedim>
3572  inline
3573  typename Scalar<dim,spacedim>::hessian_type
3574  Scalar<dim,spacedim>::hessian (const unsigned int shape_function,
3575  const unsigned int q_point) const
3576  {
3577  Assert (shape_function < fe_values->fe->dofs_per_cell,
3578  ExcIndexRange (shape_function, 0, fe_values->fe->dofs_per_cell));
3579  Assert (fe_values->update_flags & update_hessians,
3580  (typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_hessians")));
3581 
3582  // an adaptation of the
3583  // FEValuesBase::shape_hessian_component
3584  // function except that here we know the
3585  // component as fixed and we have
3586  // pre-computed and cached a bunch of
3587  // information. See the comments there.
3588  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
3589  return fe_values->finite_element_output.shape_hessians[shape_function_data[shape_function].row_index][q_point];
3590  else
3591  return hessian_type();
3592  }
3593 
3594 
3595 
3596  template <int dim, int spacedim>
3597  inline
3598  typename Scalar<dim,spacedim>::third_derivative_type
3599  Scalar<dim,spacedim>::third_derivative (const unsigned int shape_function,
3600  const unsigned int q_point) const
3601  {
3602  Assert (shape_function < fe_values->fe->dofs_per_cell,
3603  ExcIndexRange (shape_function, 0, fe_values->fe->dofs_per_cell));
3604  Assert (fe_values->update_flags & update_3rd_derivatives,
3605  (typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_3rd_derivatives")));
3606 
3607  // an adaptation of the
3608  // FEValuesBase::shape_3rdderivative_component
3609  // function except that here we know the
3610  // component as fixed and we have
3611  // pre-computed and cached a bunch of
3612  // information. See the comments there.
3613  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
3614  return fe_values->finite_element_output.shape_3rd_derivatives[shape_function_data[shape_function].row_index][q_point];
3615  else
3616  return third_derivative_type();
3617  }
3618 
3619 
3620 
3621  template <int dim, int spacedim>
3622  inline
3624  Vector<dim,spacedim>::value (const unsigned int shape_function,
3625  const unsigned int q_point) const
3626  {
3627  Assert (shape_function < fe_values->fe->dofs_per_cell,
3628  ExcIndexRange (shape_function, 0, fe_values->fe->dofs_per_cell));
3629  Assert (fe_values->update_flags & update_values,
3630  (typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_values")));
3631 
3632  // same as for the scalar case except
3633  // that we have one more index
3634  const int snc = shape_function_data[shape_function].single_nonzero_component;
3635  if (snc == -2)
3636  return value_type();
3637  else if (snc != -1)
3638  {
3639  value_type return_value;
3640  return_value[shape_function_data[shape_function].single_nonzero_component_index]
3641  = fe_values->finite_element_output.shape_values(snc,q_point);
3642  return return_value;
3643  }
3644  else
3645  {
3646  value_type return_value;
3647  for (unsigned int d=0; d<dim; ++d)
3648  if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
3649  return_value[d]
3650  = fe_values->finite_element_output.shape_values(shape_function_data[shape_function].row_index[d],q_point);
3651 
3652  return return_value;
3653  }
3654  }
3655 
3656 
3657 
3658  template <int dim, int spacedim>
3659  inline
3661  Vector<dim,spacedim>::gradient (const unsigned int shape_function,
3662  const unsigned int q_point) const
3663  {
3664  Assert (shape_function < fe_values->fe->dofs_per_cell,
3665  ExcIndexRange (shape_function, 0, fe_values->fe->dofs_per_cell));
3666  Assert (fe_values->update_flags & update_gradients,
3667  (typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_gradients")));
3668 
3669  // same as for the scalar case except
3670  // that we have one more index
3671  const int snc = shape_function_data[shape_function].single_nonzero_component;
3672  if (snc == -2)
3673  return gradient_type();
3674  else if (snc != -1)
3675  {
3676  gradient_type return_value;
3677  return_value[shape_function_data[shape_function].single_nonzero_component_index]
3678  = fe_values->finite_element_output.shape_gradients[snc][q_point];
3679  return return_value;
3680  }
3681  else
3682  {
3683  gradient_type return_value;
3684  for (unsigned int d=0; d<dim; ++d)
3685  if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
3686  return_value[d]
3687  = fe_values->finite_element_output.shape_gradients[shape_function_data[shape_function].row_index[d]][q_point];
3688 
3689  return return_value;
3690  }
3691  }
3692 
3693 
3694 
3695  template <int dim, int spacedim>
3696  inline
3698  Vector<dim,spacedim>::divergence (const unsigned int shape_function,
3699  const unsigned int q_point) const
3700  {
3701  // this function works like in
3702  // the case above
3703  Assert (shape_function < fe_values->fe->dofs_per_cell,
3704  ExcIndexRange (shape_function, 0, fe_values->fe->dofs_per_cell));
3705  Assert (fe_values->update_flags & update_gradients,
3706  (typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_gradients")));
3707 
3708  // same as for the scalar case except
3709  // that we have one more index
3710  const int snc = shape_function_data[shape_function].single_nonzero_component;
3711  if (snc == -2)
3712  return divergence_type();
3713  else if (snc != -1)
3714  return
3715  fe_values->finite_element_output.shape_gradients[snc][q_point][shape_function_data[shape_function].single_nonzero_component_index];
3716  else
3717  {
3718  divergence_type return_value = 0;
3719  for (unsigned int d=0; d<dim; ++d)
3720  if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
3721  return_value
3722  += fe_values->finite_element_output.shape_gradients[shape_function_data[shape_function].row_index[d]][q_point][d];
3723 
3724  return return_value;
3725  }
3726  }
3727 
3728 
3729 
3730  template <int dim, int spacedim>
3731  inline
3733  Vector<dim,spacedim>::curl (const unsigned int shape_function, const unsigned int q_point) const
3734  {
3735  // this function works like in the case above
3736 
3737  Assert (shape_function < fe_values->fe->dofs_per_cell,
3738  ExcIndexRange (shape_function, 0, fe_values->fe->dofs_per_cell));
3739  Assert (fe_values->update_flags & update_gradients,
3740  (typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_gradients")));
3741  // same as for the scalar case except that we have one more index
3742  const int snc = shape_function_data[shape_function].single_nonzero_component;
3743 
3744  if (snc == -2)
3745  return curl_type ();
3746 
3747  else
3748  switch (dim)
3749  {
3750  case 1:
3751  {
3752  Assert (false, ExcMessage("Computing the curl in 1d is not a useful operation"));
3753  return curl_type ();
3754  }
3755 
3756  case 2:
3757  {
3758  if (snc != -1)
3759  {
3760  curl_type return_value;
3761 
3762  // the single
3763  // nonzero component
3764  // can only be zero
3765  // or one in 2d
3766  if (shape_function_data[shape_function].single_nonzero_component_index == 0)
3767  return_value[0] = -1.0 * fe_values->finite_element_output.shape_gradients[snc][q_point][1];
3768  else
3769  return_value[0] = fe_values->finite_element_output.shape_gradients[snc][q_point][0];
3770 
3771  return return_value;
3772  }
3773 
3774  else
3775  {
3776  curl_type return_value;
3777 
3778  return_value[0] = 0.0;
3779 
3780  if (shape_function_data[shape_function].is_nonzero_shape_function_component[0])
3781  return_value[0]
3782  -= fe_values->finite_element_output.shape_gradients[shape_function_data[shape_function].row_index[0]][q_point][1];
3783 
3784  if (shape_function_data[shape_function].is_nonzero_shape_function_component[1])
3785  return_value[0]
3786  += fe_values->finite_element_output.shape_gradients[shape_function_data[shape_function].row_index[1]][q_point][0];
3787 
3788  return return_value;
3789  }
3790  }
3791 
3792  case 3:
3793  {
3794  if (snc != -1)
3795  {
3796  curl_type return_value;
3797 
3798  switch (shape_function_data[shape_function].single_nonzero_component_index)
3799  {
3800  case 0:
3801  {
3802  return_value[0] = 0;
3803  return_value[1] = fe_values->finite_element_output.shape_gradients[snc][q_point][2];
3804  return_value[2] = -1.0 * fe_values->finite_element_output.shape_gradients[snc][q_point][1];
3805  return return_value;
3806  }
3807 
3808  case 1:
3809  {
3810  return_value[0] = -1.0 * fe_values->finite_element_output.shape_gradients[snc][q_point][2];
3811  return_value[1] = 0;
3812  return_value[2] = fe_values->finite_element_output.shape_gradients[snc][q_point][0];
3813  return return_value;
3814  }
3815 
3816  default:
3817  {
3818  return_value[0] = fe_values->finite_element_output.shape_gradients[snc][q_point][1];
3819  return_value[1] = -1.0 * fe_values->finite_element_output.shape_gradients[snc][q_point][0];
3820  return_value[2] = 0;
3821  return return_value;
3822  }
3823  }
3824  }
3825 
3826  else
3827  {
3828  curl_type return_value;
3829 
3830  for (unsigned int i = 0; i < dim; ++i)
3831  return_value[i] = 0.0;
3832 
3833  if (shape_function_data[shape_function].is_nonzero_shape_function_component[0])
3834  {
3835  return_value[1]
3836  += fe_values->finite_element_output.shape_gradients[shape_function_data[shape_function].row_index[0]][q_point][2];
3837  return_value[2]
3838  -= fe_values->finite_element_output.shape_gradients[shape_function_data[shape_function].row_index[0]][q_point][1];
3839  }
3840 
3841  if (shape_function_data[shape_function].is_nonzero_shape_function_component[1])
3842  {
3843  return_value[0]
3844  -= fe_values->finite_element_output.shape_gradients[shape_function_data[shape_function].row_index[1]][q_point][2];
3845  return_value[2]
3846  += fe_values->finite_element_output.shape_gradients[shape_function_data[shape_function].row_index[1]][q_point][0];
3847  }
3848 
3849  if (shape_function_data[shape_function].is_nonzero_shape_function_component[2])
3850  {
3851  return_value[0]
3852  += fe_values->finite_element_output.shape_gradients[shape_function_data[shape_function].row_index[2]][q_point][1];
3853  return_value[1]
3854  -= fe_values->finite_element_output.shape_gradients[shape_function_data[shape_function].row_index[2]][q_point][0];
3855  }
3856 
3857  return return_value;
3858  }
3859  }
3860  }
3861  // should not end up here
3862  Assert (false, ExcInternalError());
3863  return curl_type();
3864  }
3865 
3866  template <int dim, int spacedim>
3867  inline
3869  Vector<dim,spacedim>::hessian (const unsigned int shape_function,
3870  const unsigned int q_point) const
3871  {
3872  // this function works like in
3873  // the case above
3874  Assert (shape_function < fe_values->fe->dofs_per_cell,
3875  ExcIndexRange (shape_function, 0, fe_values->fe->dofs_per_cell));
3876  Assert (fe_values->update_flags & update_hessians,
3877  (typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_hessians")));
3878 
3879  // same as for the scalar case except
3880  // that we have one more index
3881  const int snc = shape_function_data[shape_function].single_nonzero_component;
3882  if (snc == -2)
3883  return hessian_type();
3884  else if (snc != -1)
3885  {
3886  hessian_type return_value;
3887  return_value[shape_function_data[shape_function].single_nonzero_component_index]
3888  = fe_values->finite_element_output.shape_hessians[snc][q_point];
3889  return return_value;
3890  }
3891  else
3892  {
3893  hessian_type return_value;
3894  for (unsigned int d=0; d<dim; ++d)
3895  if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
3896  return_value[d]
3897  = fe_values->finite_element_output.shape_hessians[shape_function_data[shape_function].row_index[d]][q_point];
3898 
3899  return return_value;
3900  }
3901  }
3902 
3903  template <int dim, int spacedim>
3904  inline
3906  Vector<dim,spacedim>::third_derivative (const unsigned int shape_function,
3907  const unsigned int q_point) const
3908  {
3909  // this function works like in
3910  // the case above
3911  Assert (shape_function < fe_values->fe->dofs_per_cell,
3912  ExcIndexRange (shape_function, 0, fe_values->fe->dofs_per_cell));
3913  Assert (fe_values->update_flags & update_3rd_derivatives,
3914  (typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_3rd_derivatives")));
3915 
3916  // same as for the scalar case except
3917  // that we have one more index
3918  const int snc = shape_function_data[shape_function].single_nonzero_component;
3919  if (snc == -2)
3920  return third_derivative_type();
3921  else if (snc != -1)
3922  {
3923  third_derivative_type return_value;
3924  return_value[shape_function_data[shape_function].single_nonzero_component_index]
3925  = fe_values->finite_element_output.shape_3rd_derivatives[snc][q_point];
3926  return return_value;
3927  }
3928  else
3929  {
3930  third_derivative_type return_value;
3931  for (unsigned int d=0; d<dim; ++d)
3932  if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
3933  return_value[d]
3934  = fe_values->finite_element_output.shape_3rd_derivatives[shape_function_data[shape_function].row_index[d]][q_point];
3935 
3936  return return_value;
3937  }
3938  }
3939 
3940 
3941  namespace
3942  {
3947  inline
3948  ::SymmetricTensor<2,1>
3949  symmetrize_single_row (const unsigned int n,
3950  const ::Tensor<1,1> &t)
3951  {
3952  Assert (n < 1, ExcIndexRange (n, 0, 1));
3953  (void)n; // removes -Wunused-parameter warning in optimized mode
3954 
3955  const double array[1] = { t[0] };
3956  return ::SymmetricTensor<2,1>(array);
3957  }
3958 
3959 
3960  inline
3961  ::SymmetricTensor<2,2>
3962  symmetrize_single_row (const unsigned int n,
3963  const ::Tensor<1,2> &t)
3964  {
3965  switch (n)
3966  {
3967  case 0:
3968  {
3969  const double array[3] = { t[0], 0, t[1]/2 };
3970  return ::SymmetricTensor<2,2>(array);
3971  }
3972  case 1:
3973  {
3974  const double array[3] = { 0, t[1], t[0]/2 };
3975  return ::SymmetricTensor<2,2>(array);
3976  }
3977  default:
3978  {
3979  Assert (false, ExcIndexRange (n, 0, 2));
3980  return ::SymmetricTensor<2,2>();
3981  }
3982  }
3983  }
3984 
3985 
3986  inline
3987  ::SymmetricTensor<2,3>
3988  symmetrize_single_row (const unsigned int n,
3989  const ::Tensor<1,3> &t)
3990  {
3991  switch (n)
3992  {
3993  case 0:
3994  {
3995  const double array[6] = { t[0], 0, 0, t[1]/2, t[2]/2, 0 };
3996  return ::SymmetricTensor<2,3>(array);
3997  }
3998  case 1:
3999  {
4000  const double array[6] = { 0, t[1], 0, t[0]/2, 0, t[2]/2 };
4001  return ::SymmetricTensor<2,3>(array);
4002  }
4003  case 2:
4004  {
4005  const double array[6] = { 0, 0, t[2], 0, t[0]/2, t[1]/2 };
4006  return ::SymmetricTensor<2,3>(array);
4007  }
4008  default:
4009  {
4010  Assert (false, ExcIndexRange (n, 0, 3));
4011  return ::SymmetricTensor<2,3>();
4012  }
4013  }
4014  }
4015  }
4016 
4017 
4018  template <int dim, int spacedim>
4019  inline
4021  Vector<dim,spacedim>::symmetric_gradient (const unsigned int shape_function,
4022  const unsigned int q_point) const
4023  {
4024  Assert (shape_function < fe_values->fe->dofs_per_cell,
4025  ExcIndexRange (shape_function, 0, fe_values->fe->dofs_per_cell));
4026  Assert (fe_values->update_flags & update_gradients,
4027  (typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_gradients")));
4028 
4029  // same as for the scalar case except
4030  // that we have one more index
4031  const int snc = shape_function_data[shape_function].single_nonzero_component;
4032  if (snc == -2)
4033  return symmetric_gradient_type();
4034  else if (snc != -1)
4035  return symmetrize_single_row (shape_function_data[shape_function].single_nonzero_component_index,
4036  fe_values->finite_element_output.shape_gradients[snc][q_point]);
4037  else
4038  {
4039  gradient_type return_value;
4040  for (unsigned int d=0; d<dim; ++d)
4041  if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
4042  return_value[d]
4043  = fe_values->finite_element_output.shape_gradients[shape_function_data[shape_function].row_index[d]][q_point];
4044 
4045  return symmetrize(return_value);
4046  }
4047  }
4048 
4049 
4050 
4051  template <int dim, int spacedim>
4052  inline
4054  SymmetricTensor<2, dim, spacedim>::value (const unsigned int shape_function,
4055  const unsigned int q_point) const
4056  {
4057  Assert (shape_function < fe_values->fe->dofs_per_cell,
4058  ExcIndexRange (shape_function, 0, fe_values->fe->dofs_per_cell));
4059  Assert (fe_values->update_flags & update_values,
4060  (typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_values")));
4061 
4062  // similar to the vector case where we
4063  // have more then one index and we need
4064  // to convert between unrolled and
4065  // component indexing for tensors
4066  const int snc
4067  = shape_function_data[shape_function].single_nonzero_component;
4068 
4069  if (snc == -2)
4070  {
4071  // shape function is zero for the
4072  // selected components
4073  return value_type();
4074 
4075  }
4076  else if (snc != -1)
4077  {
4078  value_type return_value;
4079  const unsigned int comp =
4080  shape_function_data[shape_function].single_nonzero_component_index;
4081  return_value[value_type::unrolled_to_component_indices(comp)]
4082  = fe_values->finite_element_output.shape_values(snc,q_point);
4083  return return_value;
4084  }
4085  else
4086  {
4087  value_type return_value;
4088  for (unsigned int d = 0; d < value_type::n_independent_components; ++d)
4089  if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
4090  return_value[value_type::unrolled_to_component_indices(d)]
4091  = fe_values->finite_element_output.shape_values(shape_function_data[shape_function].row_index[d],q_point);
4092  return return_value;
4093  }
4094  }
4095 
4096 
4097  template <int dim, int spacedim>
4098  inline
4100  SymmetricTensor<2, dim, spacedim>::divergence(const unsigned int shape_function,
4101  const unsigned int q_point) const
4102  {
4103  Assert (shape_function < fe_values->fe->dofs_per_cell,
4104  ExcIndexRange (shape_function, 0, fe_values->fe->dofs_per_cell));
4105  Assert (fe_values->update_flags & update_gradients,
4106  (typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_gradients")));
4107 
4108  const int snc = shape_function_data[shape_function].single_nonzero_component;
4109 
4110  if (snc == -2)
4111  {
4112  // shape function is zero for the
4113  // selected components
4114  return divergence_type();
4115  }
4116  else if (snc != -1)
4117  {
4118  // we have a single non-zero component
4119  // when the symmetric tensor is
4120  // represented in unrolled form.
4121  // this implies we potentially have
4122  // two non-zero components when
4123  // represented in component form! we
4124  // will only have one non-zero entry
4125  // if the non-zero component lies on
4126  // the diagonal of the tensor.
4127  //
4128  // the divergence of a second-order tensor
4129  // is a first order tensor.
4130  //
4131  // assume the second-order tensor is
4132  // A with components A_{ij}. then
4133  // A_{ij} = A_{ji} and there is only
4134  // one (if diagonal) or two non-zero
4135  // entries in the tensorial
4136  // representation. define the
4137  // divergence as:
4138  // b_i := \dfrac{\partial phi_{ij}}{\partial x_j}.
4139  // (which is incidentally also
4140  // b_j := \dfrac{\partial phi_{ij}}{\partial x_i}).
4141  // In both cases, a sum is implied.
4142  //
4143  // Now, we know the nonzero component
4144  // in unrolled form: it is indicated
4145  // by 'snc'. we can figure out which
4146  // tensor components belong to this:
4147  const unsigned int comp =
4148  shape_function_data[shape_function].single_nonzero_component_index;
4149  const unsigned int ii = value_type::unrolled_to_component_indices(comp)[0];
4150  const unsigned int jj = value_type::unrolled_to_component_indices(comp)[1];
4151 
4152  // given the form of the divergence
4153  // above, if ii=jj there is only a
4154  // single nonzero component of the
4155  // full tensor and the gradient
4156  // equals
4157  // b_ii := \dfrac{\partial phi_{ii,ii}}{\partial x_ii}.
4158  // all other entries of 'b' are zero
4159  //
4160  // on the other hand, if ii!=jj, then
4161  // there are two nonzero entries in
4162  // the full tensor and
4163  // b_ii := \dfrac{\partial phi_{ii,jj}}{\partial x_ii}.
4164  // b_jj := \dfrac{\partial phi_{ii,jj}}{\partial x_jj}.
4165  // again, all other entries of 'b' are
4166  // zero
4167  const ::Tensor<1, spacedim> phi_grad = fe_values->finite_element_output.shape_gradients[snc][q_point];
4168 
4169  divergence_type return_value;
4170  return_value[ii] = phi_grad[jj];
4171 
4172  if (ii != jj)
4173  return_value[jj] = phi_grad[ii];
4174 
4175  return return_value;
4176 
4177  }
4178  else
4179  {
4180  Assert (false, ExcNotImplemented());
4181  divergence_type return_value;
4182  return return_value;
4183  }
4184  }
4185 
4186  template <int dim, int spacedim>
4187  inline
4189  Tensor<2, dim, spacedim>::value (const unsigned int shape_function,
4190  const unsigned int q_point) const
4191  {
4192  Assert (shape_function < fe_values->fe->dofs_per_cell,
4193  ExcIndexRange (shape_function, 0, fe_values->fe->dofs_per_cell));
4194  Assert (fe_values->update_flags & update_values,
4195  (typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_values")));
4196 
4197  // similar to the vector case where we
4198  // have more then one index and we need
4199  // to convert between unrolled and
4200  // component indexing for tensors
4201  const int snc
4202  = shape_function_data[shape_function].single_nonzero_component;
4203 
4204  if (snc == -2)
4205  {
4206  // shape function is zero for the
4207  // selected components
4208  return value_type();
4209 
4210  }
4211  else if (snc != -1)
4212  {
4213  value_type return_value;
4214  const unsigned int comp =
4215  shape_function_data[shape_function].single_nonzero_component_index;
4217  return_value[indices] = fe_values->finite_element_output.shape_values(snc,q_point);
4218  return return_value;
4219  }
4220  else
4221  {
4222  value_type return_value;
4223  for (unsigned int d = 0; d < dim*dim; ++d)
4224  if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
4225  {
4227  return_value[indices]
4228  = fe_values->finite_element_output.shape_values(shape_function_data[shape_function].row_index[d],q_point);
4229  }
4230  return return_value;
4231  }
4232  }
4233 
4234 
4235  template <int dim, int spacedim>
4236  inline
4238  Tensor<2, dim, spacedim>::divergence(const unsigned int shape_function,
4239  const unsigned int q_point) const
4240  {
4241  Assert (shape_function < fe_values->fe->dofs_per_cell,
4242  ExcIndexRange (shape_function, 0, fe_values->fe->dofs_per_cell));
4243  Assert (fe_values->update_flags & update_gradients,
4244  (typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_gradients")));
4245 
4246  const int snc = shape_function_data[shape_function].single_nonzero_component;
4247 
4248  if (snc == -2)
4249  {
4250  // shape function is zero for the
4251  // selected components
4252  return divergence_type();
4253  }
4254  else if (snc != -1)
4255  {
4256  // we have a single non-zero component
4257  // when the tensor is
4258  // represented in unrolled form.
4259  //
4260  // the divergence of a second-order tensor
4261  // is a first order tensor.
4262  //
4263  // assume the second-order tensor is
4264  // A with components A_{ij}.
4265  // divergence as:
4266  // b_j := \dfrac{\partial phi_{ij}}{\partial x_i}.
4267  //
4268  // Now, we know the nonzero component
4269  // in unrolled form: it is indicated
4270  // by 'snc'. we can figure out which
4271  // tensor components belong to this:
4272  const unsigned int comp =
4273  shape_function_data[shape_function].single_nonzero_component_index;
4275  const unsigned int ii = indices[0];
4276  const unsigned int jj = indices[1];
4277 
4278  const ::Tensor<1, spacedim> phi_grad = fe_values->finite_element_output.shape_gradients[snc][q_point];
4279 
4280  divergence_type return_value;
4281  return_value[jj] = phi_grad[ii];
4282 
4283  return return_value;
4284 
4285  }
4286  else
4287  {
4288  Assert (false, ExcNotImplemented());
4289  divergence_type return_value;
4290  return return_value;
4291  }
4292  }
4293 }
4294 
4295 
4296 
4297 /*------------------------ Inline functions: FEValuesBase ------------------------*/
4298 
4299 
4300 
4301 template <int dim, int spacedim>
4302 inline
4305 operator[] (const FEValuesExtractors::Scalar &scalar) const
4306 {
4307  Assert (scalar.component < fe_values_views_cache.scalars.size(),
4308  ExcIndexRange (scalar.component,
4309  0, fe_values_views_cache.scalars.size()));
4310 
4311  return fe_values_views_cache.scalars[scalar.component];
4312 }
4313 
4314 
4315 
4316 template <int dim, int spacedim>
4317 inline
4320 operator[] (const FEValuesExtractors::Vector &vector) const
4321 {
4322  Assert (vector.first_vector_component <
4323  fe_values_views_cache.vectors.size(),
4325  0, fe_values_views_cache.vectors.size()));
4326 
4327  return fe_values_views_cache.vectors[vector.first_vector_component];
4328 }
4329 
4330 template <int dim, int spacedim>
4331 inline
4335 {
4336  Assert (tensor.first_tensor_component <
4337  fe_values_views_cache.symmetric_second_order_tensors.size(),
4339  0, fe_values_views_cache.symmetric_second_order_tensors.size()));
4340 
4341  return fe_values_views_cache.symmetric_second_order_tensors[tensor.first_tensor_component];
4342 }
4343 
4344 template <int dim, int spacedim>
4345 inline
4348 operator[] (const FEValuesExtractors::Tensor<2> &tensor) const
4349 {
4350  Assert (tensor.first_tensor_component <
4351  fe_values_views_cache.second_order_tensors.size(),
4353  0, fe_values_views_cache.second_order_tensors.size()));
4354 
4355  return fe_values_views_cache.second_order_tensors[tensor.first_tensor_component];
4356 }
4357 
4358 
4359 
4360 
4361 template <int dim, int spacedim>
4362 inline
4363 const double &
4364 FEValuesBase<dim,spacedim>::shape_value (const unsigned int i,
4365  const unsigned int j) const
4366 {
4367  Assert (i < fe->dofs_per_cell,
4368  ExcIndexRange (i, 0, fe->dofs_per_cell));
4370  ExcAccessToUninitializedField("update_values"));
4371  Assert (fe->is_primitive (i),
4373 
4374  // if the entire FE is primitive,
4375  // then we can take a short-cut:
4376  if (fe->is_primitive())
4377  return this->finite_element_output.shape_values(i,j);
4378  else
4379  {
4380  // otherwise, use the mapping
4381  // between shape function
4382  // numbers and rows. note that
4383  // by the assertions above, we
4384  // know that this particular
4385  // shape function is primitive,
4386  // so we can call
4387  // system_to_component_index
4388  const unsigned int
4389  row = this->finite_element_output.shape_function_to_row_table[i * fe->n_components() + fe->system_to_component_index(i).first];
4390  return this->finite_element_output.shape_values(row, j);
4391  }
4392 }
4393 
4394 
4395 
4396 template <int dim, int spacedim>
4397 inline
4398 double
4400  const unsigned int j,
4401  const unsigned int component) const
4402 {
4403  Assert (i < fe->dofs_per_cell,
4404  ExcIndexRange (i, 0, fe->dofs_per_cell));
4405  Assert (this->update_flags & update_values,
4406  ExcAccessToUninitializedField("update_values"));
4407  Assert (component < fe->n_components(),
4408  ExcIndexRange(component, 0, fe->n_components()));
4409 
4410  // check whether the shape function
4411  // is non-zero at all within
4412  // this component:
4413  if (fe->get_nonzero_components(i)[component] == false)
4414  return 0;
4415 
4416  // look up the right row in the
4417  // table and take the data from
4418  // there
4419  const unsigned int
4420  row = this->finite_element_output.shape_function_to_row_table[i * fe->n_components() + component];
4421  return this->finite_element_output.shape_values(row, j);
4422 }
4423 
4424 
4425 
4426 template <int dim, int spacedim>
4427 inline
4428 const Tensor<1,spacedim> &
4429 FEValuesBase<dim,spacedim>::shape_grad (const unsigned int i,
4430  const unsigned int j) const
4431 {
4432  Assert (i < fe->dofs_per_cell,
4433  ExcIndexRange (i, 0, fe->dofs_per_cell));
4435  ExcAccessToUninitializedField("update_gradients"));
4436  Assert (fe->is_primitive (i),
4438 
4439  // if the entire FE is primitive,
4440  // then we can take a short-cut:
4441  if (fe->is_primitive())
4442  return this->finite_element_output.shape_gradients[i][j];
4443  else
4444  {
4445  // otherwise, use the mapping
4446  // between shape function
4447  // numbers and rows. note that
4448  // by the assertions above, we
4449  // know that this particular
4450  // shape function is primitive,
4451  // so we can call
4452  // system_to_component_index
4453  const unsigned int
4454  row = this->finite_element_output.shape_function_to_row_table[i * fe->n_components() + fe->system_to_component_index(i).first];
4455  return this->finite_element_output.shape_gradients[row][j];
4456  }
4457 }
4458 
4459 
4460 
4461 template <int dim, int spacedim>
4462 inline
4465  const unsigned int j,
4466  const unsigned int component) const
4467 {
4468  Assert (i < fe->dofs_per_cell,
4469  ExcIndexRange (i, 0, fe->dofs_per_cell));
4470  Assert (this->update_flags & update_gradients,
4471  ExcAccessToUninitializedField("update_gradients"));
4472  Assert (component < fe->n_components(),
4473  ExcIndexRange(component, 0, fe->n_components()));
4474 
4475  // check whether the shape function
4476  // is non-zero at all within
4477  // this component:
4478  if (fe->get_nonzero_components(i)[component] == false)
4479  return Tensor<1,spacedim>();
4480 
4481  // look up the right row in the
4482  // table and take the data from
4483  // there
4484  const unsigned int
4485  row = this->finite_element_output.shape_function_to_row_table[i * fe->n_components() + component];
4486  return this->finite_element_output.shape_gradients[row][j];
4487 }
4488 
4489 
4490 
4491 template <int dim, int spacedim>
4492 inline
4493 const Tensor<2,spacedim> &
4494 FEValuesBase<dim,spacedim>::shape_hessian (const unsigned int i,
4495  const unsigned int j) const
4496 {
4497  Assert (i < fe->dofs_per_cell,
4498  ExcIndexRange (i, 0, fe->dofs_per_cell));
4499  Assert (this->update_flags & update_hessians,
4500  ExcAccessToUninitializedField("update_hessians"));
4501  Assert (fe->is_primitive (i),
4502  ExcShapeFunctionNotPrimitive(i));
4503 
4504  // if the entire FE is primitive,
4505  // then we can take a short-cut:
4506  if (fe->is_primitive())
4507  return this->finite_element_output.shape_hessians[i][j];
4508  else
4509  {
4510  // otherwise, use the mapping
4511  // between shape function
4512  // numbers and rows. note that
4513  // by the assertions above, we
4514  // know that this particular
4515  // shape function is primitive,
4516  // so we can call
4517  // system_to_component_index
4518  const unsigned int
4519  row = this->finite_element_output.shape_function_to_row_table[i * fe->n_components() + fe->system_to_component_index(i).first];
4520  return this->finite_element_output.shape_hessians[row][j];
4521  }
4522 }
4523 
4524 
4525 
4526 template <int dim, int spacedim>
4527 inline
4530  const unsigned int j,
4531  const unsigned int component) const
4532 {
4533  Assert (i < fe->dofs_per_cell,
4534  ExcIndexRange (i, 0, fe->dofs_per_cell));
4535  Assert (this->update_flags & update_hessians,
4536  ExcAccessToUninitializedField("update_hessians"));
4537  Assert (component < fe->n_components(),
4538  ExcIndexRange(component, 0, fe->n_components()));
4539 
4540  // check whether the shape function
4541  // is non-zero at all within
4542  // this component:
4543  if (fe->get_nonzero_components(i)[component] == false)
4544  return Tensor<2,spacedim>();
4545 
4546  // look up the right row in the
4547  // table and take the data from
4548  // there
4549  const unsigned int
4550  row = this->finite_element_output.shape_function_to_row_table[i * fe->n_components() + component];
4551  return this->finite_element_output.shape_hessians[row][j];
4552 }
4553 
4554 
4555 
4556 template <int dim, int spacedim>
4557 inline
4558 const Tensor<3,spacedim> &
4560  const unsigned int j) const
4561 {
4562  Assert (i < fe->dofs_per_cell,
4563  ExcIndexRange (i, 0, fe->dofs_per_cell));
4564  Assert (this->update_flags & update_hessians,
4565  ExcAccessToUninitializedField("update_3rd_derivatives"));
4566  Assert (fe->is_primitive (i),
4567  ExcShapeFunctionNotPrimitive(i));
4568 
4569  // if the entire FE is primitive,
4570  // then we can take a short-cut:
4571  if (fe->is_primitive())
4572  return this->finite_element_output.shape_3rd_derivatives[i][j];
4573  else
4574  {
4575  // otherwise, use the mapping
4576  // between shape function
4577  // numbers and rows. note that
4578  // by the assertions above, we
4579  // know that this particular
4580  // shape function is primitive,
4581  // so we can call
4582  // system_to_component_index
4583  const unsigned int
4584  row = this->finite_element_output.shape_function_to_row_table[i * fe->n_components() + fe->system_to_component_index(i).first];
4585  return this->finite_element_output.shape_3rd_derivatives[row][j];
4586  }
4587 }
4588 
4589 
4590 
4591 template <int dim, int spacedim>
4592 inline
4595  const unsigned int j,
4596  const unsigned int component) const
4597 {
4598  Assert (i < fe->dofs_per_cell,
4599  ExcIndexRange (i, 0, fe->dofs_per_cell));
4600  Assert (this->update_flags & update_hessians,
4601  ExcAccessToUninitializedField("update_3rd_derivatives"));
4602  Assert (component < fe->n_components(),
4603  ExcIndexRange(component, 0, fe->n_components()));
4604 
4605  // check whether the shape function
4606  // is non-zero at all within
4607  // this component:
4608  if (fe->get_nonzero_components(i)[component] == false)
4609  return Tensor<3,spacedim>();
4610 
4611  // look up the right row in the
4612  // table and take the data from
4613  // there
4614  const unsigned int
4615  row = this->finite_element_output.shape_function_to_row_table[i * fe->n_components() + component];
4616  return this->finite_element_output.shape_3rd_derivatives[row][j];
4617 }
4618 
4619 
4620 
4621 template <int dim, int spacedim>
4622 inline
4625 {
4626  return *fe;
4627 }
4628 
4629 
4630 template <int dim, int spacedim>
4631 inline
4632 const Mapping<dim,spacedim> &
4634 {
4635  return *mapping;
4636 }
4637 
4638 
4639 
4640 template <int dim, int spacedim>
4641 inline
4644 {
4645  return this->update_flags;
4646 }
4647 
4648 
4649 
4650 template <int dim, int spacedim>
4651 inline
4652 const std::vector<Point<spacedim> > &
4654 {
4655  Assert (this->update_flags & update_quadrature_points,
4656  ExcAccessToUninitializedField("update_quadrature_points"));
4657  return this->mapping_output.quadrature_points;
4658 }
4659 
4660 
4661 
4662 template <int dim, int spacedim>
4663 inline
4664 const std::vector<double> &
4666 {
4667  Assert (this->update_flags & update_JxW_values,
4668  ExcAccessToUninitializedField("update_JxW_values"));
4669  return this->mapping_output.JxW_values;
4670 }
4671 
4672 
4673 
4674 template <int dim, int spacedim>
4675 inline
4676 const std::vector<DerivativeForm<1,dim,spacedim> > &
4678 {
4679  Assert (this->update_flags & update_jacobians,
4680  ExcAccessToUninitializedField("update_jacobians"));
4681  return this->mapping_output.jacobians;
4682 }
4683 
4684 
4685 
4686 template <int dim, int spacedim>
4687 inline
4688 const std::vector<DerivativeForm<2,dim,spacedim> > &
4690 {
4691  Assert (this->update_flags & update_jacobian_grads,
4692  ExcAccessToUninitializedField("update_jacobians_grads"));
4693  return this->mapping_output.jacobian_grads;
4694 }
4695 
4696 
4697 
4698 template <int dim, int spacedim>
4699 inline
4700 const Tensor<3,spacedim> &
4702 {
4703  Assert (this->update_flags & update_jacobian_pushed_forward_grads,
4704  ExcAccessToUninitializedField("update_jacobian_pushed_forward_grads"));
4705  return this->mapping_output.jacobian_pushed_forward_grads[i];
4706 }
4707 
4708 
4709 
4710 template <int dim, int spacedim>
4711 inline
4712 const std::vector<Tensor<3,spacedim> > &
4714 {
4715  Assert (this->update_flags & update_jacobian_pushed_forward_grads,
4716  ExcAccessToUninitializedField("update_jacobian_pushed_forward_grads"));
4717  return this->mapping_output.jacobian_pushed_forward_grads;
4718 }
4719 
4720 
4721 
4722 template <int dim, int spacedim>
4723 inline
4725 FEValuesBase<dim,spacedim>::jacobian_2nd_derivative (const unsigned int i) const
4726 {
4727  Assert (this->update_flags & update_jacobian_2nd_derivatives,
4728  ExcAccessToUninitializedField("update_jacobian_2nd_derivatives"));
4729  return this->mapping_output.jacobian_2nd_derivatives[i];
4730 }
4731 
4732 
4733 
4734 template <int dim, int spacedim>
4735 inline
4736 const std::vector<DerivativeForm<3,dim,spacedim> > &
4738 {
4739  Assert (this->update_flags & update_jacobian_2nd_derivatives,
4740  ExcAccessToUninitializedField("update_jacobian_2nd_derivatives"));
4741  return this->mapping_output.jacobian_2nd_derivatives;
4742 }
4743 
4744 
4745 
4746 template <int dim, int spacedim>
4747 inline
4748 const Tensor<4,spacedim> &
4750 {
4752  ExcAccessToUninitializedField("update_jacobian_pushed_forward_2nd_derivatives"));
4753  return this->mapping_output.jacobian_pushed_forward_2nd_derivatives[i];
4754 }
4755 
4756 
4757 
4758 template <int dim, int spacedim>
4759 inline
4760 const std::vector<Tensor<4,spacedim> > &
4762 {
4764  ExcAccessToUninitializedField("update_jacobian_pushed_forward_2nd_derivatives"));
4765  return this->mapping_output.jacobian_pushed_forward_2nd_derivatives;
4766 }
4767 
4768 
4769 
4770 template <int dim, int spacedim>
4771 inline
4773 FEValuesBase<dim,spacedim>::jacobian_3rd_derivative (const unsigned int i) const
4774 {
4775  Assert (this->update_flags & update_jacobian_3rd_derivatives,
4776  ExcAccessToUninitializedField("update_jacobian_3rd_derivatives"));
4777  return this->mapping_output.jacobian_3rd_derivatives[i];
4778 }
4779 
4780 
4781 
4782 template <int dim, int spacedim>
4783 inline
4784 const std::vector<DerivativeForm<4,dim,spacedim> > &
4786 {
4787  Assert (this->update_flags & update_jacobian_3rd_derivatives,
4788  ExcAccessToUninitializedField("update_jacobian_3rd_derivatives"));
4789  return this->mapping_output.jacobian_3rd_derivatives;
4790 }
4791 
4792 
4793 
4794 template <int dim, int spacedim>
4795 inline
4796 const Tensor<5,spacedim> &
4798 {
4800  ExcAccessToUninitializedField("update_jacobian_pushed_forward_3rd_derivatives"));
4801  return this->mapping_output.jacobian_pushed_forward_3rd_derivatives[i];
4802 }
4803 
4804 
4805 
4806 template <int dim, int spacedim>
4807 inline
4808 const std::vector<Tensor<5,spacedim> > &
4810 {
4812  ExcAccessToUninitializedField("update_jacobian_pushed_forward_3rd_derivatives"));
4813  return this->mapping_output.jacobian_pushed_forward_3rd_derivatives;
4814 }
4815 
4816 
4817 template <int dim, int spacedim>
4818 inline
4819 const std::vector<DerivativeForm<1,spacedim,dim> > &
4821 {
4822  Assert (this->update_flags & update_inverse_jacobians,
4823  ExcAccessToUninitializedField("update_inverse_jacobians"));
4824  return this->mapping_output.inverse_jacobians;
4825 }
4826 
4827 
4828 
4829 template <int dim, int spacedim>
4830 inline
4831 const Point<spacedim> &
4832 FEValuesBase<dim,spacedim>::quadrature_point (const unsigned int i) const
4833 {
4834  Assert (this->update_flags & update_quadrature_points,
4835  ExcAccessToUninitializedField("update_quadrature_points"));
4836  Assert (i<this->mapping_output.quadrature_points.size(),
4837  ExcIndexRange(i, 0, this->mapping_output.quadrature_points.size()));
4838 
4839  return this->mapping_output.quadrature_points[i];
4840 }
4841 
4842 
4843 
4844 
4845 template <int dim, int spacedim>
4846 inline
4847 double
4848 FEValuesBase<dim,spacedim>::JxW (const unsigned int i) const
4849 {
4850  Assert (this->update_flags & update_JxW_values,
4851  ExcAccessToUninitializedField("update_JxW_values"));
4852  Assert (i<this->mapping_output.JxW_values.size(),
4853  ExcIndexRange(i, 0, this->mapping_output.JxW_values.size()));
4854 
4855  return this->mapping_output.JxW_values[i];
4856 }
4857 
4858 
4859 
4860 template <int dim, int spacedim>
4861 inline
4863 FEValuesBase<dim,spacedim>::jacobian (const unsigned int i) const
4864 {
4865  Assert (this->update_flags & update_jacobians,
4866  ExcAccessToUninitializedField("update_jacobians"));
4867  Assert (i<this->mapping_output.jacobians.size(),
4868  ExcIndexRange(i, 0, this->mapping_output.jacobians.size()));
4869 
4870  return this->mapping_output.jacobians[i];
4871 }
4872 
4873 
4874 
4875 template <int dim, int spacedim>
4876 inline
4878 FEValuesBase<dim,spacedim>::jacobian_grad (const unsigned int i) const
4879 {
4880  Assert (this->update_flags & update_jacobian_grads,
4881  ExcAccessToUninitializedField("update_jacobians_grads"));
4882  Assert (i<this->mapping_output.jacobian_grads.size(),
4883  ExcIndexRange(i, 0, this->mapping_output.jacobian_grads.size()));
4884 
4885  return this->mapping_output.jacobian_grads[i];
4886 }
4887 
4888 
4889 
4890 template <int dim, int spacedim>
4891 inline
4893 FEValuesBase<dim,spacedim>::inverse_jacobian (const unsigned int i) const
4894 {
4895  Assert (this->update_flags & update_inverse_jacobians,
4896  ExcAccessToUninitializedField("update_inverse_jacobians"));
4897  Assert (i<this->mapping_output.inverse_jacobians.size(),
4898  ExcIndexRange(i, 0, this->mapping_output.inverse_jacobians.size()));
4899 
4900  return this->mapping_output.inverse_jacobians[i];
4901 }
4902 
4903 
4904 template <int dim, int spacedim>
4905 inline
4906 const Tensor<1,spacedim> &
4907 FEValuesBase<dim,spacedim>::normal_vector (const unsigned int i) const
4908 {
4909  Assert (this->update_flags & update_normal_vectors,
4910  (typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_normal_vectors")));
4911  Assert (i<this->mapping_output.normal_vectors.size(),
4912  ExcIndexRange(i, 0, this->mapping_output.normal_vectors.size()));
4913 
4914  return this->mapping_output.normal_vectors[i];
4915 }
4916 
4917 
4918 
4919 /*------------------------ Inline functions: FEValues ----------------------------*/
4920 
4921 
4922 template <int dim, int spacedim>
4923 inline
4924 const Quadrature<dim> &
4926 {
4927  return quadrature;
4928 }
4929 
4930 
4931 
4932 template <int dim, int spacedim>
4933 inline
4934 const FEValues<dim,spacedim> &
4936 {
4937  return *this;
4938 }
4939 
4940 
4941 /*------------------------ Inline functions: FEFaceValuesBase --------------------*/
4942 
4943 
4944 template <int dim, int spacedim>
4945 inline
4946 unsigned int
4948 {
4949  return present_face_index;
4950 }
4951 
4952 
4953 /*------------------------ Inline functions: FE*FaceValues --------------------*/
4954 
4955 template <int dim, int spacedim>
4956 inline
4957 const Quadrature<dim-1> &
4959 {
4960  return quadrature;
4961 }
4962 
4963 
4964 
4965 template <int dim, int spacedim>
4966 inline
4969 {
4970  return *this;
4971 }
4972 
4973 
4974 
4975 template <int dim, int spacedim>
4976 inline
4979 {
4980  return *this;
4981 }
4982 
4983 
4984 
4985 template <int dim, int spacedim>
4986 inline
4987 const Tensor<1,spacedim> &
4988 FEFaceValuesBase<dim,spacedim>::boundary_form (const unsigned int i) const
4989 {
4990  Assert (i<this->mapping_output.boundary_forms.size(),
4991  ExcIndexRange(i, 0, this->mapping_output.boundary_forms.size()));
4992  Assert (this->update_flags & update_boundary_forms,
4993  (typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_boundary_forms")));
4994 
4995  return this->mapping_output.boundary_forms[i];
4996 }
4997 
4998 #endif // DOXYGEN
4999 
5000 DEAL_II_NAMESPACE_CLOSE
5001 
5002 #endif
Transformed quadrature weights.
void get_function_hessians_from_local_dof_values(const InputVector &dof_values, std::vector< typename OutputType< typename InputVector::value_type >::hessian_type > &hessians) const
Definition: fe_values.cc:1746
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:3062
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:1678
Shape function values.
double shape_value_component(const unsigned int function_no, const unsigned int point_no, const unsigned int component) const
const Tensor< 1, spacedim > & shape_grad(const unsigned int function_no, const unsigned int quadrature_point) const
ProductType< Number, typename Scalar< dim, spacedim >::third_derivative_type >::type third_derivative_type
Definition: fe_values.h:205
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:1411
const DerivativeForm< 4, dim, spacedim > & jacobian_3rd_derivative(const unsigned int quadrature_point) const
DEAL_II_CUDA_HOST_DEV Tensor()
Definition: tensor.h:825
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:1281
hessian_type hessian(const unsigned int shape_function, const unsigned int q_point) const
Number value_type
Definition: vector.h:140
const FEValuesViews::Scalar< dim, spacedim > & operator[](const FEValuesExtractors::Scalar &scalar) const
ProductType< Number, typename Vector< dim, spacedim >::third_derivative_type >::type third_derivative_type
Definition: fe_values.h:640
CellSimilarity::Similarity cell_similarity
Definition: fe_values.h:3036
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:3520
::Tensor< 1, spacedim > divergence_type
Definition: fe_values.h:1403
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:2038
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:3266
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:1367
std::vector< ShapeFunctionData > shape_function_data
Definition: fe_values.h:1363
void get_function_third_derivatives_from_local_dof_values(const InputVector &dof_values, std::vector< typename OutputType< typename InputVector::value_type >::third_derivative_type > &third_derivatives) const
Definition: fe_values.cc:1479
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
void get_function_symmetric_gradients_from_local_dof_values(const InputVector &dof_values, std::vector< typename OutputType< typename InputVector::value_type >::symmetric_gradient_type > &symmetric_gradients) const
Definition: fe_values.cc:1612
std::vector<::FEValuesViews::Scalar< dim, spacedim > > scalars
Definition: fe_values.h:1644
const unsigned int dofs_per_cell
Definition: fe_values.h:1788
const DerivativeForm< 1, dim, spacedim > & jacobian(const unsigned int quadrature_point) const
const unsigned int component
Definition: fe_values.h:486
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:2624
::Tensor< 3, spacedim > hessian_type
Definition: fe_values.h:578
const Quadrature< dim-1 > quadrature
Definition: fe_values.h:3271
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:571
Outer normal vector, not normalized.
const FEFaceValues< dim, spacedim > & get_present_fe_values() const
std::unique_ptr< const CellIteratorBase > present_cell
Definition: fe_values.h:2931
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:155
STL namespace.
Transformed quadrature points.
void get_function_curls_from_local_dof_values(const InputVector &dof_values, std::vector< typename OutputType< typename InputVector::value_type >::curl_type > &curls) const
Definition: fe_values.cc:1702
const Tensor< 1, spacedim > & boundary_form(const unsigned int i) const
const Triangulation< dim, spacedim >::cell_iterator get_cell() const
Definition: fe_values.cc:3621
ProductType< Number, typename SymmetricTensor< 2, dim, spacedim >::divergence_type >::type divergence_type
Definition: fe_values.h:1166
value_type value(const unsigned int shape_function, const unsigned int q_point) const
::Tensor< 1, spacedim > gradient_type
Definition: fe_values.h:154
unsigned int n_components(const DoFHandler< dim, spacedim > &dh) 1
::SymmetricTensor< 2, spacedim > value_type
Definition: fe_values.h:1136
unsigned int row_index[spacedim]
Definition: fe_values.h:668
void get_function_values_from_local_dof_values(const InputVector &dof_values, std::vector< typename OutputType< typename InputVector::value_type >::value_type > &values) const
Definition: fe_values.cc:1523
void get_function_hessians(const InputVector &fe_function, std::vector< Tensor< 2, spacedim, typename InputVector::value_type > > &hessians) const
Definition: fe_values.cc:3298
void check_cell_similarity(const typename Triangulation< dim, spacedim >::cell_iterator &cell)
Definition: fe_values.cc:3754
static const unsigned int space_dimension
Definition: fe_values.h:1776
::internal::FEValuesViews::Cache< dim, spacedim > fe_values_views_cache
Definition: fe_values.h:3062
std::unique_ptr< typename FiniteElement< dim, spacedim >::InternalDataBase > fe_data
Definition: fe_values.h:3007
static::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
void get_function_laplacians_from_local_dof_values(const InputVector &dof_values, std::vector< typename OutputType< typename InputVector::value_type >::laplacian_type > &laplacians) const
Definition: fe_values.cc:1435
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:1455
ProductType< Number, typename Vector< dim, spacedim >::hessian_type >::type hessian_type
Definition: fe_values.h:634
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:1587
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:1543
const std::vector< double > & get_JxW_values() const
UpdateFlags get_update_flags() const
static::ExceptionBase & ExcInvalidUpdateFlag()
ProductType< Number, typename Vector< dim, spacedim >::gradient_type >::type gradient_type
Definition: fe_values.h:604
void get_function_laplacians_from_local_dof_values(const InputVector &dof_values, std::vector< typename OutputType< typename InputVector::value_type >::laplacian_type > &laplacians) const
Definition: fe_values.cc:1794
const Quadrature< dim > & get_quadrature() const
UpdateFlags compute_update_flags(const UpdateFlags update_flags) const
Definition: fe_values.cc:3673
void get_function_gradients_from_local_dof_values(const InputVector &dof_values, std::vector< typename OutputType< typename InputVector::value_type >::gradient_type > &gradients) const
Definition: fe_values.cc:1567
ProductType< Number, typename Vector< dim, spacedim >::value_type >::type value_type
Definition: fe_values.h:598
ProductType< Number, typename Scalar< dim, spacedim >::gradient_type >::type gradient_type
Definition: fe_values.h:187
void get_function_laplacians(const InputVector &fe_function, std::vector< typename InputVector::value_type > &laplacians) const
Definition: fe_values.cc:3399
const SmartPointer< const FiniteElement< dim, spacedim >, FEValuesBase< dim, spacedim > > fe
Definition: fe_values.h:3000
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:1722
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:534
ProductType< Number, typename Scalar< dim, spacedim >::hessian_type >::type hessian_type
Definition: fe_values.h:199
#define DeclException1(Exception1, type1, outsequence)
Definition: exceptions.h:593
std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > mapping_data
Definition: fe_values.h:2987
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:1323
std::vector< ShapeFunctionData > shape_function_data
Definition: fe_values.h:1620
void get_function_hessians_from_local_dof_values(const InputVector &dof_values, std::vector< typename OutputType< typename InputVector::value_type >::hessian_type > &hessians) const
Definition: fe_values.cc:1391
#define Assert(cond, exc)
Definition: exceptions.h:337
::Tensor< 2, spacedim > value_type
Definition: fe_values.h:1398
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:1815
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:1095
const Quadrature< dim > quadrature
Definition: fe_values.h:3170
const unsigned int first_vector_component
Definition: fe_values.h:1090
std::size_t memory_consumption() const
Definition: fe_values.cc:3653
#define DeclException0(Exception0)
Definition: exceptions.h:570
::internal::FEValues::FiniteElementRelatedData< dim, spacedim > finite_element_output
Definition: fe_values.h:3013
void get_function_gradients_from_local_dof_values(const InputVector &dof_values, std::vector< typename OutputType< typename InputVector::value_type >::gradient_type > &gradients) const
Definition: fe_values.cc:1347
void invalidate_present_cell()
Definition: fe_values.cc:3691
ProductType< Number, typename Scalar< dim, spacedim >::value_type >::type value_type
Definition: fe_values.h:181
const SmartPointer< const Mapping< dim, spacedim >, FEValuesBase< dim, spacedim > > mapping
Definition: fe_values.h:2980
const std::vector< Tensor< 3, spacedim > > & get_jacobian_pushed_forward_grads() const
Vector & operator=(const Vector< dim, spacedim > &)
Definition: fe_values.cc:250
const DerivativeForm< 3, dim, spacedim > & jacobian_2nd_derivative(const unsigned int quadrature_point) const
const std::vector< Tensor< 1, spacedim > > & get_normal_vectors() const
Definition: fe_values.cc:3641
void get_function_divergences_from_local_dof_values(const InputVector &dof_values, std::vector< typename OutputType< typename InputVector::value_type >::divergence_type > &divergences) const
Definition: fe_values.cc:1658
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:657
Tensor< 2, spacedim > shape_hessian_component(const unsigned int function_no, const unsigned int point_no, const unsigned int component) const
ProductType< Number, typename Vector< dim, spacedim >::divergence_type >::type divergence_type
Definition: fe_values.h:616
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:1139
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:1499
ProductType< Number, typename Vector< dim, spacedim >::symmetric_gradient_type >::type symmetric_gradient_type
Definition: fe_values.h:610
const DerivativeForm< 1, spacedim, dim > & inverse_jacobian(const unsigned int quadrature_point) const
const unsigned int n_quadrature_points
Definition: fe_values.h:1781
ProductType< Number, typename Vector< dim, spacedim >::value_type >::type laplacian_type
Definition: fe_values.h:622
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:1084
static const unsigned int dimension
Definition: fe_values.h:1771
::Tensor< 2, spacedim > hessian_type
Definition: fe_values.h:161
boost::signals2::connection tria_listener_mesh_transform
Definition: fe_values.h:2956
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:480
Definition: mpi.h:51
::SymmetricTensor< 2, spacedim > symmetric_gradient_type
Definition: fe_values.h:556
ProductType< Number, typename SymmetricTensor< 2, dim, spacedim >::value_type >::type value_type
Definition: fe_values.h:1160
const std::vector< Tensor< 1, spacedim > > & get_all_normal_vectors() const 1
Definition: fe_values.cc:3630
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:3709
::Tensor< 2, spacedim > gradient_type
Definition: fe_values.h:544
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:1633
const SmartPointer< const FEValuesBase< dim, spacedim > > fe_values
Definition: fe_values.h:1609
Definition: fe.h:33
const SmartPointer< const FEValuesBase< dim, spacedim > > fe_values
Definition: fe_values.h:1352
const DerivativeForm< 2, dim, spacedim > & jacobian_grad(const unsigned int quadrature_point) const
value_type value(const unsigned int shape_function, const unsigned int q_point) const
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:1766
static::ExceptionBase & ExcNotImplemented()
ProductType< Number, typename Vector< dim, spacedim >::curl_type >::type curl_type
Definition: fe_values.h:628
const FiniteElement< dim, spacedim > & get_fe() const
boost::signals2::connection tria_listener_refinement
Definition: fe_values.h:2947
::internal::FEValues::MappingRelatedData< dim, spacedim > mapping_output
Definition: fe_values.h:2993
const std::vector< DerivativeForm< 3, dim, spacedim > > & get_jacobian_2nd_derivatives() const
const Point< spacedim > & quadrature_point(const unsigned int q) const
const std::vector< DerivativeForm< 4, dim, spacedim > > & get_jacobian_3rd_derivatives() const
ProductType< Number, typename Tensor< 2, dim, spacedim >::divergence_type >::type divergence_type
Definition: fe_values.h:1422
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:168
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
ProductType< Number, typename Scalar< dim, spacedim >::value_type >::type laplacian_type
Definition: fe_values.h:193
double JxW(const unsigned int quadrature_point) const
void get_function_values_from_local_dof_values(const InputVector &dof_values, std::vector< typename OutputType< typename InputVector::value_type >::value_type > &values) const
Definition: fe_values.cc:1303
const FEValues< dim, spacedim > & get_present_fe_values() const
CellSimilarity::Similarity get_cell_similarity() const
Definition: fe_values.cc:3808
void get_function_third_derivatives_from_local_dof_values(const InputVector &dof_values, std::vector< typename OutputType< typename InputVector::value_type >::third_derivative_type > &third_derivatives) const
Definition: fe_values.cc:1839
std::vector< ShapeFunctionData > shape_function_data
Definition: fe_values.h:491
UpdateFlags update_flags
Definition: fe_values.h:3019
static::ExceptionBase & ExcInternalError()
ProductType< Number, typename Tensor< 2, dim, spacedim >::value_type >::type value_type
Definition: fe_values.h:1416
void get_function_gradients(const InputVector &fe_function, std::vector< Tensor< 1, spacedim, typename InputVector::value_type > > &gradients) const
Definition: fe_values.cc:3194
::Tensor< 4, spacedim > third_derivative_type
Definition: fe_values.h:585