include/deal.II/fe/fe_values.h

00001 //---------------------------------------------------------------------------
00002 //    @f$Id: fe_values.h 25345 2012-03-31 08:37:04Z bangerth @f$
00003 //
00004 //    Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009, 2010, 2011, 2012 by the deal.II authors
00005 //
00006 //    This file is subject to QPL and may not be  distributed
00007 //    without copyright and license information. Please refer
00008 //    to the file deal.II/doc/license.html for the  text  and
00009 //    further information on this license.
00010 //
00011 //---------------------------------------------------------------------------
00012 #ifndef __deal2__fe_values_h
00013 #define __deal2__fe_values_h
00014 
00015 
00016 #include <deal.II/base/config.h>
00017 #include <deal.II/base/exceptions.h>
00018 #include <deal.II/base/subscriptor.h>
00019 #include <deal.II/base/point.h>
00020 #include <deal.II/base/derivative_form.h>
00021 #include <deal.II/base/symmetric_tensor.h>
00022 #include <deal.II/base/vector_slice.h>
00023 #include <deal.II/base/quadrature.h>
00024 #include <deal.II/base/table.h>
00025 #include <deal.II/grid/tria.h>
00026 #include <deal.II/grid/tria_iterator.h>
00027 #include <deal.II/dofs/dof_handler.h>
00028 #include <deal.II/dofs/dof_accessor.h>
00029 #include <deal.II/hp/dof_handler.h>
00030 #include <deal.II/fe/fe.h>
00031 #include <deal.II/fe/fe_update_flags.h>
00032 #include <deal.II/fe/mapping.h>
00033 #include <deal.II/multigrid/mg_dof_handler.h>
00034 #include <deal.II/multigrid/mg_dof_accessor.h>
00035 
00036 #include <algorithm>
00037 #include <memory>
00038 
00039                                 // dummy include in order to have the
00040                                 // definition of PetscScalar available
00041                                 // without including other PETSc stuff
00042 #ifdef DEAL_II_USE_PETSC
00043 #  include <petsc.h>
00044 #endif
00045 
00046 DEAL_II_NAMESPACE_OPEN
00047 
00048 template <int dim>   class Quadrature;
00049 template <int dim, int spacedim=dim> class FEValuesBase;
00050 
00051 template <typename Number> class Vector;
00052 template <typename Number> class BlockVector;
00053 
00054 
00055 namespace internal
00056 {
00063   template <int dim>
00064   struct CurlType;
00065 
00074   template <>
00075   struct CurlType<1>{
00076       typedef Tensor<1,1>     type;
00077   };
00078 
00087   template <>
00088   struct CurlType<2>{
00089       typedef Tensor<1,1>     type;
00090   };
00091 
00100   template <>
00101   struct CurlType<3>{
00102       typedef Tensor<1,3>     type;
00103   };
00104 }
00105 
00106 
00122 namespace FEValuesExtractors
00123 {
00138   struct Scalar
00139   {
00144       unsigned int component;
00145 
00150       Scalar (const unsigned int component);
00151   };
00152 
00153 
00196   struct Vector
00197   {
00202       unsigned int first_vector_component;
00203 
00210       Vector (const unsigned int first_vector_component);
00211   };
00212 
00213 
00239   template <int rank>
00240   struct SymmetricTensor
00241   {
00246       unsigned int first_tensor_component;
00247 
00254       SymmetricTensor (const unsigned int first_tensor_component);
00255   };
00256 }
00257 
00258 
00259 
00260 
00281 namespace FEValuesViews
00282 {
00297   template <int dim, int spacedim=dim>
00298   class Scalar
00299   {
00300     public:
00308       typedef double        value_type;
00309 
00318       typedef Tensor<1,spacedim> gradient_type;
00319 
00328       typedef Tensor<2,spacedim> hessian_type;
00329 
00334       Scalar ();
00335 
00343       Scalar (const FEValuesBase<dim,spacedim> &fe_values_base,
00344               const unsigned int       component);
00345 
00352       Scalar & operator= (const Scalar<dim,spacedim> &);
00353 
00372       value_type
00373       value (const unsigned int shape_function,
00374              const unsigned int q_point) const;
00375 
00387       gradient_type
00388       gradient (const unsigned int shape_function,
00389                 const unsigned int q_point) const;
00390 
00403       hessian_type
00404       hessian (const unsigned int shape_function,
00405                const unsigned int q_point) const;
00406 
00423       template <class InputVector>
00424       void get_function_values (const InputVector& fe_function,
00425                                 std::vector<value_type>& values) const;
00426 
00443       template <class InputVector>
00444       void get_function_gradients (const InputVector& fe_function,
00445                                    std::vector<gradient_type>& gradients) const;
00446 
00463       template <class InputVector>
00464       void get_function_hessians (const InputVector& fe_function,
00465                                   std::vector<hessian_type>& hessians) const;
00466 
00488       template <class InputVector>
00489       void get_function_laplacians (const InputVector& fe_function,
00490                                     std::vector<value_type>& laplacians) const;
00491 
00492 
00493     private:
00498       const FEValuesBase<dim,spacedim> &fe_values;
00499 
00505       const unsigned int component;
00506 
00513       struct ShapeFunctionData
00514       {
00529           bool is_nonzero_shape_function_component;
00530 
00546           unsigned int row_index;
00547       };
00548 
00553       std::vector<ShapeFunctionData> shape_function_data;
00554   };
00555 
00556 
00557 
00595   template <int dim, int spacedim=dim>
00596   class Vector
00597   {
00598     public:
00606       typedef Tensor<1,spacedim>          value_type;
00607 
00616       typedef Tensor<2,spacedim>          gradient_type;
00617 
00633       typedef ::SymmetricTensor<2,spacedim> symmetric_gradient_type;
00634 
00643       typedef double                 divergence_type;
00644 
00655       typedef typename internal::CurlType<spacedim>::type   curl_type;
00656 
00665       typedef Tensor<3,spacedim>          hessian_type;
00666 
00671       Vector ();
00672 
00685       Vector (const FEValuesBase<dim,spacedim> &fe_values_base,
00686               const unsigned int first_vector_component);
00687 
00694       Vector & operator= (const Vector<dim,spacedim> &);
00695 
00720       value_type
00721       value (const unsigned int shape_function,
00722              const unsigned int q_point) const;
00723 
00735       gradient_type
00736       gradient (const unsigned int shape_function,
00737                 const unsigned int q_point) const;
00738 
00760       symmetric_gradient_type
00761       symmetric_gradient (const unsigned int shape_function,
00762                           const unsigned int q_point) const;
00763 
00775       divergence_type
00776       divergence (const unsigned int shape_function,
00777                   const unsigned int q_point) const;
00778 
00807       curl_type
00808       curl (const unsigned int shape_function,
00809                   const unsigned int q_point) const;
00810 
00823       hessian_type
00824       hessian (const unsigned int shape_function,
00825                const unsigned int q_point) const;
00826 
00843       template <class InputVector>
00844       void get_function_values (const InputVector& fe_function,
00845                                 std::vector<value_type>& values) const;
00846 
00863       template <class InputVector>
00864       void get_function_gradients (const InputVector& fe_function,
00865                                    std::vector<gradient_type>& gradients) const;
00866 
00892       template <class InputVector>
00893       void
00894       get_function_symmetric_gradients (const InputVector& fe_function,
00895                                         std::vector<symmetric_gradient_type>& symmetric_gradients) const;
00896 
00915       template <class InputVector>
00916       void get_function_divergences (const InputVector& fe_function,
00917                                      std::vector<divergence_type>& divergences) const;
00918 
00937       template <class InputVector>
00938       void get_function_curls (const InputVector& fe_function,
00939                                      std::vector<curl_type>& curls) const;
00940 
00957       template <class InputVector>
00958       void get_function_hessians (const InputVector& fe_function,
00959                                   std::vector<hessian_type>& hessians) const;
00960 
00982       template <class InputVector>
00983       void get_function_laplacians (const InputVector& fe_function,
00984                                     std::vector<value_type>& laplacians) const;
00985 
00986     private:
00991       const FEValuesBase<dim,spacedim> &fe_values;
00992 
00998       const unsigned int first_vector_component;
00999 
01006       struct ShapeFunctionData
01007       {
01023           bool is_nonzero_shape_function_component[dim];
01024 
01041           unsigned int row_index[dim];
01042 
01060           int          single_nonzero_component;
01061           unsigned int single_nonzero_component_index;
01062       };
01063 
01068       std::vector<ShapeFunctionData> shape_function_data;
01069   };
01070 
01071 
01072   template <int rank, int dim, int spacedim = dim>
01073   class SymmetricTensor;
01074 
01110   template <int dim, int spacedim>
01111   class SymmetricTensor<2,dim,spacedim>
01112   {
01113     public:
01122       typedef ::SymmetricTensor<2, spacedim> value_type;
01123 
01138       typedef Tensor<1, spacedim> divergence_type;
01139 
01144       SymmetricTensor();
01145 
01160       SymmetricTensor(const FEValuesBase<dim, spacedim> &fe_values_base,
01161                       const unsigned int first_tensor_component);
01162 
01169       SymmetricTensor & operator=(const SymmetricTensor<2, dim, spacedim> &);
01170 
01195       value_type
01196       value (const unsigned int shape_function,
01197              const unsigned int q_point) const;
01198 
01199 
01216       divergence_type
01217       divergence (const unsigned int shape_function,
01218                   const unsigned int q_point) const;
01219 
01236       template <class InputVector>
01237       void get_function_values (const InputVector& fe_function,
01238                                 std::vector<value_type>& values) const;
01239 
01263       template <class InputVector>
01264       void get_function_divergences (const InputVector& fe_function,
01265                                      std::vector<divergence_type>& divergences) const;
01266 
01267     private:
01272       const FEValuesBase<dim, spacedim> &fe_values;
01273 
01279       const unsigned int first_tensor_component;
01280 
01287       struct ShapeFunctionData
01288       {
01304           bool is_nonzero_shape_function_component[value_type::n_independent_components];
01305 
01322           unsigned int row_index[value_type::n_independent_components];
01323 
01341           int single_nonzero_component;
01342           unsigned int single_nonzero_component_index;
01343       };
01344 
01349       std::vector<ShapeFunctionData> shape_function_data;
01350   };
01351 }
01352 
01353 
01354 namespace internal
01355 {
01356   namespace FEValuesViews
01357   {
01370     template <int dim, int spacedim>
01371     struct Cache
01372     {
01378         std::vector<::FEValuesViews::Scalar<dim,spacedim> > scalars;
01379         std::vector<::FEValuesViews::Vector<dim,spacedim> > vectors;
01380         std::vector<::FEValuesViews::SymmetricTensor<2,dim,spacedim> >
01381                 symmetric_second_order_tensors;
01382 
01386         Cache (const FEValuesBase<dim,spacedim> &fe_values);
01387     };
01388   }
01389 }
01390 
01391 
01392 
01393 //TODO: Add access to mapping values to FEValuesBase
01394 //TODO: Several FEValuesBase of a system should share Mapping
01395 
01412 template <int dim, int spacedim=dim>
01413 class FEValuesData
01414 {
01415   public:
01420     void initialize (const unsigned int        n_quadrature_points,
01421                      const FiniteElement<dim,spacedim> &fe,
01422                      const UpdateFlags         flags);
01423 
01462     typedef Table<2,double> ShapeVector;
01463 
01470     typedef std::vector<std::vector<Tensor<1,spacedim> > > GradientVector;
01471 
01476     typedef std::vector<std::vector<Tensor<2,spacedim> > > HessianVector;
01477 
01485     ShapeVector shape_values;
01486 
01495     GradientVector shape_gradients;
01496 
01505     HessianVector shape_hessians;
01506 
01534     std::vector<double>       JxW_values;
01535 
01540     std::vector< DerivativeForm<1,dim,spacedim> > jacobians;
01541 
01546     std::vector<DerivativeForm<2,dim,spacedim> >  jacobian_grads;
01547 
01552     std::vector<DerivativeForm<1,spacedim,dim> > inverse_jacobians;
01553 
01561     std::vector<Point<spacedim> >  quadrature_points;
01562 
01568     std::vector<Point<spacedim> >  normal_vectors;
01569 
01575     std::vector<Tensor<1,spacedim> >  boundary_forms;
01576 
01606     std::vector<unsigned int> shape_function_to_row_table;
01607 
01613     UpdateFlags          update_flags;
01614 };
01615 
01616 
01730 template <int dim, int spacedim>
01731 class FEValuesBase : protected FEValuesData<dim,spacedim>,
01732                      public Subscriptor
01733 {
01734   public:
01739     static const unsigned int dimension = dim;
01740 
01745     static const unsigned int space_dimension = spacedim;
01746 
01750     const unsigned int n_quadrature_points;
01751 
01761     const unsigned int dofs_per_cell;
01762 
01763 
01778     FEValuesBase (const unsigned int n_q_points,
01779                   const unsigned int dofs_per_cell,
01780                   const UpdateFlags update_flags,
01781                   const Mapping<dim,spacedim> &mapping,
01782                   const FiniteElement<dim,spacedim> &fe);
01783 
01784 
01788     ~FEValuesBase ();
01790 
01791 
01825     const double & shape_value (const unsigned int function_no,
01826                                 const unsigned int point_no) const;
01827 
01860     double shape_value_component (const unsigned int function_no,
01861                                   const unsigned int point_no,
01862                                   const unsigned int component) const;
01863 
01897     const Tensor<1,spacedim> &
01898     shape_grad (const unsigned int function,
01899                 const unsigned int quadrature_point) const;
01900 
01929     Tensor<1,spacedim>
01930     shape_grad_component (const unsigned int function_no,
01931                           const unsigned int point_no,
01932                           const unsigned int component) const;
01933 
01967     const Tensor<2,spacedim> &
01968     shape_hessian (const unsigned int function_no,
01969                    const unsigned int point_no) const;
01970 
01974     const Tensor<2,spacedim> &
01975     shape_2nd_derivative (const unsigned int function_no,
01976                           const unsigned int point_no) const;
01977 
01978 
02009     Tensor<2,spacedim>
02010     shape_hessian_component (const unsigned int function_no,
02011                              const unsigned int point_no,
02012                              const unsigned int component) const;
02013 
02017     Tensor<2,spacedim>
02018     shape_2nd_derivative_component (const unsigned int function_no,
02019                                     const unsigned int point_no,
02020                                     const unsigned int component) const;
02021 
02022 
02024 
02025 
02026 
02092     template <class InputVector, typename number>
02093     void get_function_values (const InputVector& fe_function,
02094                               std::vector<number>& values) const;
02095 
02119     template <class InputVector, typename number>
02120     void get_function_values (const InputVector       &fe_function,
02121                               std::vector<Vector<number> > &values) const;
02122 
02154     template <class InputVector, typename number>
02155     void get_function_values (const InputVector& fe_function,
02156                               const VectorSlice<const std::vector<unsigned int> >& indices,
02157                               std::vector<number>& values) const;
02158 
02197     template <class InputVector, typename number>
02198     void get_function_values (const InputVector& fe_function,
02199                               const VectorSlice<const std::vector<unsigned int> >& indices,
02200                               std::vector<Vector<number> >& values) const;
02201 
02202 
02261     template <class InputVector>
02262     void get_function_values (const InputVector& fe_function,
02263                               const VectorSlice<const std::vector<unsigned int> >& indices,
02264                               VectorSlice<std::vector<std::vector<double> > > values,
02265                               const bool quadrature_points_fastest) const;
02266 
02268 
02269 
02270 
02339     template <class InputVector>
02340     void get_function_gradients (const InputVector      &fe_function,
02341                                  std::vector<Tensor<1,spacedim> > &gradients) const;
02342 
02375     template <class InputVector>
02376     void get_function_gradients (const InputVector               &fe_function,
02377                                  std::vector<std::vector<Tensor<1,spacedim> > > &gradients) const;
02378 
02385     template <class InputVector>
02386     void get_function_gradients (const InputVector& fe_function,
02387                                  const VectorSlice<const std::vector<unsigned int> >& indices,
02388                                  std::vector<Tensor<1,spacedim> >& gradients) const;
02389 
02396     template <class InputVector>
02397     void get_function_gradients (const InputVector& fe_function,
02398                                  const VectorSlice<const std::vector<unsigned int> >& indices,
02399                                  VectorSlice<std::vector<std::vector<Tensor<1,spacedim> > > > gradients,
02400                                  bool quadrature_points_fastest = false) const;
02401 
02406     template <class InputVector>
02407     void get_function_grads (const InputVector      &fe_function,
02408                              std::vector<Tensor<1,spacedim> > &gradients) const;
02409 
02414     template <class InputVector>
02415     void get_function_grads (const InputVector               &fe_function,
02416                              std::vector<std::vector<Tensor<1,spacedim> > > &gradients) const;
02417 
02422     template <class InputVector>
02423     void get_function_grads (const InputVector& fe_function,
02424                              const VectorSlice<const std::vector<unsigned int> >& indices,
02425                              std::vector<Tensor<1,spacedim> >& gradients) const;
02426 
02431     template <class InputVector>
02432     void get_function_grads (const InputVector& fe_function,
02433                              const VectorSlice<const std::vector<unsigned int> >& indices,
02434                              std::vector<std::vector<Tensor<1,spacedim> > >& gradients,
02435                              bool quadrature_points_fastest = false) const;
02436 
02438 
02439 
02440 
02511     template <class InputVector>
02512     void
02513     get_function_hessians (const InputVector& fe_function,
02514                            std::vector<Tensor<2,spacedim> >& hessians) const;
02515 
02548     template <class InputVector>
02549     void
02550     get_function_hessians (const InputVector      &fe_function,
02551                            std::vector<std::vector<Tensor<2,spacedim> > > &hessians,
02552                            bool quadrature_points_fastest = false) const;
02553 
02561     template <class InputVector>
02562     void get_function_hessians (
02563       const InputVector& fe_function,
02564       const VectorSlice<const std::vector<unsigned int> >& indices,
02565       std::vector<Tensor<2,spacedim> >& hessians) const;
02566 
02574     template <class InputVector>
02575     void get_function_hessians (
02576       const InputVector& fe_function,
02577       const VectorSlice<const std::vector<unsigned int> >& indices,
02578       VectorSlice<std::vector<std::vector<Tensor<2,spacedim> > > > hessians,
02579       bool quadrature_points_fastest = false) const;
02580 
02584     template <class InputVector>
02585     void
02586     get_function_2nd_derivatives (const InputVector&,
02587                                   std::vector<Tensor<2,spacedim> >&) const;
02588 
02592     template <class InputVector>
02593     void
02594     get_function_2nd_derivatives (const InputVector&,
02595                                   std::vector<std::vector<Tensor<2,spacedim> > >&,
02596                                   bool = false) const;
02597 
02676     template <class InputVector, typename number>
02677     void
02678     get_function_laplacians (const InputVector& fe_function,
02679                              std::vector<number>& laplacians) const;
02680 
02714     template <class InputVector, typename number>
02715     void
02716     get_function_laplacians (const InputVector      &fe_function,
02717                              std::vector<Vector<number> > &laplacians) const;
02718 
02726     template <class InputVector, typename number>
02727     void get_function_laplacians (
02728       const InputVector& fe_function,
02729       const VectorSlice<const std::vector<unsigned int> >& indices,
02730       std::vector<number>& laplacians) const;
02731 
02739     template <class InputVector, typename number>
02740     void get_function_laplacians (
02741       const InputVector& fe_function,
02742       const VectorSlice<const std::vector<unsigned int> >& indices,
02743       std::vector<Vector<number> >& laplacians) const;
02744 
02752     template <class InputVector, typename number>
02753     void get_function_laplacians (
02754       const InputVector& fe_function,
02755       const VectorSlice<const std::vector<unsigned int> >& indices,
02756       std::vector<std::vector<number> >& laplacians,
02757       bool quadrature_points_fastest = false) const;
02759 
02761 
02762 
02767     const Point<spacedim> & quadrature_point (const unsigned int i) const;
02768 
02773     const std::vector<Point<spacedim> > & get_quadrature_points () const;
02774 
02798     double JxW (const unsigned int quadrature_point) const;
02799 
02804     const std::vector<double> & get_JxW_values () const;
02805 
02812     const DerivativeForm<1,dim,spacedim> & jacobian (const unsigned int quadrature_point) const;
02813 
02818     const std::vector<DerivativeForm<1,dim,spacedim> > & get_jacobians () const;
02819 
02827     const DerivativeForm<2,dim,spacedim> & jacobian_grad (const unsigned int quadrature_point) const;
02828 
02834     const std::vector<DerivativeForm<2,dim,spacedim> > & get_jacobian_grads () const;
02835 
02842     const DerivativeForm<1,spacedim,dim> & inverse_jacobian (const unsigned int quadrature_point) const;
02843 
02849     const std::vector<DerivativeForm<1,spacedim,dim> > & get_inverse_jacobians () const;
02864     const Point<spacedim> & normal_vector (const unsigned int i) const;
02865 
02876     const std::vector<Point<spacedim> > & get_normal_vectors () const;
02877 
02885     void transform (std::vector<Tensor<1,spacedim> >& transformed,
02886                     const std::vector<Tensor<1,dim> >& original,
02887                     MappingType mapping) const;
02888 
02898     const Point<spacedim> & cell_normal_vector (const unsigned int i) const;
02899 
02908     const std::vector<Point<spacedim> > & get_cell_normal_vectors () const;
02909 
02911 
02913 
02914 
02925     const FEValuesViews::Scalar<dim,spacedim> &
02926     operator[] (const FEValuesExtractors::Scalar &scalar) const;
02927 
02938     const FEValuesViews::Vector<dim,spacedim> &
02939     operator[] (const FEValuesExtractors::Vector &vector) const;
02940 
02952     const FEValuesViews::SymmetricTensor<2,dim,spacedim> &
02953     operator[] (const FEValuesExtractors::SymmetricTensor<2> &tensor) const;
02954 
02956 
02958 
02959 
02964     const Mapping<dim,spacedim> & get_mapping () const;
02965 
02971     const FiniteElement<dim,spacedim> & get_fe () const;
02972 
02977     UpdateFlags get_update_flags () const;
02978 
02983     const typename Triangulation<dim,spacedim>::cell_iterator get_cell () const;
02984 
02994     CellSimilarity::Similarity get_cell_similarity () const;
02995 
03001     std::size_t memory_consumption () const;
03003 
03004 
03015     DeclException0 (ExcAccessToUninitializedField);
03021     DeclException0 (ExcCannotInitializeField);
03027     DeclException0 (ExcInvalidUpdateFlag);
03033     DeclException0 (ExcFEDontMatch);
03039     DeclException1 (ExcShapeFunctionNotPrimitive,
03040                     int,
03041                     << "The shape function with index " << arg1
03042                     << " is not primitive, i.e. it is vector-valued and "
03043                     << "has more than one non-zero vector component. This "
03044                     << "function cannot be called for these shape functions. "
03045                     << "Maybe you want to use the same function with the "
03046                     << "_component suffix?");
03052     DeclException0 (ExcFENotPrimitive);
03053 
03054   protected:
03114     class CellIteratorBase;
03115 
03122     template <typename CI> class CellIterator;
03123     class TriaCellIterator;
03124 
03133     std::auto_ptr<const CellIteratorBase> present_cell;
03134 
03142     boost::signals2::connection tria_listener;
03143 
03150     void invalidate_present_cell ();
03151 
03161     void
03162     maybe_invalidate_previous_present_cell (const typename Triangulation<dim,spacedim>::cell_iterator &cell);
03163 
03167     const SmartPointer<const Mapping<dim,spacedim>,FEValuesBase<dim,spacedim> > mapping;
03168 
03172     const SmartPointer<const FiniteElement<dim,spacedim>,FEValuesBase<dim,spacedim> > fe;
03173 
03174 
03178     SmartPointer<typename Mapping<dim,spacedim>::InternalDataBase,FEValuesBase<dim,spacedim> > mapping_data;
03179 
03183     SmartPointer<typename Mapping<dim,spacedim>::InternalDataBase,FEValuesBase<dim,spacedim> > fe_data;
03184 
03202     UpdateFlags compute_update_flags (const UpdateFlags update_flags) const;
03203 
03213     CellSimilarity::Similarity cell_similarity;
03214 
03223     void
03224     check_cell_similarity (const typename Triangulation<dim,spacedim>::cell_iterator &cell);
03225 
03226   private:
03233     FEValuesBase (const FEValuesBase &);
03234 
03241     FEValuesBase & operator= (const FEValuesBase &);
03242 
03247     internal::FEValuesViews::Cache<dim,spacedim> fe_values_views_cache;
03248 
03254     template <int, int> friend class FEValuesViews::Scalar;
03255     template <int, int> friend class FEValuesViews::Vector;
03256     template <int, int, int> friend class FEValuesViews::SymmetricTensor;
03257 };
03258 
03259 
03260 
03271 template <int dim, int spacedim=dim>
03272 class FEValues : public FEValuesBase<dim,spacedim>
03273 {
03274   public:
03281     static const unsigned int integral_dimension = dim;
03282 
03290     FEValues (const Mapping<dim,spacedim>       &mapping,
03291               const FiniteElement<dim,spacedim> &fe,
03292               const Quadrature<dim>             &quadrature,
03293               const UpdateFlags                  update_flags);
03294 
03299     FEValues (const FiniteElement<dim,spacedim> &fe,
03300               const Quadrature<dim>             &quadrature,
03301               const UpdateFlags                  update_flags);
03302 
03316     void reinit (const typename DoFHandler<dim,spacedim>::cell_iterator &cell);
03317 
03331     void reinit (const typename hp::DoFHandler<dim,spacedim>::cell_iterator &cell);
03332 
03346     void reinit (const typename MGDoFHandler<dim,spacedim>::cell_iterator &cell);
03347 
03374     void reinit (const typename Triangulation<dim,spacedim>::cell_iterator &cell);
03375 
03381     const Quadrature<dim> & get_quadrature () const;
03382 
03388     std::size_t memory_consumption () const;
03389 
03413     const FEValues<dim,spacedim> & get_present_fe_values () const;
03414 
03415   private:
03420     const Quadrature<dim> quadrature;
03421 
03426     void initialize (const UpdateFlags update_flags);
03427 
03440     void do_reinit ();
03441 };
03442 
03443 
03454 template <int dim, int spacedim=dim>
03455 class FEFaceValuesBase : public FEValuesBase<dim,spacedim>
03456 {
03457   public:
03464     static const unsigned int integral_dimension = dim-1;
03465 
03484     FEFaceValuesBase (const unsigned int                 n_q_points,
03485                       const unsigned int                 dofs_per_cell,
03486                       const UpdateFlags                  update_flags,
03487                       const Mapping<dim,spacedim>       &mapping,
03488                       const FiniteElement<dim,spacedim> &fe,
03489                       const Quadrature<dim-1>&           quadrature);
03490 
03497     const Tensor<1,spacedim> & boundary_form (const unsigned int i) const;
03498 
03505     const std::vector<Tensor<1,spacedim> > & get_boundary_forms () const;
03506 
03512     unsigned int get_face_index() const;
03513 
03519     const Quadrature<dim-1> & get_quadrature () const;
03520 
03526     std::size_t memory_consumption () const;
03527 
03528   protected:
03529 
03535     unsigned int present_face_index;
03536 
03541     const Quadrature<dim-1> quadrature;
03542 };
03543 
03544 
03545 
03560 template <int dim, int spacedim=dim>
03561 class FEFaceValues : public FEFaceValuesBase<dim,spacedim>
03562 {
03563   public:
03569     static const unsigned int dimension = dim;
03570 
03571     static const unsigned int space_dimension = spacedim;
03572 
03579     static const unsigned int integral_dimension = dim-1;
03580 
03588     FEFaceValues (const Mapping<dim,spacedim>       &mapping,
03589                   const FiniteElement<dim,spacedim> &fe,
03590                   const Quadrature<dim-1>           &quadrature,
03591                   const UpdateFlags                  update_flags);
03592 
03597     FEFaceValues (const FiniteElement<dim,spacedim> &fe,
03598                   const Quadrature<dim-1>           &quadrature,
03599                   const UpdateFlags                  update_flags);
03600 
03607     void reinit (const typename DoFHandler<dim,spacedim>::cell_iterator &cell,
03608                  const unsigned int                                      face_no);
03609 
03623     void reinit (const typename hp::DoFHandler<dim,spacedim>::cell_iterator &cell,
03624                  const unsigned int                                          face_no);
03625 
03639     void reinit (const typename MGDoFHandler<dim,spacedim>::cell_iterator &cell,
03640                  const unsigned int                                        face_no);
03641 
03668     void reinit (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
03669                  const unsigned int                                         face_no);
03670 
03694     const FEFaceValues<dim,spacedim> & get_present_fe_values () const;
03695   private:
03696 
03701     void initialize (const UpdateFlags update_flags);
03702 
03715     void do_reinit (const unsigned int face_no);
03716 };
03717 
03718 
03737 template <int dim, int spacedim=dim>
03738 class FESubfaceValues : public FEFaceValuesBase<dim,spacedim>
03739 {
03740   public:
03745     static const unsigned int dimension = dim;
03746 
03751     static const unsigned int space_dimension = spacedim;
03752 
03759     static const unsigned int integral_dimension = dim-1;
03760 
03768     FESubfaceValues (const Mapping<dim,spacedim>       &mapping,
03769                      const FiniteElement<dim,spacedim> &fe,
03770                      const Quadrature<dim-1>  &face_quadrature,
03771                      const UpdateFlags         update_flags);
03772 
03777     FESubfaceValues (const FiniteElement<dim,spacedim> &fe,
03778                      const Quadrature<dim-1>  &face_quadrature,
03779                      const UpdateFlags         update_flags);
03780 
03794     void reinit (const typename DoFHandler<dim,spacedim>::cell_iterator &cell,
03795                  const unsigned int                    face_no,
03796                  const unsigned int                    subface_no);
03797 
03811     void reinit (const typename hp::DoFHandler<dim,spacedim>::cell_iterator &cell,
03812                  const unsigned int                    face_no,
03813                  const unsigned int                    subface_no);
03814 
03828     void reinit (const typename MGDoFHandler<dim,spacedim>::cell_iterator &cell,
03829                  const unsigned int                    face_no,
03830                  const unsigned int                    subface_no);
03831 
03858     void reinit (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
03859                  const unsigned int                    face_no,
03860                  const unsigned int                    subface_no);
03861 
03885     const FESubfaceValues<dim,spacedim> & get_present_fe_values () const;
03886 
03892     DeclException0 (ExcReinitCalledWithBoundaryFace);
03893 
03899     DeclException0 (ExcFaceHasNoSubfaces);
03900 
03901   private:
03902 
03907     void initialize (const UpdateFlags update_flags);
03908 
03921     void do_reinit (const unsigned int face_no,
03922                     const unsigned int subface_no);
03923 };
03924 
03925 
03926 #ifndef DOXYGEN
03927 
03928 
03929 /*------------------------ Inline functions: namespace FEValuesExtractors --------*/
03930 
03931 namespace FEValuesExtractors
03932 {
03933   inline
03934   Scalar::Scalar (const unsigned int component)
03935                   :
03936                   component (component)
03937   {}
03938 
03939 
03940 
03941   inline
03942   Vector::Vector (const unsigned int first_vector_component)
03943                   :
03944                   first_vector_component (first_vector_component)
03945   {}
03946 
03947   template <int rank>
03948   inline
03949   SymmetricTensor<rank>::SymmetricTensor (const unsigned int first_tensor_component)
03950                   :
03951                   first_tensor_component (first_tensor_component)
03952   {}
03953 }
03954 
03955 
03956 /*------------------------ Inline functions: namespace FEValuesViews --------*/
03957 
03958 namespace FEValuesViews
03959 {
03960   template <int dim, int spacedim>
03961   inline
03962   typename Scalar<dim,spacedim>::value_type
03963   Scalar<dim,spacedim>::value (const unsigned int shape_function,
03964                                const unsigned int q_point) const
03965   {
03966     typedef FEValuesBase<dim,spacedim> FVB;
03967     Assert (shape_function < fe_values.fe->dofs_per_cell,
03968             ExcIndexRange (shape_function, 0, fe_values.fe->dofs_per_cell));
03969     Assert (fe_values.update_flags & update_values,
03970             typename FVB::ExcAccessToUninitializedField());
03971 
03972                                      // an adaptation of the
03973                                      // FEValuesBase::shape_value_component
03974                                      // function except that here we know the
03975                                      // component as fixed and we have
03976                                      // pre-computed and cached a bunch of
03977                                      // information. see the comments there
03978     if (shape_function_data[shape_function].is_nonzero_shape_function_component)
03979       return fe_values.shape_values(shape_function_data[shape_function]
03980                                     .row_index,
03981                                     q_point);
03982     else
03983       return 0;
03984   }
03985 
03986 
03987 
03988 
03989   template <int dim, int spacedim>
03990   inline
03991   typename Scalar<dim,spacedim>::gradient_type
03992   Scalar<dim,spacedim>::gradient (const unsigned int shape_function,
03993                                   const unsigned int q_point) const
03994   {
03995     typedef FEValuesBase<dim,spacedim> FVB;
03996     Assert (shape_function < fe_values.fe->dofs_per_cell,
03997             ExcIndexRange (shape_function, 0, fe_values.fe->dofs_per_cell));
03998     Assert (fe_values.update_flags & update_gradients,
03999             typename FVB::ExcAccessToUninitializedField());
04000 
04001                                      // an adaptation of the
04002                                      // FEValuesBase::shape_grad_component
04003                                      // function except that here we know the
04004                                      // component as fixed and we have
04005                                      // pre-computed and cached a bunch of
04006                                      // information. see the comments there
04007     if (shape_function_data[shape_function].is_nonzero_shape_function_component)
04008       return fe_values.shape_gradients[shape_function_data[shape_function]
04009                                        .row_index][q_point];
04010     else
04011       return gradient_type();
04012   }
04013 
04014 
04015 
04016   template <int dim, int spacedim>
04017   inline
04018   typename Scalar<dim,spacedim>::hessian_type
04019   Scalar<dim,spacedim>::hessian (const unsigned int shape_function,
04020                                  const unsigned int q_point) const
04021   {
04022     typedef FEValuesBase<dim,spacedim> FVB;
04023     Assert (shape_function < fe_values.fe->dofs_per_cell,
04024             ExcIndexRange (shape_function, 0, fe_values.fe->dofs_per_cell));
04025     Assert (fe_values.update_flags & update_hessians,
04026             typename FVB::ExcAccessToUninitializedField());
04027 
04028                                      // an adaptation of the
04029                                      // FEValuesBase::shape_grad_component
04030                                      // function except that here we know the
04031                                      // component as fixed and we have
04032                                      // pre-computed and cached a bunch of
04033                                      // information. see the comments there
04034     if (shape_function_data[shape_function].is_nonzero_shape_function_component)
04035       return fe_values.shape_hessians[shape_function_data[shape_function].row_index][q_point];
04036     else
04037       return hessian_type();
04038   }
04039 
04040 
04041 
04042   template <int dim, int spacedim>
04043   inline
04044   typename Vector<dim,spacedim>::value_type
04045   Vector<dim,spacedim>::value (const unsigned int shape_function,
04046                                const unsigned int q_point) const
04047   {
04048     typedef FEValuesBase<dim,spacedim> FVB;
04049     Assert (shape_function < fe_values.fe->dofs_per_cell,
04050             ExcIndexRange (shape_function, 0, fe_values.fe->dofs_per_cell));
04051     Assert (fe_values.update_flags & update_values,
04052             typename FVB::ExcAccessToUninitializedField());
04053 
04054                                      // same as for the scalar case except
04055                                      // that we have one more index
04056     const int snc = shape_function_data[shape_function].single_nonzero_component;
04057     if (snc == -2)
04058       return value_type();
04059     else if (snc != -1)
04060       {
04061         value_type return_value;
04062         return_value[shape_function_data[shape_function].single_nonzero_component_index]
04063           = fe_values.shape_values(snc,q_point);
04064         return return_value;
04065       }
04066     else
04067       {
04068         value_type return_value;
04069         for (unsigned int d=0; d<dim; ++d)
04070           if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
04071             return_value[d]
04072               = fe_values.shape_values(shape_function_data[shape_function].row_index[d],q_point);
04073 
04074         return return_value;
04075       }
04076   }
04077 
04078 
04079 
04080   template <int dim, int spacedim>
04081   inline
04082   typename Vector<dim,spacedim>::gradient_type
04083   Vector<dim,spacedim>::gradient (const unsigned int shape_function,
04084                                   const unsigned int q_point) const
04085   {
04086     typedef FEValuesBase<dim,spacedim> FVB;
04087     Assert (shape_function < fe_values.fe->dofs_per_cell,
04088             ExcIndexRange (shape_function, 0, fe_values.fe->dofs_per_cell));
04089     Assert (fe_values.update_flags & update_gradients,
04090             typename FVB::ExcAccessToUninitializedField());
04091 
04092                                      // same as for the scalar case except
04093                                      // that we have one more index
04094     const int snc = shape_function_data[shape_function].single_nonzero_component;
04095     if (snc == -2)
04096       return gradient_type();
04097     else if (snc != -1)
04098       {
04099         gradient_type return_value;
04100         return_value[shape_function_data[shape_function].single_nonzero_component_index]
04101           = fe_values.shape_gradients[snc][q_point];
04102         return return_value;
04103       }
04104     else
04105       {
04106         gradient_type return_value;
04107         for (unsigned int d=0; d<dim; ++d)
04108           if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
04109             return_value[d]
04110               = fe_values.shape_gradients[shape_function_data[shape_function].row_index[d]][q_point];
04111 
04112         return return_value;
04113       }
04114   }
04115 
04116 
04117 
04118   template <int dim, int spacedim>
04119   inline
04120   typename Vector<dim,spacedim>::divergence_type
04121   Vector<dim,spacedim>::divergence (const unsigned int shape_function,
04122                                     const unsigned int q_point) const
04123   {
04124                                      // this function works like in
04125                                      // the case above
04126     typedef FEValuesBase<dim,spacedim> FVB;
04127     Assert (shape_function < fe_values.fe->dofs_per_cell,
04128             ExcIndexRange (shape_function, 0, fe_values.fe->dofs_per_cell));
04129     Assert (fe_values.update_flags & update_gradients,
04130             typename FVB::ExcAccessToUninitializedField());
04131 
04132                                      // same as for the scalar case except
04133                                      // that we have one more index
04134     const int snc = shape_function_data[shape_function].single_nonzero_component;
04135     if (snc == -2)
04136       return divergence_type();
04137     else if (snc != -1)
04138       return
04139         fe_values.shape_gradients[snc][q_point][shape_function_data[shape_function].single_nonzero_component_index];
04140     else
04141       {
04142         divergence_type return_value = 0;
04143         for (unsigned int d=0; d<dim; ++d)
04144           if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
04145             return_value
04146               += fe_values.shape_gradients[shape_function_data[shape_function].row_index[d]][q_point][d];
04147 
04148         return return_value;
04149       }
04150   }
04151 
04152 
04153 
04154   template <int dim, int spacedim>
04155   inline
04156   typename Vector<dim,spacedim>::curl_type
04157   Vector<dim,spacedim>::curl (const unsigned int shape_function, const unsigned int q_point) const {
04158      // this function works like in the case above
04159      typedef FEValuesBase<dim,spacedim> FVB;
04160 
04161      Assert (shape_function < fe_values.fe->dofs_per_cell,
04162          ExcIndexRange (shape_function, 0, fe_values.fe->dofs_per_cell));
04163      Assert (fe_values.update_flags & update_gradients,
04164          typename FVB::ExcAccessToUninitializedField());
04165      // same as for the scalar case except that we have one more index
04166      const int snc = shape_function_data[shape_function].single_nonzero_component;
04167 
04168      if (snc == -2)
04169         return curl_type ();
04170 
04171      else
04172         switch (dim) {
04173            case 1: {
04174                   Assert (false, ExcMessage("Computing the curl in 1d is not a useful operation"));
04175               return curl_type ();
04176            }
04177 
04178            case 2: {
04179               if (snc != -1) {
04180                  curl_type return_value;
04181 
04182                                                   // the single
04183                                                   // nonzero component
04184                                                   // can only be zero
04185                                                   // or one in 2d
04186                  if (shape_function_data[shape_function].single_nonzero_component_index == 0)
04187                    return_value[0] = -1.0 * fe_values.shape_gradients[snc][q_point][1];
04188                  else
04189                    return_value[0] = fe_values.shape_gradients[snc][q_point][0];
04190 
04191                  return return_value;
04192               }
04193 
04194               else {
04195                  curl_type return_value;
04196 
04197                  return_value[0] = 0.0;
04198 
04199                  if (shape_function_data[shape_function].is_nonzero_shape_function_component[0])
04200                     return_value[0]
04201                       -= fe_values.shape_gradients[shape_function_data[shape_function].row_index[0]][q_point][1];
04202 
04203                  if (shape_function_data[shape_function].is_nonzero_shape_function_component[1])
04204                     return_value[0]
04205                       += fe_values.shape_gradients[shape_function_data[shape_function].row_index[1]][q_point][0];
04206 
04207                  return return_value;
04208               }
04209            }
04210 
04211            case 3: {
04212               if (snc != -1) {
04213                  curl_type return_value;
04214 
04215                  switch (shape_function_data[shape_function].single_nonzero_component_index) {
04216                     case 0: {
04217                        return_value[0] = 0;
04218                        return_value[1] = fe_values.shape_gradients[snc][q_point][2];
04219                        return_value[2] = -1.0 * fe_values.shape_gradients[snc][q_point][1];
04220                        return return_value;
04221                     }
04222 
04223                     case 1: {
04224                        return_value[0] = -1.0 * fe_values.shape_gradients[snc][q_point][2];
04225                        return_value[1] = 0;
04226                        return_value[2] = fe_values.shape_gradients[snc][q_point][0];
04227                        return return_value;
04228                     }
04229 
04230                     default: {
04231                        return_value[0] = fe_values.shape_gradients[snc][q_point][1];
04232                        return_value[1] = -1.0 * fe_values.shape_gradients[snc][q_point][0];
04233                        return_value[2] = 0;
04234                        return return_value;
04235                     }
04236                  }
04237               }
04238 
04239               else {
04240                  curl_type return_value;
04241 
04242                  for (unsigned int i = 0; i < dim; ++i)
04243                     return_value[i] = 0.0;
04244 
04245                  if (shape_function_data[shape_function].is_nonzero_shape_function_component[0]) {
04246                     return_value[1]
04247                       += fe_values.shape_gradients[shape_function_data[shape_function].row_index[0]][q_point][2];
04248                     return_value[2]
04249                       -= fe_values.shape_gradients[shape_function_data[shape_function].row_index[0]][q_point][1];
04250                  }
04251 
04252                  if (shape_function_data[shape_function].is_nonzero_shape_function_component[1]) {
04253                     return_value[0]
04254                       -= fe_values.shape_gradients[shape_function_data[shape_function].row_index[1]][q_point][2];
04255                     return_value[2]
04256                       += fe_values.shape_gradients[shape_function_data[shape_function].row_index[1]][q_point][0];
04257                  }
04258 
04259                  if (shape_function_data[shape_function].is_nonzero_shape_function_component[2]) {
04260                     return_value[0]
04261                       += fe_values.shape_gradients[shape_function_data[shape_function].row_index[2]][q_point][1];
04262                     return_value[1]
04263                       -= fe_values.shape_gradients[shape_function_data[shape_function].row_index[2]][q_point][0];
04264                  }
04265 
04266                  return return_value;
04267               }
04268            }
04269         }
04270                                 // should not end up here
04271      Assert (false, ExcInternalError());
04272      return curl_type();
04273   }
04274 
04275   template <int dim, int spacedim>
04276   inline
04277   typename Vector<dim,spacedim>::hessian_type
04278   Vector<dim,spacedim>::hessian (const unsigned int shape_function,
04279                                  const unsigned int q_point) const
04280   {
04281                                      // this function works like in
04282                                      // the case above
04283     typedef FEValuesBase<dim,spacedim> FVB;
04284     Assert (shape_function < fe_values.fe->dofs_per_cell,
04285             ExcIndexRange (shape_function, 0, fe_values.fe->dofs_per_cell));
04286     Assert (fe_values.update_flags & update_hessians,
04287             typename FVB::ExcAccessToUninitializedField());
04288 
04289                                      // same as for the scalar case except
04290                                      // that we have one more index
04291     const int snc = shape_function_data[shape_function].single_nonzero_component;
04292     if (snc == -2)
04293       return hessian_type();
04294     else if (snc != -1)
04295       {
04296         hessian_type return_value;
04297         return_value[shape_function_data[shape_function].single_nonzero_component_index]
04298           = fe_values.shape_hessians[snc][q_point];
04299         return return_value;
04300       }
04301     else
04302       {
04303         hessian_type return_value;
04304         for (unsigned int d=0; d<dim; ++d)
04305           if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
04306             return_value[d]
04307               = fe_values.shape_hessians[shape_function_data[shape_function].row_index[d]][q_point];
04308 
04309         return return_value;
04310       }
04311   }
04312 
04313 
04314   namespace
04315   {
04322     inline
04323     ::SymmetricTensor<2,1>
04324     symmetrize_single_row (const unsigned int n,
04325                            const Tensor<1,1> &t)
04326     {
04327       Assert (n < 1, ExcIndexRange (n, 0, 1));
04328 
04329       const double array[1] = { t[0] };
04330       return ::SymmetricTensor<2,1>(array);
04331     }
04332 
04333 
04334     inline
04335     ::SymmetricTensor<2,2>
04336     symmetrize_single_row (const unsigned int n,
04337                            const Tensor<1,2> &t)
04338     {
04339       switch (n)
04340         {
04341           case 0:
04342           {
04343             const double array[3] = { t[0], 0, t[1]/2 };
04344             return ::SymmetricTensor<2,2>(array);
04345           }
04346           case 1:
04347           {
04348             const double array[3] = { 0, t[1], t[0]/2 };
04349             return ::SymmetricTensor<2,2>(array);
04350           }
04351           default:
04352           {
04353             Assert (false, ExcIndexRange (n, 0, 2));
04354             return ::SymmetricTensor<2,2>();
04355           }
04356         }
04357     }
04358 
04359 
04360     inline
04361     ::SymmetricTensor<2,3>
04362     symmetrize_single_row (const unsigned int n,
04363                            const Tensor<1,3> &t)
04364     {
04365       switch (n)
04366         {
04367           case 0:
04368           {
04369             const double array[6] = { t[0], 0, 0, t[1]/2, t[2]/2, 0 };
04370             return ::SymmetricTensor<2,3>(array);
04371           }
04372           case 1:
04373           {
04374             const double array[6] = { 0, t[1], 0, t[0]/2, 0, t[2]/2 };
04375             return ::SymmetricTensor<2,3>(array);
04376           }
04377           case 2:
04378           {
04379             const double array[6] = { 0, 0, t[2], 0, t[0]/2, t[1]/2 };
04380             return ::SymmetricTensor<2,3>(array);
04381           }
04382           default:
04383           {
04384             Assert (false, ExcIndexRange (n, 0, 3));
04385             return ::SymmetricTensor<2,3>();
04386           }
04387         }
04388     }
04389   }
04390 
04391 
04392   template <int dim, int spacedim>
04393   inline
04394   typename Vector<dim,spacedim>::symmetric_gradient_type
04395   Vector<dim,spacedim>::symmetric_gradient (const unsigned int shape_function,
04396                                             const unsigned int q_point) const
04397   {
04398     typedef FEValuesBase<dim,spacedim> FVB;
04399     Assert (shape_function < fe_values.fe->dofs_per_cell,
04400             ExcIndexRange (shape_function, 0, fe_values.fe->dofs_per_cell));
04401     Assert (fe_values.update_flags & update_gradients,
04402             typename FVB::ExcAccessToUninitializedField());
04403 
04404                                      // same as for the scalar case except
04405                                      // that we have one more index
04406     const int snc = shape_function_data[shape_function].single_nonzero_component;
04407     if (snc == -2)
04408       return symmetric_gradient_type();
04409     else if (snc != -1)
04410       return symmetrize_single_row (shape_function_data[shape_function].single_nonzero_component_index,
04411                                     fe_values.shape_gradients[snc][q_point]);
04412     else
04413       {
04414         gradient_type return_value;
04415         for (unsigned int d=0; d<dim; ++d)
04416           if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
04417             return_value[d]
04418               = fe_values.shape_gradients[shape_function_data[shape_function].row_index[d]][q_point];
04419 
04420         return symmetrize(return_value);
04421       }
04422   }
04423 
04424 
04425 
04426   template <int dim, int spacedim>
04427   inline
04428   typename SymmetricTensor<2, dim, spacedim>::value_type
04429   SymmetricTensor<2, dim, spacedim>::value (const unsigned int shape_function,
04430                                             const unsigned int q_point) const
04431   {
04432     typedef FEValuesBase<dim,spacedim> FVB;
04433     Assert (shape_function < fe_values.fe->dofs_per_cell,
04434             ExcIndexRange (shape_function, 0, fe_values.fe->dofs_per_cell));
04435     Assert (fe_values.update_flags & update_values,
04436             typename FVB::ExcAccessToUninitializedField());
04437 
04438                                      // similar to the vector case where we
04439                                      // have more then one index and we need
04440                                      // to convert between unrolled and
04441                                      // component indexing for tensors
04442     const int snc
04443       = shape_function_data[shape_function].single_nonzero_component;
04444 
04445     if (snc == -2)
04446       {
04447                                          // shape function is zero for the
04448                                          // selected components
04449         return value_type();
04450 
04451       }
04452     else if (snc != -1)
04453       {
04454         value_type return_value;
04455         const unsigned int comp =
04456           shape_function_data[shape_function].single_nonzero_component_index;
04457         return_value[value_type::unrolled_to_component_indices(comp)]
04458           = fe_values.shape_values(snc,q_point);
04459         return return_value;
04460       }
04461     else
04462       {
04463         value_type return_value;
04464         for (unsigned int d = 0; d < value_type::n_independent_components; ++d)
04465           if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
04466             return_value[value_type::unrolled_to_component_indices(d)]
04467               = fe_values.shape_values(shape_function_data[shape_function].row_index[d],q_point);
04468         return return_value;
04469       }
04470   }
04471 
04472 
04473   template <int dim, int spacedim>
04474   inline
04475   typename SymmetricTensor<2, dim, spacedim>::divergence_type
04476   SymmetricTensor<2, dim, spacedim>::divergence(const unsigned int shape_function,
04477                                                 const unsigned int q_point) const
04478   {
04479     typedef FEValuesBase<dim,spacedim> FVB;
04480     Assert (shape_function < fe_values.fe->dofs_per_cell,
04481             ExcIndexRange (shape_function, 0, fe_values.fe->dofs_per_cell));
04482     Assert (fe_values.update_flags & update_gradients,
04483             typename FVB::ExcAccessToUninitializedField());
04484 
04485     const int snc = shape_function_data[shape_function].single_nonzero_component;
04486 
04487     if (snc == -2)
04488       {
04489                                          // shape function is zero for the
04490                                          // selected components
04491         return divergence_type();
04492       }
04493     else if (snc != -1)
04494       {
04495                                          // we have a single non-zero component
04496                                          // when the symmetric tensor is
04497                                          // represented in unrolled form.
04498                                          // this implies we potentially have
04499                                          // two non-zero components when
04500                                          // represented in component form!  we
04501                                          // will only have one non-zero entry
04502                                          // if the non-zero component lies on
04503                                          // the diagonal of the tensor.
04504                                          //
04505                                          // the divergence of a second-order tensor
04506                                          // is a first order tensor.
04507                                          //
04508                                          // assume the second-order tensor is
04509                                          // A with components A_{ij}.  then
04510                                          // A_{ij} = A_{ji} and there is only
04511                                          // one (if diagonal) or two non-zero
04512                                          // entries in the tensorial
04513                                          // representation.  define the
04514                                          // divergence as:
04515                                          // b_i := \dfrac{\partial phi_{ij}}{\partial x_j}.
04516                                          // (which is incidentally also
04517                                          // b_j := \dfrac{\partial phi_{ij}}{\partial x_i}).
04518                                          // In both cases, a sum is implied.
04519                                          //
04520                                          // Now, we know the nonzero component
04521                                          // in unrolled form: it is indicated
04522                                          // by 'snc'. we can figure out which
04523                                          // tensor components belong to this:
04524         const unsigned int comp =
04525           shape_function_data[shape_function].single_nonzero_component_index;
04526         const unsigned int ii = value_type::unrolled_to_component_indices(comp)[0];
04527         const unsigned int jj = value_type::unrolled_to_component_indices(comp)[1];
04528 
04529                                          // given the form of the divergence
04530                                          // above, if ii=jj there is only a
04531                                          // single nonzero component of the
04532                                          // full tensor and the gradient
04533                                          // equals
04534                                          // b_ii := \dfrac{\partial phi_{ii,ii}}{\partial x_ii}.
04535                                          // all other entries of 'b' are zero
04536                                          //
04537                                          // on the other hand, if ii!=jj, then
04538                                          // there are two nonzero entries in
04539                                          // the full tensor and
04540                                          // b_ii := \dfrac{\partial phi_{ii,jj}}{\partial x_ii}.
04541                                          // b_jj := \dfrac{\partial phi_{ii,jj}}{\partial x_jj}.
04542                                          // again, all other entries of 'b' are
04543                                          // zero
04544         const Tensor<1, spacedim> phi_grad = fe_values.shape_gradients[snc][q_point];
04545 
04546         divergence_type return_value;
04547         return_value[ii] = phi_grad[jj];
04548 
04549         if (ii != jj)
04550           return_value[jj] = phi_grad[ii];
04551 
04552         return return_value;
04553 
04554       }
04555     else
04556       {
04557         Assert (false, ExcNotImplemented());
04558         divergence_type return_value;
04559         return return_value;
04560       }
04561   }
04562 
04563 }
04564 
04565 
04566 
04567 /*------------------------ Inline functions: FEValuesBase ------------------------*/
04568 
04569 
04570 
04571 template <int dim, int spacedim>
04572 inline
04573 const FEValuesViews::Scalar<dim,spacedim> &
04574 FEValuesBase<dim,spacedim>::
04575 operator[] (const FEValuesExtractors::Scalar &scalar) const
04576 {
04577   Assert (scalar.component < fe_values_views_cache.scalars.size(),
04578           ExcIndexRange (scalar.component,
04579                          0, fe_values_views_cache.scalars.size()));
04580 
04581   return fe_values_views_cache.scalars[scalar.component];
04582 }
04583 
04584 
04585 
04586 template <int dim, int spacedim>
04587 inline
04588 const FEValuesViews::Vector<dim,spacedim> &
04589 FEValuesBase<dim,spacedim>::
04590 operator[] (const FEValuesExtractors::Vector &vector) const
04591 {
04592   Assert (vector.first_vector_component <
04593           fe_values_views_cache.vectors.size(),
04594           ExcIndexRange (vector.first_vector_component,
04595                          0, fe_values_views_cache.vectors.size()));
04596 
04597   return fe_values_views_cache.vectors[vector.first_vector_component];
04598 }
04599 
04600 template <int dim, int spacedim>
04601 inline
04602 const FEValuesViews::SymmetricTensor<2,dim,spacedim> &
04603 FEValuesBase<dim,spacedim>::
04604 operator[] (const FEValuesExtractors::SymmetricTensor<2> &tensor) const
04605 {
04606   Assert (tensor.first_tensor_component <
04607           fe_values_views_cache.symmetric_second_order_tensors.size(),
04608           ExcIndexRange (tensor.first_tensor_component,
04609                          0, fe_values_views_cache.symmetric_second_order_tensors.size()));
04610 
04611   return fe_values_views_cache.symmetric_second_order_tensors[tensor.first_tensor_component];
04612 }
04613 
04614 template <int dim, int spacedim>
04615 inline
04616 const double &
04617 FEValuesBase<dim,spacedim>::shape_value (const unsigned int i,
04618                                          const unsigned int j) const
04619 {
04620   Assert (i < fe->dofs_per_cell,
04621           ExcIndexRange (i, 0, fe->dofs_per_cell));
04622   Assert (this->update_flags & update_values,
04623           ExcAccessToUninitializedField());
04624   Assert (fe->is_primitive (i),
04625           ExcShapeFunctionNotPrimitive(i));
04626 
04627                                    // if the entire FE is primitive,
04628                                    // then we can take a short-cut:
04629   if (fe->is_primitive())
04630     return this->shape_values(i,j);
04631   else
04632                                      // otherwise, use the mapping
04633                                      // between shape function numbers
04634                                      // and rows. note that by the
04635                                      // assertions above, we know that
04636                                      // this particular shape function
04637                                      // is primitive, so there is no
04638                                      // question to which vector
04639                                      // component the call of this
04640                                      // function refers
04641     return this->shape_values(this->shape_function_to_row_table[i], j);
04642 }
04643 
04644 
04645 
04646 template <int dim, int spacedim>
04647 inline
04648 double
04649 FEValuesBase<dim,spacedim>::shape_value_component (const unsigned int i,
04650                                                    const unsigned int j,
04651                                                    const unsigned int component) const
04652 {
04653   Assert (i < fe->dofs_per_cell,
04654           ExcIndexRange (i, 0, fe->dofs_per_cell));
04655   Assert (this->update_flags & update_values,
04656           ExcAccessToUninitializedField());
04657   Assert (component < fe->n_components(),
04658           ExcIndexRange(component, 0, fe->n_components()));
04659 
04660                                    // if this particular shape
04661                                    // function is primitive, then we
04662                                    // can take a short-cut by checking
04663                                    // whether the requested component
04664                                    // is the only non-zero one (note
04665                                    // that calling
04666                                    // system_to_component_table only
04667                                    // works if the shape function is
04668                                    // primitive):
04669   if (fe->is_primitive(i))
04670     {
04671       if (component == fe->system_to_component_index(i).first)
04672         return this->shape_values(this->shape_function_to_row_table[i],j);
04673       else
04674         return 0;
04675     }
04676   else
04677     {
04678                                        // no, this shape function is
04679                                        // not primitive. then we have
04680                                        // to loop over its components
04681                                        // to find the corresponding
04682                                        // row in the arrays of this
04683                                        // object. before that check
04684                                        // whether the shape function
04685                                        // is non-zero at all within
04686                                        // this component:
04687       if (fe->get_nonzero_components(i)[component] == false)
04688         return 0;
04689 
04690                                        // count how many non-zero
04691                                        // component the shape function
04692                                        // has before the one we are
04693                                        // looking for, and add this to
04694                                        // the offset of the first
04695                                        // non-zero component of this
04696                                        // shape function in the arrays
04697                                        // we index presently:
04698       const unsigned int
04699         row = (this->shape_function_to_row_table[i]
04700                +
04701                std::count (fe->get_nonzero_components(i).begin(),
04702                            fe->get_nonzero_components(i).begin()+component,
04703                            true));
04704       return this->shape_values(row, j);
04705     }
04706 }
04707 
04708 
04709 
04710 template <int dim, int spacedim>
04711 inline
04712 const Tensor<1,spacedim> &
04713 FEValuesBase<dim,spacedim>::shape_grad (const unsigned int i,
04714                                         const unsigned int j) const
04715 {
04716   Assert (i < fe->dofs_per_cell,
04717           ExcIndexRange (i, 0, fe->dofs_per_cell));
04718   Assert (this->update_flags & update_gradients,
04719           ExcAccessToUninitializedField());
04720   Assert (fe->is_primitive (i),
04721           ExcShapeFunctionNotPrimitive(i));
04722   Assert (i<this->shape_gradients.size(),
04723           ExcIndexRange (i, 0, this->shape_gradients.size()));
04724   Assert (j<this->shape_gradients[0].size(),
04725           ExcIndexRange (j, 0, this->shape_gradients[0].size()));
04726 
04727                                    // if the entire FE is primitive,
04728                                    // then we can take a short-cut:
04729   if (fe->is_primitive())
04730     return this->shape_gradients[i][j];
04731   else
04732                                      // otherwise, use the mapping
04733                                      // between shape function numbers
04734                                      // and rows. note that by the
04735                                      // assertions above, we know that
04736                                      // this particular shape function
04737                                      // is primitive, so there is no
04738                                      // question to which vector
04739                                      // component the call of this
04740                                      // function refers
04741     return this->shape_gradients[this->shape_function_to_row_table[i]][j];
04742 }
04743 
04744 
04745 
04746 template <int dim, int spacedim>
04747 inline
04748 Tensor<1,spacedim>
04749 FEValuesBase<dim,spacedim>::shape_grad_component (const unsigned int i,
04750                                                   const unsigned int j,
04751                                                   const unsigned int component) const
04752 {
04753   Assert (i < fe->dofs_per_cell,
04754           ExcIndexRange (i, 0, fe->dofs_per_cell));
04755   Assert (this->update_flags & update_gradients,
04756           ExcAccessToUninitializedField());
04757   Assert (component < fe->n_components(),
04758           ExcIndexRange(component, 0, fe->n_components()));
04759 
04760                                    // if this particular shape
04761                                    // function is primitive, then we
04762                                    // can take a short-cut by checking
04763                                    // whether the requested component
04764                                    // is the only non-zero one (note
04765                                    // that calling
04766                                    // system_to_component_table only
04767                                    // works if the shape function is
04768                                    // primitive):
04769   if (fe->is_primitive(i))
04770     {
04771       if (component == fe->system_to_component_index(i).first)
04772         return this->shape_gradients[this->shape_function_to_row_table[i]][j];
04773       else
04774         return Tensor<1,spacedim>();
04775     }
04776   else
04777     {
04778                                        // no, this shape function is
04779                                        // not primitive. then we have
04780                                        // to loop over its components
04781                                        // to find the corresponding
04782                                        // row in the arrays of this
04783                                        // object. before that check
04784                                        // whether the shape function
04785                                        // is non-zero at all within
04786                                        // this component:
04787       if (fe->get_nonzero_components(i)[component] == false)
04788         return Tensor<1,spacedim>();
04789 
04790                                        // count how many non-zero
04791                                        // component the shape function
04792                                        // has before the one we are
04793                                        // looking for, and add this to
04794                                        // the offset of the first
04795                                        // non-zero component of this
04796                                        // shape function in the arrays
04797                                        // we index presently:
04798       const unsigned int
04799         row = (this->shape_function_to_row_table[i]
04800                +
04801                std::count (fe->get_nonzero_components(i).begin(),
04802                            fe->get_nonzero_components(i).begin()+component,
04803                            true));
04804       return this->shape_gradients[row][j];
04805     }
04806 }
04807 
04808 
04809 
04810 template <int dim, int spacedim>
04811 inline
04812 const Tensor<2,spacedim> &
04813 FEValuesBase<dim,spacedim>::shape_hessian (const unsigned int i,
04814                                            const unsigned int j) const
04815 {
04816   Assert (i < fe->dofs_per_cell,
04817           ExcIndexRange (i, 0, fe->dofs_per_cell));
04818   Assert (this->update_flags & update_hessians,
04819           ExcAccessToUninitializedField());
04820   Assert (fe->is_primitive (i),
04821           ExcShapeFunctionNotPrimitive(i));
04822   Assert (i<this->shape_hessians.size(),
04823           ExcIndexRange (i, 0, this->shape_hessians.size()));
04824   Assert (j<this->shape_hessians[0].size(),
04825           ExcIndexRange (j, 0, this->shape_hessians[0].size()));
04826 
04827                                    // if the entire FE is primitive,
04828                                    // then we can take a short-cut:
04829   if (fe->is_primitive())
04830     return this->shape_hessians[i][j];
04831   else
04832                                      // otherwise, use the mapping
04833                                      // between shape function numbers
04834                                      // and rows. note that by the
04835                                      // assertions above, we know that
04836                                      // this particular shape function
04837                                      // is primitive, so there is no
04838                                      // question to which vector
04839                                      // component the call of this
04840                                      // function refers
04841     return this->shape_hessians[this->shape_function_to_row_table[i]][j];
04842 }
04843 
04844 
04845 
04846 template <int dim, int spacedim>
04847 inline
04848 const Tensor<2,spacedim> &
04849 FEValuesBase<dim,spacedim>::shape_2nd_derivative (const unsigned int i,
04850                                                   const unsigned int j) const
04851 {
04852   return shape_hessian(i,j);
04853 }
04854 
04855 
04856 
04857 template <int dim, int spacedim>
04858 inline
04859 Tensor<2,spacedim>
04860 FEValuesBase<dim,spacedim>::shape_hessian_component (const unsigned int i,
04861                                                      const unsigned int j,
04862                                                      const unsigned int component) const
04863 {
04864   Assert (i < fe->dofs_per_cell,
04865           ExcIndexRange (i, 0, fe->dofs_per_cell));
04866   Assert (this->update_flags & update_hessians,
04867           ExcAccessToUninitializedField());
04868   Assert (component < fe->n_components(),
04869           ExcIndexRange(component, 0, fe->n_components()));
04870 
04871                                    // if this particular shape
04872                                    // function is primitive, then we
04873                                    // can take a short-cut by checking
04874                                    // whether the requested component
04875                                    // is the only non-zero one (note
04876                                    // that calling
04877                                    // system_to_component_table only
04878                                    // works if the shape function is
04879                                    // primitive):
04880   if (fe->is_primitive(i))
04881     {
04882       if (component == fe->system_to_component_index(i).first)
04883         return this->shape_hessians[this->shape_function_to_row_table[i]][j];
04884       else
04885         return Tensor<2,spacedim>();
04886     }
04887   else
04888     {
04889                                        // no, this shape function is
04890                                        // not primitive. then we have
04891                                        // to loop over its components
04892                                        // to find the corresponding
04893                                        // row in the arrays of this
04894                                        // object. before that check
04895                                        // whether the shape function
04896                                        // is non-zero at all within
04897                                        // this component:
04898       if (fe->get_nonzero_components(i)[component] == false)
04899         return Tensor<2,spacedim>();
04900 
04901                                        // count how many non-zero
04902                                        // component the shape function
04903                                        // has before the one we are
04904                                        // looking for, and add this to
04905                                        // the offset of the first
04906                                        // non-zero component of this
04907                                        // shape function in the arrays
04908                                        // we index presently:
04909       const unsigned int
04910         row = (this->shape_function_to_row_table[i]
04911                +
04912                std::count (fe->get_nonzero_components(i).begin(),
04913                            fe->get_nonzero_components(i).begin()+component,
04914                            true));
04915       return this->shape_hessians[row][j];
04916     }
04917 }
04918 
04919 
04920 
04921 template <int dim, int spacedim>
04922 inline
04923 Tensor<2,spacedim>
04924 FEValuesBase<dim,spacedim>::shape_2nd_derivative_component (const unsigned int i,
04925                                                             const unsigned int j,
04926                                                             const unsigned int component) const
04927 {
04928   return shape_hessian_component(i,j,component);
04929 }
04930 
04931 
04932 
04933 template <int dim, int spacedim>
04934 inline
04935 const FiniteElement<dim,spacedim> &
04936 FEValuesBase<dim,spacedim>::get_fe () const
04937 {
04938   return *fe;
04939 }
04940 
04941 
04942 template <int dim, int spacedim>
04943 inline
04944 const Mapping<dim,spacedim> &
04945 FEValuesBase<dim,spacedim>::get_mapping () const
04946 {
04947   return *mapping;
04948 }
04949 
04950 
04951 
04952 template <int dim, int spacedim>
04953 inline
04954 UpdateFlags
04955 FEValuesBase<dim,spacedim>::get_update_flags () const
04956 {
04957   return this->update_flags;
04958 }
04959 
04960 
04961 
04962 template <int dim, int spacedim>
04963 inline
04964 const std::vector<Point<spacedim> > &
04965 FEValuesBase<dim,spacedim>::get_quadrature_points () const
04966 {
04967   Assert (this->update_flags & update_quadrature_points, ExcAccessToUninitializedField());
04968   return this->quadrature_points;
04969 }
04970 
04971 
04972 
04973 template <int dim, int spacedim>
04974 inline
04975 const std::vector<double> &
04976 FEValuesBase<dim,spacedim>::get_JxW_values () const
04977 {
04978   Assert (this->update_flags & update_JxW_values, ExcAccessToUninitializedField());
04979   return this->JxW_values;
04980 }
04981 
04982 
04983 
04984 template <int dim, int spacedim>
04985 inline
04986 const std::vector<DerivativeForm<1,dim,spacedim> > &
04987 FEValuesBase<dim,spacedim>::get_jacobians () const
04988 {
04989   Assert (this->update_flags & update_jacobians, ExcAccessToUninitializedField());
04990   return this->jacobians;
04991 }
04992 
04993 
04994 
04995 template <int dim, int spacedim>
04996 inline
04997 const std::vector<DerivativeForm<2,dim,spacedim> > &
04998 FEValuesBase<dim,spacedim>::get_jacobian_grads () const
04999 {
05000   Assert (this->update_flags & update_jacobian_grads, ExcAccessToUninitializedField());
05001   return this->jacobian_grads;
05002 }
05003 
05004 
05005 
05006 template <int dim, int spacedim>
05007 inline
05008 const std::vector<DerivativeForm<1,spacedim,dim> > &
05009 FEValuesBase<dim,spacedim>::get_inverse_jacobians () const
05010 {
05011   Assert (this->update_flags & update_inverse_jacobians, ExcAccessToUninitializedField());
05012   return this->inverse_jacobians;
05013 }
05014 
05015 
05016 
05017 template <int dim, int spacedim>
05018 inline
05019 const Point<spacedim> &
05020 FEValuesBase<dim,spacedim>::quadrature_point (const unsigned int i) const
05021 {
05022   Assert (this->update_flags & update_quadrature_points, ExcAccessToUninitializedField());
05023   Assert (i<this->quadrature_points.size(), ExcIndexRange(i, 0, this->quadrature_points.size()));
05024 
05025   return this->quadrature_points[i];
05026 }
05027 
05028 
05029 
05030 
05031 template <int dim, int spacedim>
05032 inline
05033 double
05034 FEValuesBase<dim,spacedim>::JxW (const unsigned int i) const
05035 {
05036   Assert (this->update_flags & update_JxW_values, ExcAccessToUninitializedField());
05037   Assert (i<this->JxW_values.size(), ExcIndexRange(i, 0, this->JxW_values.size()));
05038 
05039   return this->JxW_values[i];
05040 }
05041 
05042 
05043 
05044 template <int dim, int spacedim>
05045 inline
05046 const DerivativeForm<1,dim,spacedim> &
05047 FEValuesBase<dim,spacedim>::jacobian (const unsigned int i) const
05048 {
05049   Assert (this->update_flags & update_jacobians, ExcAccessToUninitializedField());
05050   Assert (i<this->jacobians.size(), ExcIndexRange(i, 0, this->jacobians.size()));
05051 
05052   return this->jacobians[i];
05053 }
05054 
05055 
05056 
05057 template <int dim, int spacedim>
05058 inline
05059 const DerivativeForm<2,dim,spacedim> &
05060 FEValuesBase<dim,spacedim>::jacobian_grad (const unsigned int i) const
05061 {
05062   Assert (this->update_flags & update_jacobian_grads, ExcAccessToUninitializedField());
05063   Assert (i<this->jacobian_grads.size(), ExcIndexRange(i, 0, this->jacobian_grads.size()));
05064 
05065   return this->jacobian_grads[i];
05066 }
05067 
05068 
05069 
05070 template <int dim, int spacedim>
05071 inline
05072 const DerivativeForm<1,spacedim,dim> &
05073 FEValuesBase<dim,spacedim>::inverse_jacobian (const unsigned int i) const
05074 {
05075   Assert (this->update_flags & update_inverse_jacobians, ExcAccessToUninitializedField());
05076   Assert (i<this->inverse_jacobians.size(), ExcIndexRange(i, 0, this->inverse_jacobians.size()));
05077 
05078   return this->inverse_jacobians[i];
05079 }
05080 
05081 
05082 template <int dim, int spacedim>
05083 template <class InputVector>
05084 inline
05085 void
05086 FEValuesBase<dim,spacedim>::get_function_grads (const InputVector           &fe_function,
05087                                                 std::vector<Tensor<1,spacedim> > &gradients) const
05088 {
05089   get_function_gradients(fe_function, gradients);
05090 }
05091 
05092 
05093 
05094 template <int dim, int spacedim>
05095 template <class InputVector>
05096 inline
05097 void
05098 FEValuesBase<dim,spacedim>::get_function_grads (
05099   const InputVector& fe_function,
05100   const VectorSlice<const std::vector<unsigned int> >& indices,
05101   std::vector<Tensor<1,spacedim> > &values) const
05102 {
05103   get_function_gradients(fe_function, indices, values);
05104 }
05105 
05106 
05107 
05108 template <int dim, int spacedim>
05109 template <class InputVector>
05110 inline
05111 void
05112 FEValuesBase<dim,spacedim>::
05113 get_function_grads (const InputVector                         &fe_function,
05114                     std::vector<std::vector<Tensor<1,spacedim> > > &gradients) const
05115 {
05116   get_function_gradients(fe_function, gradients);
05117 }
05118 
05119 
05120 
05121 template <int dim, int spacedim>
05122 template <class InputVector>
05123 inline
05124 void
05125 FEValuesBase<dim,spacedim>::get_function_grads (
05126   const InputVector& fe_function,
05127   const VectorSlice<const std::vector<unsigned int> >& indices,
05128   std::vector<std::vector<Tensor<1,spacedim> > >& values,
05129   bool q_points_fastest) const
05130 {
05131   get_function_gradients(fe_function, indices, values, q_points_fastest);
05132 }
05133 
05134 
05135 
05136 template <int dim, int spacedim>
05137 template <class InputVector>
05138 inline
05139 void
05140 FEValuesBase<dim,spacedim>::
05141 get_function_2nd_derivatives (const InputVector           &fe_function,
05142                               std::vector<Tensor<2,spacedim> > &hessians) const
05143 {
05144   get_function_hessians(fe_function, hessians);
05145 }
05146 
05147 
05148 
05149 template <int dim, int spacedim>
05150 template <class InputVector>
05151 inline
05152 void
05153 FEValuesBase<dim,spacedim>::
05154 get_function_2nd_derivatives (const InputVector                         &fe_function,
05155                               std::vector<std::vector<Tensor<2,spacedim> > > &hessians,
05156                               bool quadrature_points_fastest) const
05157 {
05158   get_function_hessians(fe_function, hessians, quadrature_points_fastest);
05159 }
05160 
05161 
05162 
05163 template <int dim, int spacedim>
05164 inline
05165 const Point<spacedim> &
05166 FEValuesBase<dim,spacedim>::normal_vector (const unsigned int i) const
05167 {
05168   typedef FEValuesBase<dim,spacedim> FVB;
05169   Assert (this->update_flags & update_normal_vectors,
05170           typename FVB::ExcAccessToUninitializedField());
05171   Assert (i<this->normal_vectors.size(),
05172           ExcIndexRange(i, 0, this->normal_vectors.size()));
05173 
05174   return this->normal_vectors[i];
05175 }
05176 
05177 
05178 
05179 template <int dim, int spacedim>
05180 inline
05181 const Point<spacedim> &
05182 FEValuesBase<dim,spacedim>::cell_normal_vector (const unsigned int i) const
05183 {
05184   return this->normal_vector(i);
05185 }
05186 
05187 
05188 
05189 
05190 /*------------------------ Inline functions: FEValues ----------------------------*/
05191 
05192 
05193 template <int dim, int spacedim>
05194 inline
05195 const Quadrature<dim> &
05196 FEValues<dim,spacedim>::get_quadrature () const
05197 {
05198   return quadrature;
05199 }
05200 
05201 
05202 
05203 template <int dim, int spacedim>
05204 inline
05205 const FEValues<dim,spacedim> &
05206 FEValues<dim,spacedim>::get_present_fe_values () const
05207 {
05208   return *this;
05209 }
05210 
05211 
05212 /*------------------------ Inline functions: FEFaceValuesBase --------------------*/
05213 
05214 
05215 template <int dim, int spacedim>
05216 inline
05217 unsigned int
05218 FEFaceValuesBase<dim,spacedim>::get_face_index () const
05219 {
05220   return present_face_index;
05221 }
05222 
05223 
05224 /*------------------------ Inline functions: FE*FaceValues --------------------*/
05225 
05226 template <int dim, int spacedim>
05227 inline
05228 const Quadrature<dim-1> &
05229 FEFaceValuesBase<dim,spacedim>::get_quadrature () const
05230 {
05231   return quadrature;
05232 }
05233 
05234 
05235 
05236 template <int dim, int spacedim>
05237 inline
05238 const FEFaceValues<dim,spacedim> &
05239 FEFaceValues<dim,spacedim>::get_present_fe_values () const
05240 {
05241   return *this;
05242 }
05243 
05244 
05245 
05246 template <int dim, int spacedim>
05247 inline
05248 const FESubfaceValues<dim,spacedim> &
05249 FESubfaceValues<dim,spacedim>::get_present_fe_values () const
05250 {
05251   return *this;
05252 }
05253 
05254 
05255 
05256 template <int dim, int spacedim>
05257 inline
05258 const Tensor<1,spacedim> &
05259 FEFaceValuesBase<dim,spacedim>::boundary_form (const unsigned int i) const
05260 {
05261   typedef FEValuesBase<dim,spacedim> FVB;
05262   Assert (i<this->boundary_forms.size(),
05263           ExcIndexRange(i, 0, this->boundary_forms.size()));
05264   Assert (this->update_flags & update_boundary_forms,
05265           typename FVB::ExcAccessToUninitializedField());
05266 
05267   return this->boundary_forms[i];
05268 }
05269 
05270 #endif // DOXYGEN
05271 
05272 DEAL_II_NAMESPACE_CLOSE
05273 
05274 #endif
05275 
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines

deal.II documentation generated on Tue May 22 2012 12:06:08 by doxygen 1.7.3