Reference documentation for deal.II version 8.1.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Modules Pages
fe_poly_tensor.h
1 // ---------------------------------------------------------------------
2 // @f$Id: fe_poly_tensor.h 30036 2013-07-18 16:55:32Z maier @f$
3 //
4 // Copyright (C) 2005 - 2013 by the deal.II authors
5 //
6 // This file is part of the deal.II library.
7 //
8 // The deal.II library is free software; you can use it, redistribute
9 // it, and/or modify it under the terms of the GNU Lesser General
10 // Public License as published by the Free Software Foundation; either
11 // version 2.1 of the License, or (at your option) any later version.
12 // The full text of the license can be found in the file LICENSE at
13 // the top level of the deal.II distribution.
14 //
15 // ---------------------------------------------------------------------
16 
17 #ifndef __deal2__fe_poly_tensor_h
18 #define __deal2__fe_poly_tensor_h
19 
20 
21 #include <deal.II/lac/full_matrix.h>
22 #include <deal.II/fe/fe.h>
23 #include <deal.II/base/derivative_form.h>
24 
26 
119 template <class POLY, int dim, int spacedim=dim>
120 class FE_PolyTensor : public FiniteElement<dim,spacedim>
121 {
122 public:
131  FE_PolyTensor (const unsigned int degree,
132  const FiniteElementData<dim> &fe_data,
133  const std::vector<bool> &restriction_is_additive_flags,
134  const std::vector<ComponentMask> &nonzero_components);
135 
141  virtual double shape_value (const unsigned int i,
142  const Point<dim> &p) const;
143 
144  virtual double shape_value_component (const unsigned int i,
145  const Point<dim> &p,
146  const unsigned int component) const;
147 
153  virtual Tensor<1,dim> shape_grad (const unsigned int i,
154  const Point<dim> &p) const;
155 
156  virtual Tensor<1,dim> shape_grad_component (const unsigned int i,
157  const Point<dim> &p,
158  const unsigned int component) const;
159 
165  virtual Tensor<2,dim> shape_grad_grad (const unsigned int i,
166  const Point<dim> &p) const;
167 
168  virtual Tensor<2,dim> shape_grad_grad_component (const unsigned int i,
169  const Point<dim> &p,
170  const unsigned int component) const;
171 
182  virtual UpdateFlags update_once (const UpdateFlags flags) const;
193  virtual UpdateFlags update_each (const UpdateFlags flags) const;
194 
195 protected:
203 
204  virtual
206  get_data (const UpdateFlags,
207  const Mapping<dim,spacedim> &mapping,
208  const Quadrature<dim> &quadrature) const ;
209 
210  virtual void
211  fill_fe_values (const Mapping<dim,spacedim> &mapping,
212  const typename Triangulation<dim,spacedim>::cell_iterator &cell,
213  const Quadrature<dim> &quadrature,
214  typename Mapping<dim,spacedim>::InternalDataBase &mapping_internal,
215  typename Mapping<dim,spacedim>::InternalDataBase &fe_internal,
217  CellSimilarity::Similarity &cell_similarity) const;
218 
219  virtual void
221  const typename Triangulation<dim,spacedim>::cell_iterator &cell,
222  const unsigned int face_no,
223  const Quadrature<dim-1> &quadrature,
224  typename Mapping<dim,spacedim>::InternalDataBase &mapping_internal,
225  typename Mapping<dim,spacedim>::InternalDataBase &fe_internal,
226  FEValuesData<dim,spacedim> &data) const ;
227 
228  virtual void
230  const typename Triangulation<dim,spacedim>::cell_iterator &cell,
231  const unsigned int face_no,
232  const unsigned int sub_no,
233  const Quadrature<dim-1> &quadrature,
234  typename Mapping<dim,spacedim>::InternalDataBase &mapping_internal,
235  typename Mapping<dim,spacedim>::InternalDataBase &fe_internal,
236  FEValuesData<dim,spacedim> &data) const ;
237 
253  class InternalData : public FiniteElement<dim,spacedim>::InternalDataBase
254  {
255  public:
265  std::vector<std::vector<Tensor<1,dim> > > shape_values;
266 
276  std::vector< std::vector< DerivativeForm<1, dim, spacedim> > > shape_grads;
277  };
278 
285 
307 
319 
325  mutable std::vector<Tensor<1,dim> > cached_values;
326 
332  mutable std::vector<Tensor<2,dim> > cached_grads;
333 
339  mutable std::vector<Tensor<3,dim> > cached_grad_grads;
340 };
341 
342 DEAL_II_NAMESPACE_CLOSE
343 
344 #endif
virtual UpdateFlags update_each(const UpdateFlags flags) const
MappingType
Definition: mapping.h:58
Point< dim > cached_point
const std::vector< bool > restriction_is_additive_flags
Definition: fe.h:1909
std::vector< Tensor< 2, dim > > cached_grads
virtual double shape_value_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
virtual UpdateFlags update_once(const UpdateFlags flags) const
std::vector< std::vector< Tensor< 1, dim > > > shape_values
virtual Tensor< 2, dim > shape_grad_grad(const unsigned int i, const Point< dim > &p) const
virtual Mapping< dim, spacedim >::InternalDataBase * get_data(const UpdateFlags, const Mapping< dim, spacedim > &mapping, const Quadrature< dim > &quadrature) const
virtual Tensor< 1, dim > shape_grad(const unsigned int i, const Point< dim > &p) const
std::vector< Tensor< 3, dim > > cached_grad_grads
MappingType mapping_type
std::vector< std::vector< DerivativeForm< 1, dim, spacedim > > > shape_grads
UpdateFlags
const unsigned int degree
Definition: fe_base.h:287
const std::vector< ComponentMask > nonzero_components
Definition: fe.h:1918
virtual double shape_value(const unsigned int i, const Point< dim > &p) const
virtual void fill_fe_face_values(const Mapping< dim, spacedim > &mapping, const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const Quadrature< dim-1 > &quadrature, typename Mapping< dim, spacedim >::InternalDataBase &mapping_internal, typename Mapping< dim, spacedim >::InternalDataBase &fe_internal, FEValuesData< dim, spacedim > &data) const
virtual void fill_fe_values(const Mapping< dim, spacedim > &mapping, const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Quadrature< dim > &quadrature, typename Mapping< dim, spacedim >::InternalDataBase &mapping_internal, typename Mapping< dim, spacedim >::InternalDataBase &fe_internal, FEValuesData< dim, spacedim > &data, CellSimilarity::Similarity &cell_similarity) const
FE_PolyTensor(const unsigned int degree, const FiniteElementData< dim > &fe_data, const std::vector< bool > &restriction_is_additive_flags, const std::vector< ComponentMask > &nonzero_components)
virtual void fill_fe_subface_values(const Mapping< dim, spacedim > &mapping, const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const unsigned int sub_no, const Quadrature< dim-1 > &quadrature, typename Mapping< dim, spacedim >::InternalDataBase &mapping_internal, typename Mapping< dim, spacedim >::InternalDataBase &fe_internal, FEValuesData< dim, spacedim > &data) const
FullMatrix< double > inverse_node_matrix
std::vector< Tensor< 1, dim > > cached_values
virtual Tensor< 2, dim > shape_grad_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
virtual Tensor< 1, dim > shape_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const