Reference documentation for deal.II version Git 2618e0f 2017-11-23 17:25:26 +0100
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  DEAL_II_DEPRECATED
2742  const std::vector<Tensor<1,spacedim> > &get_all_normal_vectors () const;
2743 
2751  const std::vector<Tensor<1,spacedim> > &get_normal_vectors () const;
2752 
2754 
2756 
2757 
2767  operator[] (const FEValuesExtractors::Scalar &scalar) const;
2768 
2778  operator[] (const FEValuesExtractors::Vector &vector) const;
2779 
2791 
2792 
2802  operator[] (const FEValuesExtractors::Tensor<2> &tensor) const;
2803 
2805 
2807 
2808 
2812  const Mapping<dim,spacedim> &get_mapping () const;
2813 
2817  const FiniteElement<dim,spacedim> &get_fe () const;
2818 
2822  UpdateFlags get_update_flags () const;
2823 
2827  const typename Triangulation<dim,spacedim>::cell_iterator get_cell () const;
2828 
2835 
2840  std::size_t memory_consumption () const;
2842 
2843 
2851  char *,
2852  << "You are requesting information from an FEValues/FEFaceValues/FESubfaceValues "
2853  << "object for which this kind of information has not been computed. What "
2854  << "information these objects compute is determined by the update_* flags you "
2855  << "pass to the constructor. Here, the operation you are attempting requires "
2856  << "the <" << arg1 << "> flag to be set, but it was apparently not specified "
2857  << "upon construction.");
2858 
2865  "The FiniteElement you provided to FEValues and the FiniteElement that belongs "
2866  "to the DoFHandler that provided the cell iterator do not match.");
2873  int,
2874  << "The shape function with index " << arg1
2875  << " is not primitive, i.e. it is vector-valued and "
2876  << "has more than one non-zero vector component. This "
2877  << "function cannot be called for these shape functions. "
2878  << "Maybe you want to use the same function with the "
2879  << "_component suffix?");
2880 
2887  "The given FiniteElement is not a primitive element but the requested operation "
2888  "only works for those. See FiniteElement::is_primitive() for more information.");
2889 
2890 protected:
2919  class CellIteratorBase;
2920 
2925  template <typename CI> class CellIterator;
2927 
2933  std::unique_ptr<const CellIteratorBase> present_cell;
2934 
2942  boost::signals2::connection tria_listener_refinement;
2943 
2951  boost::signals2::connection tria_listener_mesh_transform;
2952 
2958  void invalidate_present_cell ();
2959 
2969  void
2971 
2976 
2982  std::unique_ptr<typename Mapping<dim,spacedim>::InternalDataBase> mapping_data;
2983 
2989 
2990 
2996 
3002  std::unique_ptr<typename FiniteElement<dim,spacedim>::InternalDataBase> fe_data;
3003 
3009 
3010 
3015 
3024  UpdateFlags compute_update_flags (const UpdateFlags update_flags) const;
3025 
3032 
3038  void
3040 
3041 private:
3046  FEValuesBase (const FEValuesBase &);
3047 
3053 
3058 
3063  template <int, int> friend class FEValuesViews::Scalar;
3064  template <int, int> friend class FEValuesViews::Vector;
3065  template <int, int, int> friend class FEValuesViews::SymmetricTensor;
3066  template <int, int, int> friend class FEValuesViews::Tensor;
3067 };
3068 
3069 
3070 
3081 template <int dim, int spacedim=dim>
3082 class FEValues : public FEValuesBase<dim,spacedim>
3083 {
3084 public:
3089  static const unsigned int integral_dimension = dim;
3090 
3097  const Quadrature<dim> &quadrature,
3098  const UpdateFlags update_flags);
3099 
3106  const Quadrature<dim> &quadrature,
3107  const UpdateFlags update_flags);
3108 
3115  template <template <int, int> class DoFHandlerType, bool level_dof_access>
3116  void reinit (const TriaIterator<DoFCellAccessor<DoFHandlerType<dim,spacedim>,level_dof_access> > &cell);
3117 
3131  void reinit (const typename Triangulation<dim,spacedim>::cell_iterator &cell);
3132 
3137  const Quadrature<dim> &get_quadrature () const;
3138 
3143  std::size_t memory_consumption () const;
3144 
3160 
3161 private:
3166 
3170  void initialize (const UpdateFlags update_flags);
3171 
3178  void do_reinit ();
3179 };
3180 
3181 
3192 template <int dim, int spacedim=dim>
3193 class FEFaceValuesBase : public FEValuesBase<dim,spacedim>
3194 {
3195 public:
3200  static const unsigned int integral_dimension = dim-1;
3201 
3213  FEFaceValuesBase (const unsigned int n_q_points,
3214  const unsigned int dofs_per_cell,
3215  const UpdateFlags update_flags,
3218  const Quadrature<dim-1>& quadrature);
3219 
3227  const Tensor<1,spacedim> &boundary_form (const unsigned int i) const;
3228 
3235  const std::vector<Tensor<1,spacedim> > &get_boundary_forms () const;
3236 
3241  unsigned int get_face_index() const;
3242 
3247  const Quadrature<dim-1> & get_quadrature () const;
3248 
3253  std::size_t memory_consumption () const;
3254 
3255 protected:
3256 
3261  unsigned int present_face_index;
3262 
3266  const Quadrature<dim-1> quadrature;
3267 };
3268 
3269 
3270 
3285 template <int dim, int spacedim=dim>
3286 class FEFaceValues : public FEFaceValuesBase<dim,spacedim>
3287 {
3288 public:
3293  static const unsigned int dimension = dim;
3294 
3295  static const unsigned int space_dimension = spacedim;
3296 
3301  static const unsigned int integral_dimension = dim-1;
3302 
3310  const UpdateFlags update_flags);
3311 
3319  const UpdateFlags update_flags);
3320 
3325  template <template <int, int> class DoFHandlerType, bool level_dof_access>
3326  void reinit (const TriaIterator<DoFCellAccessor<DoFHandlerType<dim,spacedim>,level_dof_access> > &cell,
3327  const unsigned int face_no);
3328 
3342  void reinit (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
3343  const unsigned int face_no);
3344 
3360 private:
3361 
3365  void initialize (const UpdateFlags update_flags);
3366 
3373  void do_reinit (const unsigned int face_no);
3374 };
3375 
3376 
3394 template <int dim, int spacedim=dim>
3395 class FESubfaceValues : public FEFaceValuesBase<dim,spacedim>
3396 {
3397 public:
3401  static const unsigned int dimension = dim;
3402 
3406  static const unsigned int space_dimension = spacedim;
3407 
3412  static const unsigned int integral_dimension = dim-1;
3413 
3420  const Quadrature<dim-1> &face_quadrature,
3421  const UpdateFlags update_flags);
3422 
3429  const Quadrature<dim-1> &face_quadrature,
3430  const UpdateFlags update_flags);
3431 
3438  template <template <int, int> class DoFHandlerType, bool level_dof_access>
3439  void reinit (const TriaIterator<DoFCellAccessor<DoFHandlerType<dim,spacedim>,level_dof_access> > &cell,
3440  const unsigned int face_no,
3441  const unsigned int subface_no);
3442 
3456  void reinit (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
3457  const unsigned int face_no,
3458  const unsigned int subface_no);
3459 
3475 
3482 
3489 
3490 private:
3491 
3495  void initialize (const UpdateFlags update_flags);
3496 
3503  void do_reinit (const unsigned int face_no,
3504  const unsigned int subface_no);
3505 };
3506 
3507 
3508 #ifndef DOXYGEN
3509 
3510 
3511 /*------------------------ Inline functions: namespace FEValuesViews --------*/
3512 
3513 namespace FEValuesViews
3514 {
3515  template <int dim, int spacedim>
3516  inline
3517  typename Scalar<dim,spacedim>::value_type
3518  Scalar<dim,spacedim>::value (const unsigned int shape_function,
3519  const unsigned int q_point) const
3520  {
3521  Assert (shape_function < fe_values->fe->dofs_per_cell,
3522  ExcIndexRange (shape_function, 0, fe_values->fe->dofs_per_cell));
3523  Assert (fe_values->update_flags & update_values,
3524  ((typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_values"))));
3525 
3526  // an adaptation of the FEValuesBase::shape_value_component function
3527  // except that here we know the component as fixed and we have
3528  // pre-computed and cached a bunch of information. See the comments there.
3529  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
3530  return fe_values->finite_element_output.shape_values(shape_function_data[shape_function]
3531  .row_index,
3532  q_point);
3533  else
3534  return 0;
3535  }
3536 
3537 
3538 
3539 
3540  template <int dim, int spacedim>
3541  inline
3542  typename Scalar<dim,spacedim>::gradient_type
3543  Scalar<dim,spacedim>::gradient (const unsigned int shape_function,
3544  const unsigned int q_point) const
3545  {
3546  Assert (shape_function < fe_values->fe->dofs_per_cell,
3547  ExcIndexRange (shape_function, 0, fe_values->fe->dofs_per_cell));
3548  Assert (fe_values->update_flags & update_gradients,
3549  (typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_gradients")));
3550 
3551  // an adaptation of the
3552  // FEValuesBase::shape_grad_component
3553  // function except that here we know the
3554  // component as fixed and we have
3555  // pre-computed and cached a bunch of
3556  // information. See the comments there.
3557  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
3558  return fe_values->finite_element_output.shape_gradients[shape_function_data[shape_function]
3559  .row_index][q_point];
3560  else
3561  return gradient_type();
3562  }
3563 
3564 
3565 
3566  template <int dim, int spacedim>
3567  inline
3568  typename Scalar<dim,spacedim>::hessian_type
3569  Scalar<dim,spacedim>::hessian (const unsigned int shape_function,
3570  const unsigned int q_point) const
3571  {
3572  Assert (shape_function < fe_values->fe->dofs_per_cell,
3573  ExcIndexRange (shape_function, 0, fe_values->fe->dofs_per_cell));
3574  Assert (fe_values->update_flags & update_hessians,
3575  (typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_hessians")));
3576 
3577  // an adaptation of the
3578  // FEValuesBase::shape_hessian_component
3579  // function except that here we know the
3580  // component as fixed and we have
3581  // pre-computed and cached a bunch of
3582  // information. See the comments there.
3583  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
3584  return fe_values->finite_element_output.shape_hessians[shape_function_data[shape_function].row_index][q_point];
3585  else
3586  return hessian_type();
3587  }
3588 
3589 
3590 
3591  template <int dim, int spacedim>
3592  inline
3593  typename Scalar<dim,spacedim>::third_derivative_type
3594  Scalar<dim,spacedim>::third_derivative (const unsigned int shape_function,
3595  const unsigned int q_point) const
3596  {
3597  Assert (shape_function < fe_values->fe->dofs_per_cell,
3598  ExcIndexRange (shape_function, 0, fe_values->fe->dofs_per_cell));
3599  Assert (fe_values->update_flags & update_3rd_derivatives,
3600  (typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_3rd_derivatives")));
3601 
3602  // an adaptation of the
3603  // FEValuesBase::shape_3rdderivative_component
3604  // function except that here we know the
3605  // component as fixed and we have
3606  // pre-computed and cached a bunch of
3607  // information. See the comments there.
3608  if (shape_function_data[shape_function].is_nonzero_shape_function_component)
3609  return fe_values->finite_element_output.shape_3rd_derivatives[shape_function_data[shape_function].row_index][q_point];
3610  else
3611  return third_derivative_type();
3612  }
3613 
3614 
3615 
3616  template <int dim, int spacedim>
3617  inline
3619  Vector<dim,spacedim>::value (const unsigned int shape_function,
3620  const unsigned int q_point) const
3621  {
3622  Assert (shape_function < fe_values->fe->dofs_per_cell,
3623  ExcIndexRange (shape_function, 0, fe_values->fe->dofs_per_cell));
3624  Assert (fe_values->update_flags & update_values,
3625  (typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_values")));
3626 
3627  // same as for the scalar case except
3628  // that we have one more index
3629  const int snc = shape_function_data[shape_function].single_nonzero_component;
3630  if (snc == -2)
3631  return value_type();
3632  else if (snc != -1)
3633  {
3634  value_type return_value;
3635  return_value[shape_function_data[shape_function].single_nonzero_component_index]
3636  = fe_values->finite_element_output.shape_values(snc,q_point);
3637  return return_value;
3638  }
3639  else
3640  {
3641  value_type return_value;
3642  for (unsigned int d=0; d<dim; ++d)
3643  if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
3644  return_value[d]
3645  = fe_values->finite_element_output.shape_values(shape_function_data[shape_function].row_index[d],q_point);
3646 
3647  return return_value;
3648  }
3649  }
3650 
3651 
3652 
3653  template <int dim, int spacedim>
3654  inline
3656  Vector<dim,spacedim>::gradient (const unsigned int shape_function,
3657  const unsigned int q_point) const
3658  {
3659  Assert (shape_function < fe_values->fe->dofs_per_cell,
3660  ExcIndexRange (shape_function, 0, fe_values->fe->dofs_per_cell));
3661  Assert (fe_values->update_flags & update_gradients,
3662  (typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_gradients")));
3663 
3664  // same as for the scalar case except
3665  // that we have one more index
3666  const int snc = shape_function_data[shape_function].single_nonzero_component;
3667  if (snc == -2)
3668  return gradient_type();
3669  else if (snc != -1)
3670  {
3671  gradient_type return_value;
3672  return_value[shape_function_data[shape_function].single_nonzero_component_index]
3673  = fe_values->finite_element_output.shape_gradients[snc][q_point];
3674  return return_value;
3675  }
3676  else
3677  {
3678  gradient_type return_value;
3679  for (unsigned int d=0; d<dim; ++d)
3680  if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
3681  return_value[d]
3682  = fe_values->finite_element_output.shape_gradients[shape_function_data[shape_function].row_index[d]][q_point];
3683 
3684  return return_value;
3685  }
3686  }
3687 
3688 
3689 
3690  template <int dim, int spacedim>
3691  inline
3693  Vector<dim,spacedim>::divergence (const unsigned int shape_function,
3694  const unsigned int q_point) const
3695  {
3696  // this function works like in
3697  // the case above
3698  Assert (shape_function < fe_values->fe->dofs_per_cell,
3699  ExcIndexRange (shape_function, 0, fe_values->fe->dofs_per_cell));
3700  Assert (fe_values->update_flags & update_gradients,
3701  (typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_gradients")));
3702 
3703  // same as for the scalar case except
3704  // that we have one more index
3705  const int snc = shape_function_data[shape_function].single_nonzero_component;
3706  if (snc == -2)
3707  return divergence_type();
3708  else if (snc != -1)
3709  return
3710  fe_values->finite_element_output.shape_gradients[snc][q_point][shape_function_data[shape_function].single_nonzero_component_index];
3711  else
3712  {
3713  divergence_type return_value = 0;
3714  for (unsigned int d=0; d<dim; ++d)
3715  if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
3716  return_value
3717  += fe_values->finite_element_output.shape_gradients[shape_function_data[shape_function].row_index[d]][q_point][d];
3718 
3719  return return_value;
3720  }
3721  }
3722 
3723 
3724 
3725  template <int dim, int spacedim>
3726  inline
3728  Vector<dim,spacedim>::curl (const unsigned int shape_function, const unsigned int q_point) const
3729  {
3730  // this function works like in the case above
3731 
3732  Assert (shape_function < fe_values->fe->dofs_per_cell,
3733  ExcIndexRange (shape_function, 0, fe_values->fe->dofs_per_cell));
3734  Assert (fe_values->update_flags & update_gradients,
3735  (typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_gradients")));
3736  // same as for the scalar case except that we have one more index
3737  const int snc = shape_function_data[shape_function].single_nonzero_component;
3738 
3739  if (snc == -2)
3740  return curl_type ();
3741 
3742  else
3743  switch (dim)
3744  {
3745  case 1:
3746  {
3747  Assert (false, ExcMessage("Computing the curl in 1d is not a useful operation"));
3748  return curl_type ();
3749  }
3750 
3751  case 2:
3752  {
3753  if (snc != -1)
3754  {
3755  curl_type return_value;
3756 
3757  // the single
3758  // nonzero component
3759  // can only be zero
3760  // or one in 2d
3761  if (shape_function_data[shape_function].single_nonzero_component_index == 0)
3762  return_value[0] = -1.0 * fe_values->finite_element_output.shape_gradients[snc][q_point][1];
3763  else
3764  return_value[0] = fe_values->finite_element_output.shape_gradients[snc][q_point][0];
3765 
3766  return return_value;
3767  }
3768 
3769  else
3770  {
3771  curl_type return_value;
3772 
3773  return_value[0] = 0.0;
3774 
3775  if (shape_function_data[shape_function].is_nonzero_shape_function_component[0])
3776  return_value[0]
3777  -= fe_values->finite_element_output.shape_gradients[shape_function_data[shape_function].row_index[0]][q_point][1];
3778 
3779  if (shape_function_data[shape_function].is_nonzero_shape_function_component[1])
3780  return_value[0]
3781  += fe_values->finite_element_output.shape_gradients[shape_function_data[shape_function].row_index[1]][q_point][0];
3782 
3783  return return_value;
3784  }
3785  }
3786 
3787  case 3:
3788  {
3789  if (snc != -1)
3790  {
3791  curl_type return_value;
3792 
3793  switch (shape_function_data[shape_function].single_nonzero_component_index)
3794  {
3795  case 0:
3796  {
3797  return_value[0] = 0;
3798  return_value[1] = fe_values->finite_element_output.shape_gradients[snc][q_point][2];
3799  return_value[2] = -1.0 * fe_values->finite_element_output.shape_gradients[snc][q_point][1];
3800  return return_value;
3801  }
3802 
3803  case 1:
3804  {
3805  return_value[0] = -1.0 * fe_values->finite_element_output.shape_gradients[snc][q_point][2];
3806  return_value[1] = 0;
3807  return_value[2] = fe_values->finite_element_output.shape_gradients[snc][q_point][0];
3808  return return_value;
3809  }
3810 
3811  default:
3812  {
3813  return_value[0] = fe_values->finite_element_output.shape_gradients[snc][q_point][1];
3814  return_value[1] = -1.0 * fe_values->finite_element_output.shape_gradients[snc][q_point][0];
3815  return_value[2] = 0;
3816  return return_value;
3817  }
3818  }
3819  }
3820 
3821  else
3822  {
3823  curl_type return_value;
3824 
3825  for (unsigned int i = 0; i < dim; ++i)
3826  return_value[i] = 0.0;
3827 
3828  if (shape_function_data[shape_function].is_nonzero_shape_function_component[0])
3829  {
3830  return_value[1]
3831  += fe_values->finite_element_output.shape_gradients[shape_function_data[shape_function].row_index[0]][q_point][2];
3832  return_value[2]
3833  -= fe_values->finite_element_output.shape_gradients[shape_function_data[shape_function].row_index[0]][q_point][1];
3834  }
3835 
3836  if (shape_function_data[shape_function].is_nonzero_shape_function_component[1])
3837  {
3838  return_value[0]
3839  -= fe_values->finite_element_output.shape_gradients[shape_function_data[shape_function].row_index[1]][q_point][2];
3840  return_value[2]
3841  += fe_values->finite_element_output.shape_gradients[shape_function_data[shape_function].row_index[1]][q_point][0];
3842  }
3843 
3844  if (shape_function_data[shape_function].is_nonzero_shape_function_component[2])
3845  {
3846  return_value[0]
3847  += fe_values->finite_element_output.shape_gradients[shape_function_data[shape_function].row_index[2]][q_point][1];
3848  return_value[1]
3849  -= fe_values->finite_element_output.shape_gradients[shape_function_data[shape_function].row_index[2]][q_point][0];
3850  }
3851 
3852  return return_value;
3853  }
3854  }
3855  }
3856  // should not end up here
3857  Assert (false, ExcInternalError());
3858  return curl_type();
3859  }
3860 
3861  template <int dim, int spacedim>
3862  inline
3864  Vector<dim,spacedim>::hessian (const unsigned int shape_function,
3865  const unsigned int q_point) const
3866  {
3867  // this function works like in
3868  // the case above
3869  Assert (shape_function < fe_values->fe->dofs_per_cell,
3870  ExcIndexRange (shape_function, 0, fe_values->fe->dofs_per_cell));
3871  Assert (fe_values->update_flags & update_hessians,
3872  (typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_hessians")));
3873 
3874  // same as for the scalar case except
3875  // that we have one more index
3876  const int snc = shape_function_data[shape_function].single_nonzero_component;
3877  if (snc == -2)
3878  return hessian_type();
3879  else if (snc != -1)
3880  {
3881  hessian_type return_value;
3882  return_value[shape_function_data[shape_function].single_nonzero_component_index]
3883  = fe_values->finite_element_output.shape_hessians[snc][q_point];
3884  return return_value;
3885  }
3886  else
3887  {
3888  hessian_type return_value;
3889  for (unsigned int d=0; d<dim; ++d)
3890  if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
3891  return_value[d]
3892  = fe_values->finite_element_output.shape_hessians[shape_function_data[shape_function].row_index[d]][q_point];
3893 
3894  return return_value;
3895  }
3896  }
3897 
3898  template <int dim, int spacedim>
3899  inline
3901  Vector<dim,spacedim>::third_derivative (const unsigned int shape_function,
3902  const unsigned int q_point) const
3903  {
3904  // this function works like in
3905  // the case above
3906  Assert (shape_function < fe_values->fe->dofs_per_cell,
3907  ExcIndexRange (shape_function, 0, fe_values->fe->dofs_per_cell));
3908  Assert (fe_values->update_flags & update_3rd_derivatives,
3909  (typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_3rd_derivatives")));
3910 
3911  // same as for the scalar case except
3912  // that we have one more index
3913  const int snc = shape_function_data[shape_function].single_nonzero_component;
3914  if (snc == -2)
3915  return third_derivative_type();
3916  else if (snc != -1)
3917  {
3918  third_derivative_type return_value;
3919  return_value[shape_function_data[shape_function].single_nonzero_component_index]
3920  = fe_values->finite_element_output.shape_3rd_derivatives[snc][q_point];
3921  return return_value;
3922  }
3923  else
3924  {
3925  third_derivative_type return_value;
3926  for (unsigned int d=0; d<dim; ++d)
3927  if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
3928  return_value[d]
3929  = fe_values->finite_element_output.shape_3rd_derivatives[shape_function_data[shape_function].row_index[d]][q_point];
3930 
3931  return return_value;
3932  }
3933  }
3934 
3935 
3936  namespace
3937  {
3942  inline
3943  ::SymmetricTensor<2,1>
3944  symmetrize_single_row (const unsigned int n,
3945  const ::Tensor<1,1> &t)
3946  {
3947  Assert (n < 1, ExcIndexRange (n, 0, 1));
3948  (void)n; // removes -Wunused-parameter warning in optimized mode
3949 
3950  const double array[1] = { t[0] };
3951  return ::SymmetricTensor<2,1>(array);
3952  }
3953 
3954 
3955  inline
3956  ::SymmetricTensor<2,2>
3957  symmetrize_single_row (const unsigned int n,
3958  const ::Tensor<1,2> &t)
3959  {
3960  switch (n)
3961  {
3962  case 0:
3963  {
3964  const double array[3] = { t[0], 0, t[1]/2 };
3965  return ::SymmetricTensor<2,2>(array);
3966  }
3967  case 1:
3968  {
3969  const double array[3] = { 0, t[1], t[0]/2 };
3970  return ::SymmetricTensor<2,2>(array);
3971  }
3972  default:
3973  {
3974  Assert (false, ExcIndexRange (n, 0, 2));
3975  return ::SymmetricTensor<2,2>();
3976  }
3977  }
3978  }
3979 
3980 
3981  inline
3982  ::SymmetricTensor<2,3>
3983  symmetrize_single_row (const unsigned int n,
3984  const ::Tensor<1,3> &t)
3985  {
3986  switch (n)
3987  {
3988  case 0:
3989  {
3990  const double array[6] = { t[0], 0, 0, t[1]/2, t[2]/2, 0 };
3991  return ::SymmetricTensor<2,3>(array);
3992  }
3993  case 1:
3994  {
3995  const double array[6] = { 0, t[1], 0, t[0]/2, 0, t[2]/2 };
3996  return ::SymmetricTensor<2,3>(array);
3997  }
3998  case 2:
3999  {
4000  const double array[6] = { 0, 0, t[2], 0, t[0]/2, t[1]/2 };
4001  return ::SymmetricTensor<2,3>(array);
4002  }
4003  default:
4004  {
4005  Assert (false, ExcIndexRange (n, 0, 3));
4006  return ::SymmetricTensor<2,3>();
4007  }
4008  }
4009  }
4010  }
4011 
4012 
4013  template <int dim, int spacedim>
4014  inline
4016  Vector<dim,spacedim>::symmetric_gradient (const unsigned int shape_function,
4017  const unsigned int q_point) const
4018  {
4019  Assert (shape_function < fe_values->fe->dofs_per_cell,
4020  ExcIndexRange (shape_function, 0, fe_values->fe->dofs_per_cell));
4021  Assert (fe_values->update_flags & update_gradients,
4022  (typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_gradients")));
4023 
4024  // same as for the scalar case except
4025  // that we have one more index
4026  const int snc = shape_function_data[shape_function].single_nonzero_component;
4027  if (snc == -2)
4028  return symmetric_gradient_type();
4029  else if (snc != -1)
4030  return symmetrize_single_row (shape_function_data[shape_function].single_nonzero_component_index,
4031  fe_values->finite_element_output.shape_gradients[snc][q_point]);
4032  else
4033  {
4034  gradient_type return_value;
4035  for (unsigned int d=0; d<dim; ++d)
4036  if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
4037  return_value[d]
4038  = fe_values->finite_element_output.shape_gradients[shape_function_data[shape_function].row_index[d]][q_point];
4039 
4040  return symmetrize(return_value);
4041  }
4042  }
4043 
4044 
4045 
4046  template <int dim, int spacedim>
4047  inline
4049  SymmetricTensor<2, dim, spacedim>::value (const unsigned int shape_function,
4050  const unsigned int q_point) const
4051  {
4052  Assert (shape_function < fe_values->fe->dofs_per_cell,
4053  ExcIndexRange (shape_function, 0, fe_values->fe->dofs_per_cell));
4054  Assert (fe_values->update_flags & update_values,
4055  (typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_values")));
4056 
4057  // similar to the vector case where we
4058  // have more then one index and we need
4059  // to convert between unrolled and
4060  // component indexing for tensors
4061  const int snc
4062  = shape_function_data[shape_function].single_nonzero_component;
4063 
4064  if (snc == -2)
4065  {
4066  // shape function is zero for the
4067  // selected components
4068  return value_type();
4069 
4070  }
4071  else if (snc != -1)
4072  {
4073  value_type return_value;
4074  const unsigned int comp =
4075  shape_function_data[shape_function].single_nonzero_component_index;
4076  return_value[value_type::unrolled_to_component_indices(comp)]
4077  = fe_values->finite_element_output.shape_values(snc,q_point);
4078  return return_value;
4079  }
4080  else
4081  {
4082  value_type return_value;
4083  for (unsigned int d = 0; d < value_type::n_independent_components; ++d)
4084  if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
4085  return_value[value_type::unrolled_to_component_indices(d)]
4086  = fe_values->finite_element_output.shape_values(shape_function_data[shape_function].row_index[d],q_point);
4087  return return_value;
4088  }
4089  }
4090 
4091 
4092  template <int dim, int spacedim>
4093  inline
4095  SymmetricTensor<2, dim, spacedim>::divergence(const unsigned int shape_function,
4096  const unsigned int q_point) const
4097  {
4098  Assert (shape_function < fe_values->fe->dofs_per_cell,
4099  ExcIndexRange (shape_function, 0, fe_values->fe->dofs_per_cell));
4100  Assert (fe_values->update_flags & update_gradients,
4101  (typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_gradients")));
4102 
4103  const int snc = shape_function_data[shape_function].single_nonzero_component;
4104 
4105  if (snc == -2)
4106  {
4107  // shape function is zero for the
4108  // selected components
4109  return divergence_type();
4110  }
4111  else if (snc != -1)
4112  {
4113  // we have a single non-zero component
4114  // when the symmetric tensor is
4115  // represented in unrolled form.
4116  // this implies we potentially have
4117  // two non-zero components when
4118  // represented in component form! we
4119  // will only have one non-zero entry
4120  // if the non-zero component lies on
4121  // the diagonal of the tensor.
4122  //
4123  // the divergence of a second-order tensor
4124  // is a first order tensor.
4125  //
4126  // assume the second-order tensor is
4127  // A with components A_{ij}. then
4128  // A_{ij} = A_{ji} and there is only
4129  // one (if diagonal) or two non-zero
4130  // entries in the tensorial
4131  // representation. define the
4132  // divergence as:
4133  // b_i := \dfrac{\partial phi_{ij}}{\partial x_j}.
4134  // (which is incidentally also
4135  // b_j := \dfrac{\partial phi_{ij}}{\partial x_i}).
4136  // In both cases, a sum is implied.
4137  //
4138  // Now, we know the nonzero component
4139  // in unrolled form: it is indicated
4140  // by 'snc'. we can figure out which
4141  // tensor components belong to this:
4142  const unsigned int comp =
4143  shape_function_data[shape_function].single_nonzero_component_index;
4144  const unsigned int ii = value_type::unrolled_to_component_indices(comp)[0];
4145  const unsigned int jj = value_type::unrolled_to_component_indices(comp)[1];
4146 
4147  // given the form of the divergence
4148  // above, if ii=jj there is only a
4149  // single nonzero component of the
4150  // full tensor and the gradient
4151  // equals
4152  // b_ii := \dfrac{\partial phi_{ii,ii}}{\partial x_ii}.
4153  // all other entries of 'b' are zero
4154  //
4155  // on the other hand, if ii!=jj, then
4156  // there are two nonzero entries in
4157  // the full tensor and
4158  // b_ii := \dfrac{\partial phi_{ii,jj}}{\partial x_ii}.
4159  // b_jj := \dfrac{\partial phi_{ii,jj}}{\partial x_jj}.
4160  // again, all other entries of 'b' are
4161  // zero
4162  const ::Tensor<1, spacedim> phi_grad = fe_values->finite_element_output.shape_gradients[snc][q_point];
4163 
4164  divergence_type return_value;
4165  return_value[ii] = phi_grad[jj];
4166 
4167  if (ii != jj)
4168  return_value[jj] = phi_grad[ii];
4169 
4170  return return_value;
4171 
4172  }
4173  else
4174  {
4175  Assert (false, ExcNotImplemented());
4176  divergence_type return_value;
4177  return return_value;
4178  }
4179  }
4180 
4181  template <int dim, int spacedim>
4182  inline
4184  Tensor<2, dim, spacedim>::value (const unsigned int shape_function,
4185  const unsigned int q_point) const
4186  {
4187  Assert (shape_function < fe_values->fe->dofs_per_cell,
4188  ExcIndexRange (shape_function, 0, fe_values->fe->dofs_per_cell));
4189  Assert (fe_values->update_flags & update_values,
4190  (typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_values")));
4191 
4192  // similar to the vector case where we
4193  // have more then one index and we need
4194  // to convert between unrolled and
4195  // component indexing for tensors
4196  const int snc
4197  = shape_function_data[shape_function].single_nonzero_component;
4198 
4199  if (snc == -2)
4200  {
4201  // shape function is zero for the
4202  // selected components
4203  return value_type();
4204 
4205  }
4206  else if (snc != -1)
4207  {
4208  value_type return_value;
4209  const unsigned int comp =
4210  shape_function_data[shape_function].single_nonzero_component_index;
4212  return_value[indices] = fe_values->finite_element_output.shape_values(snc,q_point);
4213  return return_value;
4214  }
4215  else
4216  {
4217  value_type return_value;
4218  for (unsigned int d = 0; d < dim*dim; ++d)
4219  if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
4220  {
4222  return_value[indices]
4223  = fe_values->finite_element_output.shape_values(shape_function_data[shape_function].row_index[d],q_point);
4224  }
4225  return return_value;
4226  }
4227  }
4228 
4229 
4230  template <int dim, int spacedim>
4231  inline
4233  Tensor<2, dim, spacedim>::divergence(const unsigned int shape_function,
4234  const unsigned int q_point) const
4235  {
4236  Assert (shape_function < fe_values->fe->dofs_per_cell,
4237  ExcIndexRange (shape_function, 0, fe_values->fe->dofs_per_cell));
4238  Assert (fe_values->update_flags & update_gradients,
4239  (typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_gradients")));
4240 
4241  const int snc = shape_function_data[shape_function].single_nonzero_component;
4242 
4243  if (snc == -2)
4244  {
4245  // shape function is zero for the
4246  // selected components
4247  return divergence_type();
4248  }
4249  else if (snc != -1)
4250  {
4251  // we have a single non-zero component
4252  // when the tensor is
4253  // represented in unrolled form.
4254  //
4255  // the divergence of a second-order tensor
4256  // is a first order tensor.
4257  //
4258  // assume the second-order tensor is
4259  // A with components A_{ij}.
4260  // divergence as:
4261  // b_j := \dfrac{\partial phi_{ij}}{\partial x_i}.
4262  //
4263  // Now, we know the nonzero component
4264  // in unrolled form: it is indicated
4265  // by 'snc'. we can figure out which
4266  // tensor components belong to this:
4267  const unsigned int comp =
4268  shape_function_data[shape_function].single_nonzero_component_index;
4270  const unsigned int ii = indices[0];
4271  const unsigned int jj = indices[1];
4272 
4273  const ::Tensor<1, spacedim> phi_grad = fe_values->finite_element_output.shape_gradients[snc][q_point];
4274 
4275  divergence_type return_value;
4276  return_value[jj] = phi_grad[ii];
4277 
4278  return return_value;
4279 
4280  }
4281  else
4282  {
4283  Assert (false, ExcNotImplemented());
4284  divergence_type return_value;
4285  return return_value;
4286  }
4287  }
4288 }
4289 
4290 
4291 
4292 /*------------------------ Inline functions: FEValuesBase ------------------------*/
4293 
4294 
4295 
4296 template <int dim, int spacedim>
4297 inline
4300 operator[] (const FEValuesExtractors::Scalar &scalar) const
4301 {
4302  Assert (scalar.component < fe_values_views_cache.scalars.size(),
4303  ExcIndexRange (scalar.component,
4304  0, fe_values_views_cache.scalars.size()));
4305 
4306  return fe_values_views_cache.scalars[scalar.component];
4307 }
4308 
4309 
4310 
4311 template <int dim, int spacedim>
4312 inline
4315 operator[] (const FEValuesExtractors::Vector &vector) const
4316 {
4317  Assert (vector.first_vector_component <
4318  fe_values_views_cache.vectors.size(),
4320  0, fe_values_views_cache.vectors.size()));
4321 
4322  return fe_values_views_cache.vectors[vector.first_vector_component];
4323 }
4324 
4325 template <int dim, int spacedim>
4326 inline
4330 {
4331  Assert (tensor.first_tensor_component <
4332  fe_values_views_cache.symmetric_second_order_tensors.size(),
4334  0, fe_values_views_cache.symmetric_second_order_tensors.size()));
4335 
4336  return fe_values_views_cache.symmetric_second_order_tensors[tensor.first_tensor_component];
4337 }
4338 
4339 template <int dim, int spacedim>
4340 inline
4343 operator[] (const FEValuesExtractors::Tensor<2> &tensor) const
4344 {
4345  Assert (tensor.first_tensor_component <
4346  fe_values_views_cache.second_order_tensors.size(),
4348  0, fe_values_views_cache.second_order_tensors.size()));
4349 
4350  return fe_values_views_cache.second_order_tensors[tensor.first_tensor_component];
4351 }
4352 
4353 
4354 
4355 
4356 template <int dim, int spacedim>
4357 inline
4358 const double &
4359 FEValuesBase<dim,spacedim>::shape_value (const unsigned int i,
4360  const unsigned int j) const
4361 {
4362  Assert (i < fe->dofs_per_cell,
4363  ExcIndexRange (i, 0, fe->dofs_per_cell));
4364  Assert (this->update_flags & update_values,
4365  ExcAccessToUninitializedField("update_values"));
4366  Assert (fe->is_primitive (i),
4367  ExcShapeFunctionNotPrimitive(i));
4368  Assert (present_cell.get() != 0,
4369  ExcMessage ("FEValues object is not reinit'ed to any cell"));
4370  // if the entire FE is primitive,
4371  // then we can take a short-cut:
4372  if (fe->is_primitive())
4373  return this->finite_element_output.shape_values(i,j);
4374  else
4375  {
4376  // otherwise, use the mapping
4377  // between shape function
4378  // numbers and rows. note that
4379  // by the assertions above, we
4380  // know that this particular
4381  // shape function is primitive,
4382  // so we can call
4383  // system_to_component_index
4384  const unsigned int
4385  row = this->finite_element_output.shape_function_to_row_table[i * fe->n_components() + fe->system_to_component_index(i).first];
4386  return this->finite_element_output.shape_values(row, j);
4387  }
4388 }
4389 
4390 
4391 
4392 template <int dim, int spacedim>
4393 inline
4394 double
4396  const unsigned int j,
4397  const unsigned int component) const
4398 {
4399  Assert (i < fe->dofs_per_cell,
4400  ExcIndexRange (i, 0, fe->dofs_per_cell));
4401  Assert (this->update_flags & update_values,
4402  ExcAccessToUninitializedField("update_values"));
4403  Assert (component < fe->n_components(),
4404  ExcIndexRange(component, 0, fe->n_components()));
4405  Assert (present_cell.get() != 0,
4406  ExcMessage ("FEValues object is not reinit'ed to any cell"));
4407 
4408  // check whether the shape function
4409  // is non-zero at all within
4410  // this component:
4411  if (fe->get_nonzero_components(i)[component] == false)
4412  return 0;
4413 
4414  // look up the right row in the
4415  // table and take the data from
4416  // there
4417  const unsigned int
4418  row = this->finite_element_output.shape_function_to_row_table[i * fe->n_components() + component];
4419  return this->finite_element_output.shape_values(row, j);
4420 }
4421 
4422 
4423 
4424 template <int dim, int spacedim>
4425 inline
4426 const Tensor<1,spacedim> &
4427 FEValuesBase<dim,spacedim>::shape_grad (const unsigned int i,
4428  const unsigned int j) const
4429 {
4430  Assert (i < fe->dofs_per_cell,
4431  ExcIndexRange (i, 0, fe->dofs_per_cell));
4432  Assert (this->update_flags & update_gradients,
4433  ExcAccessToUninitializedField("update_gradients"));
4434  Assert (fe->is_primitive (i),
4435  ExcShapeFunctionNotPrimitive(i));
4436  Assert (present_cell.get() != 0,
4437  ExcMessage ("FEValues object is not reinit'ed to any cell"));
4438  // if the entire FE is primitive,
4439  // then we can take a short-cut:
4440  if (fe->is_primitive())
4441  return this->finite_element_output.shape_gradients[i][j];
4442  else
4443  {
4444  // otherwise, use the mapping
4445  // between shape function
4446  // numbers and rows. note that
4447  // by the assertions above, we
4448  // know that this particular
4449  // shape function is primitive,
4450  // so we can call
4451  // system_to_component_index
4452  const unsigned int
4453  row = this->finite_element_output.shape_function_to_row_table[i * fe->n_components() + fe->system_to_component_index(i).first];
4454  return this->finite_element_output.shape_gradients[row][j];
4455  }
4456 }
4457 
4458 
4459 
4460 template <int dim, int spacedim>
4461 inline
4464  const unsigned int j,
4465  const unsigned int component) const
4466 {
4467  Assert (i < fe->dofs_per_cell,
4468  ExcIndexRange (i, 0, fe->dofs_per_cell));
4469  Assert (this->update_flags & update_gradients,
4470  ExcAccessToUninitializedField("update_gradients"));
4471  Assert (component < fe->n_components(),
4472  ExcIndexRange(component, 0, fe->n_components()));
4473  Assert (present_cell.get() != 0,
4474  ExcMessage ("FEValues object is not reinit'ed to any cell"));
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  Assert (present_cell.get() != 0,
4504  ExcMessage ("FEValues object is not reinit'ed to any cell"));
4505  // if the entire FE is primitive,
4506  // then we can take a short-cut:
4507  if (fe->is_primitive())
4508  return this->finite_element_output.shape_hessians[i][j];
4509  else
4510  {
4511  // otherwise, use the mapping
4512  // between shape function
4513  // numbers and rows. note that
4514  // by the assertions above, we
4515  // know that this particular
4516  // shape function is primitive,
4517  // so we can call
4518  // system_to_component_index
4519  const unsigned int
4520  row = this->finite_element_output.shape_function_to_row_table[i * fe->n_components() + fe->system_to_component_index(i).first];
4521  return this->finite_element_output.shape_hessians[row][j];
4522  }
4523 }
4524 
4525 
4526 
4527 template <int dim, int spacedim>
4528 inline
4531  const unsigned int j,
4532  const unsigned int component) const
4533 {
4534  Assert (i < fe->dofs_per_cell,
4535  ExcIndexRange (i, 0, fe->dofs_per_cell));
4536  Assert (this->update_flags & update_hessians,
4537  ExcAccessToUninitializedField("update_hessians"));
4538  Assert (component < fe->n_components(),
4539  ExcIndexRange(component, 0, fe->n_components()));
4540  Assert (present_cell.get() != 0,
4541  ExcMessage ("FEValues object is not reinit'ed to any cell"));
4542  // check whether the shape function
4543  // is non-zero at all within
4544  // this component:
4545  if (fe->get_nonzero_components(i)[component] == false)
4546  return Tensor<2,spacedim>();
4547 
4548  // look up the right row in the
4549  // table and take the data from
4550  // there
4551  const unsigned int
4552  row = this->finite_element_output.shape_function_to_row_table[i * fe->n_components() + component];
4553  return this->finite_element_output.shape_hessians[row][j];
4554 }
4555 
4556 
4557 
4558 template <int dim, int spacedim>
4559 inline
4560 const Tensor<3,spacedim> &
4562  const unsigned int j) const
4563 {
4564  Assert (i < fe->dofs_per_cell,
4565  ExcIndexRange (i, 0, fe->dofs_per_cell));
4566  Assert (this->update_flags & update_hessians,
4567  ExcAccessToUninitializedField("update_3rd_derivatives"));
4568  Assert (fe->is_primitive (i),
4569  ExcShapeFunctionNotPrimitive(i));
4570  Assert (present_cell.get() != 0,
4571  ExcMessage ("FEValues object is not reinit'ed to any cell"));
4572  // if the entire FE is primitive,
4573  // then we can take a short-cut:
4574  if (fe->is_primitive())
4575  return this->finite_element_output.shape_3rd_derivatives[i][j];
4576  else
4577  {
4578  // otherwise, use the mapping
4579  // between shape function
4580  // numbers and rows. note that
4581  // by the assertions above, we
4582  // know that this particular
4583  // shape function is primitive,
4584  // so we can call
4585  // system_to_component_index
4586  const unsigned int
4587  row = this->finite_element_output.shape_function_to_row_table[i * fe->n_components() + fe->system_to_component_index(i).first];
4588  return this->finite_element_output.shape_3rd_derivatives[row][j];
4589  }
4590 }
4591 
4592 
4593 
4594 template <int dim, int spacedim>
4595 inline
4598  const unsigned int j,
4599  const unsigned int component) const
4600 {
4601  Assert (i < fe->dofs_per_cell,
4602  ExcIndexRange (i, 0, fe->dofs_per_cell));
4603  Assert (this->update_flags & update_hessians,
4604  ExcAccessToUninitializedField("update_3rd_derivatives"));
4605  Assert (component < fe->n_components(),
4606  ExcIndexRange(component, 0, fe->n_components()));
4607  Assert (present_cell.get() != 0,
4608  ExcMessage ("FEValues object is not reinit'ed to any cell"));
4609  // check whether the shape function
4610  // is non-zero at all within
4611  // this component:
4612  if (fe->get_nonzero_components(i)[component] == false)
4613  return Tensor<3,spacedim>();
4614 
4615  // look up the right row in the
4616  // table and take the data from
4617  // there
4618  const unsigned int
4619  row = this->finite_element_output.shape_function_to_row_table[i * fe->n_components() + component];
4620  return this->finite_element_output.shape_3rd_derivatives[row][j];
4621 }
4622 
4623 
4624 
4625 template <int dim, int spacedim>
4626 inline
4629 {
4630  return *fe;
4631 }
4632 
4633 
4634 template <int dim, int spacedim>
4635 inline
4636 const Mapping<dim,spacedim> &
4638 {
4639  return *mapping;
4640 }
4641 
4642 
4643 
4644 template <int dim, int spacedim>
4645 inline
4648 {
4649  return this->update_flags;
4650 }
4651 
4652 
4653 
4654 template <int dim, int spacedim>
4655 inline
4656 const std::vector<Point<spacedim> > &
4658 {
4659  Assert (this->update_flags & update_quadrature_points,
4660  ExcAccessToUninitializedField("update_quadrature_points"));
4661  Assert (present_cell.get() != 0,
4662  ExcMessage ("FEValues object is not reinit'ed to any cell"));
4663  return this->mapping_output.quadrature_points;
4664 }
4665 
4666 
4667 
4668 template <int dim, int spacedim>
4669 inline
4670 const std::vector<double> &
4672 {
4673  Assert (this->update_flags & update_JxW_values,
4674  ExcAccessToUninitializedField("update_JxW_values"));
4675  Assert (present_cell.get() != 0,
4676  ExcMessage ("FEValues object is not reinit'ed to any cell"));
4677  return this->mapping_output.JxW_values;
4678 }
4679 
4680 
4681 
4682 template <int dim, int spacedim>
4683 inline
4684 const std::vector<DerivativeForm<1,dim,spacedim> > &
4686 {
4687  Assert (this->update_flags & update_jacobians,
4688  ExcAccessToUninitializedField("update_jacobians"));
4689  Assert (present_cell.get() != 0,
4690  ExcMessage ("FEValues object is not reinit'ed to any cell"));
4691  return this->mapping_output.jacobians;
4692 }
4693 
4694 
4695 
4696 template <int dim, int spacedim>
4697 inline
4698 const std::vector<DerivativeForm<2,dim,spacedim> > &
4700 {
4701  Assert (this->update_flags & update_jacobian_grads,
4702  ExcAccessToUninitializedField("update_jacobians_grads"));
4703  Assert (present_cell.get() != 0,
4704  ExcMessage ("FEValues object is not reinit'ed to any cell"));
4705  return this->mapping_output.jacobian_grads;
4706 }
4707 
4708 
4709 
4710 template <int dim, int spacedim>
4711 inline
4712 const Tensor<3,spacedim> &
4714 {
4715  Assert (this->update_flags & update_jacobian_pushed_forward_grads,
4716  ExcAccessToUninitializedField("update_jacobian_pushed_forward_grads"));
4717  Assert (present_cell.get() != 0,
4718  ExcMessage ("FEValues object is not reinit'ed to any cell"));
4719  return this->mapping_output.jacobian_pushed_forward_grads[i];
4720 }
4721 
4722 
4723 
4724 template <int dim, int spacedim>
4725 inline
4726 const std::vector<Tensor<3,spacedim> > &
4728 {
4729  Assert (this->update_flags & update_jacobian_pushed_forward_grads,
4730  ExcAccessToUninitializedField("update_jacobian_pushed_forward_grads"));
4731  Assert (present_cell.get() != 0,
4732  ExcMessage ("FEValues object is not reinit'ed to any cell"));
4733  return this->mapping_output.jacobian_pushed_forward_grads;
4734 }
4735 
4736 
4737 
4738 template <int dim, int spacedim>
4739 inline
4741 FEValuesBase<dim,spacedim>::jacobian_2nd_derivative (const unsigned int i) const
4742 {
4743  Assert (this->update_flags & update_jacobian_2nd_derivatives,
4744  ExcAccessToUninitializedField("update_jacobian_2nd_derivatives"));
4745  Assert (present_cell.get() != 0,
4746  ExcMessage ("FEValues object is not reinit'ed to any cell"));
4747  return this->mapping_output.jacobian_2nd_derivatives[i];
4748 }
4749 
4750 
4751 
4752 template <int dim, int spacedim>
4753 inline
4754 const std::vector<DerivativeForm<3,dim,spacedim> > &
4756 {
4757  Assert (this->update_flags & update_jacobian_2nd_derivatives,
4758  ExcAccessToUninitializedField("update_jacobian_2nd_derivatives"));
4759  Assert (present_cell.get() != 0,
4760  ExcMessage ("FEValues object is not reinit'ed to any cell"));
4761  return this->mapping_output.jacobian_2nd_derivatives;
4762 }
4763 
4764 
4765 
4766 template <int dim, int spacedim>
4767 inline
4768 const Tensor<4,spacedim> &
4770 {
4772  ExcAccessToUninitializedField("update_jacobian_pushed_forward_2nd_derivatives"));
4773  Assert (present_cell.get() != 0,
4774  ExcMessage ("FEValues object is not reinit'ed to any cell"));
4775  return this->mapping_output.jacobian_pushed_forward_2nd_derivatives[i];
4776 }
4777 
4778 
4779 
4780 template <int dim, int spacedim>
4781 inline
4782 const std::vector<Tensor<4,spacedim> > &
4784 {
4785  Assert (this->update_flags & update_jacobian_pushed_forward_2nd_derivatives,
4786  ExcAccessToUninitializedField("update_jacobian_pushed_forward_2nd_derivatives"));
4787  Assert (present_cell.get() != 0,
4788  ExcMessage ("FEValues object is not reinit'ed to any cell"));
4789  return this->mapping_output.jacobian_pushed_forward_2nd_derivatives;
4790 }
4791 
4792 
4793 
4794 template <int dim, int spacedim>
4795 inline
4797 FEValuesBase<dim,spacedim>::jacobian_3rd_derivative (const unsigned int i) const
4798 {
4799  Assert (this->update_flags & update_jacobian_3rd_derivatives,
4800  ExcAccessToUninitializedField("update_jacobian_3rd_derivatives"));
4801  Assert (present_cell.get() != 0,
4802  ExcMessage ("FEValues object is not reinit'ed to any cell"));
4803  return this->mapping_output.jacobian_3rd_derivatives[i];
4804 }
4805 
4806 
4807 
4808 template <int dim, int spacedim>
4809 inline
4810 const std::vector<DerivativeForm<4,dim,spacedim> > &
4812 {
4813  Assert (this->update_flags & update_jacobian_3rd_derivatives,
4814  ExcAccessToUninitializedField("update_jacobian_3rd_derivatives"));
4815  Assert (present_cell.get() != 0,
4816  ExcMessage ("FEValues object is not reinit'ed to any cell"));
4817  return this->mapping_output.jacobian_3rd_derivatives;
4818 }
4819 
4820 
4821 
4822 template <int dim, int spacedim>
4823 inline
4824 const Tensor<5,spacedim> &
4826 {
4828  ExcAccessToUninitializedField("update_jacobian_pushed_forward_3rd_derivatives"));
4829  Assert (present_cell.get() != 0,
4830  ExcMessage ("FEValues object is not reinit'ed to any cell"));
4831  return this->mapping_output.jacobian_pushed_forward_3rd_derivatives[i];
4832 }
4833 
4834 
4835 
4836 template <int dim, int spacedim>
4837 inline
4838 const std::vector<Tensor<5,spacedim> > &
4840 {
4841  Assert (this->update_flags & update_jacobian_pushed_forward_3rd_derivatives,
4842  ExcAccessToUninitializedField("update_jacobian_pushed_forward_3rd_derivatives"));
4843  Assert (present_cell.get() != 0,
4844  ExcMessage ("FEValues object is not reinit'ed to any cell"));
4845  return this->mapping_output.jacobian_pushed_forward_3rd_derivatives;
4846 }
4847 
4848 
4849 template <int dim, int spacedim>
4850 inline
4851 const std::vector<DerivativeForm<1,spacedim,dim> > &
4853 {
4854  Assert (this->update_flags & update_inverse_jacobians,
4855  ExcAccessToUninitializedField("update_inverse_jacobians"));
4856  Assert (present_cell.get() != 0,
4857  ExcMessage ("FEValues object is not reinit'ed to any cell"));
4858  return this->mapping_output.inverse_jacobians;
4859 }
4860 
4861 
4862 
4863 template <int dim, int spacedim>
4864 inline
4865 const Point<spacedim> &
4866 FEValuesBase<dim,spacedim>::quadrature_point (const unsigned int i) const
4867 {
4868  Assert (this->update_flags & update_quadrature_points,
4869  ExcAccessToUninitializedField("update_quadrature_points"));
4870  Assert (i<this->mapping_output.quadrature_points.size(),
4871  ExcIndexRange(i, 0, this->mapping_output.quadrature_points.size()));
4872  Assert (present_cell.get() != 0,
4873  ExcMessage ("FEValues object is not reinit'ed to any cell"));
4874 
4875  return this->mapping_output.quadrature_points[i];
4876 }
4877 
4878 
4879 
4880 
4881 template <int dim, int spacedim>
4882 inline
4883 double
4884 FEValuesBase<dim,spacedim>::JxW (const unsigned int i) const
4885 {
4886  Assert (this->update_flags & update_JxW_values,
4887  ExcAccessToUninitializedField("update_JxW_values"));
4888  Assert (i<this->mapping_output.JxW_values.size(),
4889  ExcIndexRange(i, 0, this->mapping_output.JxW_values.size()));
4890  Assert (present_cell.get() != 0,
4891  ExcMessage ("FEValues object is not reinit'ed to any cell"));
4892 
4893  return this->mapping_output.JxW_values[i];
4894 }
4895 
4896 
4897 
4898 template <int dim, int spacedim>
4899 inline
4901 FEValuesBase<dim,spacedim>::jacobian (const unsigned int i) const
4902 {
4903  Assert (this->update_flags & update_jacobians,
4904  ExcAccessToUninitializedField("update_jacobians"));
4905  Assert (i<this->mapping_output.jacobians.size(),
4906  ExcIndexRange(i, 0, this->mapping_output.jacobians.size()));
4907  Assert (present_cell.get() != 0,
4908  ExcMessage ("FEValues object is not reinit'ed to any cell"));
4909 
4910  return this->mapping_output.jacobians[i];
4911 }
4912 
4913 
4914 
4915 template <int dim, int spacedim>
4916 inline
4918 FEValuesBase<dim,spacedim>::jacobian_grad (const unsigned int i) const
4919 {
4920  Assert (this->update_flags & update_jacobian_grads,
4921  ExcAccessToUninitializedField("update_jacobians_grads"));
4922  Assert (i<this->mapping_output.jacobian_grads.size(),
4923  ExcIndexRange(i, 0, this->mapping_output.jacobian_grads.size()));
4924  Assert (present_cell.get() != 0,
4925  ExcMessage ("FEValues object is not reinit'ed to any cell"));
4926 
4927  return this->mapping_output.jacobian_grads[i];
4928 }
4929 
4930 
4931 
4932 template <int dim, int spacedim>
4933 inline
4935 FEValuesBase<dim,spacedim>::inverse_jacobian (const unsigned int i) const
4936 {
4937  Assert (this->update_flags & update_inverse_jacobians,
4938  ExcAccessToUninitializedField("update_inverse_jacobians"));
4939  Assert (i<this->mapping_output.inverse_jacobians.size(),
4940  ExcIndexRange(i, 0, this->mapping_output.inverse_jacobians.size()));
4941  Assert (present_cell.get() != 0,
4942  ExcMessage ("FEValues object is not reinit'ed to any cell"));
4943 
4944  return this->mapping_output.inverse_jacobians[i];
4945 }
4946 
4947 
4948 template <int dim, int spacedim>
4949 inline
4950 const Tensor<1,spacedim> &
4951 FEValuesBase<dim,spacedim>::normal_vector (const unsigned int i) const
4952 {
4953  Assert (this->update_flags & update_normal_vectors,
4954  (typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_normal_vectors")));
4955  Assert (i<this->mapping_output.normal_vectors.size(),
4956  ExcIndexRange(i, 0, this->mapping_output.normal_vectors.size()));
4957  Assert (present_cell.get() != 0,
4958  ExcMessage ("FEValues object is not reinit'ed to any cell"));
4959 
4960  return this->mapping_output.normal_vectors[i];
4961 }
4962 
4963 
4964 
4965 /*------------------------ Inline functions: FEValues ----------------------------*/
4966 
4967 
4968 template <int dim, int spacedim>
4969 inline
4970 const Quadrature<dim> &
4972 {
4973  return quadrature;
4974 }
4975 
4976 
4977 
4978 template <int dim, int spacedim>
4979 inline
4980 const FEValues<dim,spacedim> &
4982 {
4983  return *this;
4984 }
4985 
4986 
4987 /*------------------------ Inline functions: FEFaceValuesBase --------------------*/
4988 
4989 
4990 template <int dim, int spacedim>
4991 inline
4992 unsigned int
4994 {
4995  return present_face_index;
4996 }
4997 
4998 
4999 /*------------------------ Inline functions: FE*FaceValues --------------------*/
5000 
5001 template <int dim, int spacedim>
5002 inline
5003 const Quadrature<dim-1> &
5005 {
5006  return quadrature;
5007 }
5008 
5009 
5010 
5011 template <int dim, int spacedim>
5012 inline
5015 {
5016  return *this;
5017 }
5018 
5019 
5020 
5021 template <int dim, int spacedim>
5022 inline
5025 {
5026  return *this;
5027 }
5028 
5029 
5030 
5031 template <int dim, int spacedim>
5032 inline
5033 const Tensor<1,spacedim> &
5034 FEFaceValuesBase<dim,spacedim>::boundary_form (const unsigned int i) const
5035 {
5036  Assert (i<this->mapping_output.boundary_forms.size(),
5037  ExcIndexRange(i, 0, this->mapping_output.boundary_forms.size()));
5038  Assert (this->update_flags & update_boundary_forms,
5039  (typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_boundary_forms")));
5040 
5041  return this->mapping_output.boundary_forms[i];
5042 }
5043 
5044 #endif // DOXYGEN
5045 
5046 DEAL_II_NAMESPACE_CLOSE
5047 
5048 #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:1801
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:3141
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:1733
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:1466
const DerivativeForm< 4, dim, spacedim > & jacobian_3rd_derivative(const unsigned int quadrature_point) const
DEAL_II_CUDA_HOST_DEV Tensor()
Definition: tensor.h:825
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:1336
hessian_type hessian(const unsigned int shape_function, const unsigned int q_point) const
Number value_type
Definition: vector.h:110
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:3031
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:3599
::Tensor< 1, spacedim > divergence_type
Definition: fe_values.h:1403
const Tensor< 3, spacedim > & jacobian_pushed_forward_grad(const unsigned int quadrature_point) const
static const unsigned int dimension
Definition: fe_values.h:3401
Cache(const FEValuesBase< dim, spacedim > &fe_values)
Definition: fe_values.cc:2093
std::size_t memory_consumption() const
Definition: fe_values.cc:4155
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:3261
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:1422
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:1534
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:1667
static::ExceptionBase & ExcFaceHasNoSubfaces()
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:2679
::Tensor< 3, spacedim > hessian_type
Definition: fe_values.h:578
const Quadrature< dim-1 > quadrature
Definition: fe_values.h:3266
Volume element.
static const unsigned int dimension
Definition: fe_values.h:3293
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
static::ExceptionBase & ExcReinitCalledWithBoundaryFace()
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:2926
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:156
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:1757
const Tensor< 1, spacedim > & boundary_form(const unsigned int i) const
void do_reinit(const unsigned int face_no)
Definition: fe_values.cc:4301
const Triangulation< dim, spacedim >::cell_iterator get_cell() const
Definition: fe_values.cc:3700
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
::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:1578
void get_function_hessians(const InputVector &fe_function, std::vector< Tensor< 2, spacedim, typename InputVector::value_type > > &hessians) const
Definition: fe_values.cc:3377
void check_cell_similarity(const typename Triangulation< dim, spacedim >::cell_iterator &cell)
Definition: fe_values.cc:3833
static const unsigned int space_dimension
Definition: fe_values.h:1776
::internal::FEValuesViews::Cache< dim, spacedim > fe_values_views_cache
Definition: fe_values.h:3057
std::unique_ptr< typename FiniteElement< dim, spacedim >::InternalDataBase > fe_data
Definition: fe_values.h:3002
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:1490
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:1510
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:1642
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:1598
const std::vector< double > & get_JxW_values() const
void do_reinit(const unsigned int face_no, const unsigned int subface_no)
Definition: fe_values.cc:4487
UpdateFlags get_update_flags() const
static const unsigned int integral_dimension
Definition: fe_values.h:3412
void reinit(const TriaIterator< DoFCellAccessor< DoFHandlerType< dim, spacedim >, level_dof_access > > &cell, const unsigned int face_no, const unsigned int subface_no)
Definition: fe_values.cc:4414
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:1849
const Quadrature< dim > & get_quadrature() const
UpdateFlags compute_update_flags(const UpdateFlags update_flags) const
Definition: fe_values.cc:3752
static const unsigned int integral_dimension
Definition: fe_values.h:3301
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:1622
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:3478
const SmartPointer< const FiniteElement< dim, spacedim >, FEValuesBase< dim, spacedim > > fe
Definition: fe_values.h:2995
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:1777
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
const std::vector< Tensor< 1, spacedim > > & get_boundary_forms() const
Definition: fe_values.cc:4144
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:2982
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:1378
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:1446
#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:1870
gradient_type gradient(const unsigned int shape_function, const unsigned int q_point) const
void reinit(const TriaIterator< DoFCellAccessor< DoFHandlerType< dim, spacedim >, level_dof_access > > &cell, const unsigned int face_no)
Definition: fe_values.cc:4248
Abstract base class for mapping classes.
Definition: dof_tools.h:46
std::vector< ShapeFunctionData > shape_function_data
Definition: fe_values.h:1095
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:582
const Quadrature< dim > quadrature
Definition: fe_values.h:3165
const unsigned int first_vector_component
Definition: fe_values.h:1090
std::size_t memory_consumption() const
Definition: fe_values.cc:3732
#define DeclException0(Exception0)
Definition: exceptions.h:570
::internal::FEValues::FiniteElementRelatedData< dim, spacedim > finite_element_output
Definition: fe_values.h:3008
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:1402
void invalidate_present_cell()
Definition: fe_values.cc:3770
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:2975
const std::vector< Tensor< 3, spacedim > > & get_jacobian_pushed_forward_grads() const
static const unsigned int integral_dimension
Definition: fe_values.h:3089
Vector & operator=(const Vector< dim, spacedim > &)
Definition: fe_values.cc:251
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:3720
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:1713
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
static const unsigned int integral_dimension
Definition: fe_values.h:3200
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.
FEValuesBase & operator=(const FEValuesBase &)
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:1554
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:2951
static const unsigned int space_dimension
Definition: fe_values.h:3406
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:53
::SymmetricTensor< 2, spacedim > symmetric_gradient_type
Definition: fe_values.h:556
FEFaceValuesBase(const unsigned int n_q_points, const unsigned int dofs_per_cell, const UpdateFlags update_flags, const Mapping< dim, spacedim > &mapping, const FiniteElement< dim, spacedim > &fe, const Quadrature< dim-1 > &quadrature)
Definition: fe_values.cc:4124
ProductType< Number, typename SymmetricTensor< 2, dim, spacedim >::value_type >::type value_type
Definition: fe_values.h:1160
const std::vector< Tensor< 4, spacedim > > & get_jacobian_pushed_forward_2nd_derivatives() const
void initialize(const UpdateFlags update_flags)
Definition: fe_values.cc:3945
Shape function gradients.
Normal vectors.
void maybe_invalidate_previous_present_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell)
Definition: fe_values.cc:3788
::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:1688
const SmartPointer< const FEValuesBase< dim, spacedim > > fe_values
Definition: fe_values.h:1609
Definition: fe.h:33
void initialize(const UpdateFlags update_flags)
Definition: fe_values.cc:4373
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
const std::vector< Tensor< 1, spacedim > > & get_all_normal_vectors() const
Definition: fe_values.cc:3709
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:1821
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:2942
::internal::FEValues::MappingRelatedData< dim, spacedim > mapping_output
Definition: fe_values.h:2988
const std::vector< DerivativeForm< 3, dim, spacedim > > & get_jacobian_2nd_derivatives() const
unsigned int n_components(const DoFHandler< dim, spacedim > &dh)
const Point< spacedim > & quadrature_point(const unsigned int q) const
void initialize(const UpdateFlags update_flags)
Definition: fe_values.cc:4206
FEValues(const Mapping< dim, spacedim > &mapping, const FiniteElement< dim, spacedim > &fe, const Quadrature< dim > &quadrature, const UpdateFlags update_flags)
Definition: fe_values.cc:3909
void do_reinit()
Definition: fe_values.cc:4079
const std::vector< DerivativeForm< 4, dim, spacedim > > & get_jacobian_3rd_derivatives() const
void reinit(const TriaIterator< DoFCellAccessor< DoFHandlerType< dim, spacedim >, level_dof_access > > &cell)
Definition: fe_values.cc:4050
ProductType< Number, typename Tensor< 2, dim, spacedim >::divergence_type >::type divergence_type
Definition: fe_values.h:1422
FESubfaceValues(const Mapping< dim, spacedim > &mapping, const FiniteElement< dim, spacedim > &fe, const Quadrature< dim-1 > &face_quadrature, const UpdateFlags update_flags)
Definition: fe_values.cc:4339
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:1358
FEFaceValues(const Mapping< dim, spacedim > &mapping, const FiniteElement< dim, spacedim > &fe, const Quadrature< dim-1 > &quadrature, const UpdateFlags update_flags)
Definition: fe_values.cc:4172
const FEValues< dim, spacedim > & get_present_fe_values() const
CellSimilarity::Similarity get_cell_similarity() const
Definition: fe_values.cc:3887
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:1894
std::vector< ShapeFunctionData > shape_function_data
Definition: fe_values.h:491
UpdateFlags update_flags
Definition: fe_values.h:3014
std::size_t memory_consumption() const
Definition: fe_values.cc:4113
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:3273
::Tensor< 4, spacedim > third_derivative_type
Definition: fe_values.h:585