Reference documentation for deal.II version 8.1.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
fe_poly_face.templates.h
1 // ---------------------------------------------------------------------
2 // @f$Id: fe_poly_face.templates.h 30036 2013-07-18 16:55:32Z maier @f$
3 //
4 // Copyright (C) 2009 - 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 
18 #include <deal.II/base/qprojector.h>
19 #include <deal.II/base/polynomial_space.h>
20 #include <deal.II/fe/fe_values.h>
21 #include <deal.II/fe/fe_poly_face.h>
22 
23 
24 DEAL_II_NAMESPACE_OPEN
25 
26 template <class POLY, int dim, int spacedim>
28  const POLY &poly_space,
29  const FiniteElementData<dim> &fe_data,
30  const std::vector<bool> &restriction_is_additive_flags):
31  FiniteElement<dim,spacedim> (fe_data,
32  restriction_is_additive_flags,
33  std::vector<ComponentMask> (1, ComponentMask(1,true))),
34  poly_space(poly_space)
35 {
36  AssertDimension(dim, POLY::dimension+1);
37 }
38 
39 
40 template <class POLY, int dim, int spacedim>
41 unsigned int
43 {
44  return this->degree;
45 }
46 
47 
48 //---------------------------------------------------------------------------
49 // Auxiliary functions
50 //---------------------------------------------------------------------------
51 
52 
53 
54 
55 template <class POLY, int dim, int spacedim>
58 {
59  // for this kind of elements, only
60  // the values can be precomputed
61  // once and for all. set this flag
62  // if the values are requested at
63  // all
64  return update_default;
65 }
66 
67 
68 
69 template <class POLY, int dim, int spacedim>
72 {
73  UpdateFlags out = flags & update_values;
74  if (flags & update_gradients)
75  out |= update_gradients | update_covariant_transformation;
76  if (flags & update_hessians)
77  out |= update_hessians | update_covariant_transformation;
78  if (flags & update_cell_normal_vectors)
79  out |= update_cell_normal_vectors | update_JxW_values;
80 
81  return out;
82 }
83 
84 
85 
86 //---------------------------------------------------------------------------
87 // Data field initialization
88 //---------------------------------------------------------------------------
89 
90 template <class POLY, int dim, int spacedim>
93  const UpdateFlags,
94  const Mapping<dim,spacedim> &,
95  const Quadrature<dim> &) const
96 {
97  InternalData *data = new InternalData;
98  return data;
99 }
100 
101 
102 template <class POLY, int dim, int spacedim>
105  const UpdateFlags update_flags,
106  const Mapping<dim,spacedim> &,
107  const Quadrature<dim-1>& quadrature) const
108 {
109  // generate a new data object and
110  // initialize some fields
111  InternalData *data = new InternalData;
112 
113  // check what needs to be
114  // initialized only once and what
115  // on every cell/face/subface we
116  // visit
117  data->update_once = update_once(update_flags);
118  data->update_each = update_each(update_flags);
119  data->update_flags = data->update_once | data->update_each;
120 
121  const UpdateFlags flags(data->update_flags);
122  const unsigned int n_q_points = quadrature.size();
123 
124  // some scratch arrays
125  std::vector<double> values(0);
126  std::vector<Tensor<1,dim-1> > grads(0);
127  std::vector<Tensor<2,dim-1> > grad_grads(0);
128 
129  // initialize fields only if really
130  // necessary. otherwise, don't
131  // allocate memory
132  if (flags & update_values)
133  {
134  values.resize (poly_space.n());
135  data->shape_values.resize (poly_space.n(),
136  std::vector<double> (n_q_points));
137  for (unsigned int i=0; i<n_q_points; ++i)
138  {
139  poly_space.compute(quadrature.point(i),
140  values, grads, grad_grads);
141 
142  for (unsigned int k=0; k<poly_space.n(); ++k)
143  data->shape_values[k][i] = values[k];
144  }
145  }
146  // No derivatives of this element
147  // are implemented.
148  if (flags & update_gradients || flags & update_hessians)
149  {
150  Assert(false, ExcNotImplemented());
151  }
152 
153  return data;
154 }
155 
156 
157 
158 template <class POLY, int dim, int spacedim>
161  const UpdateFlags flags,
162  const Mapping<dim,spacedim> &mapping,
163  const Quadrature<dim-1>& quadrature) const
164 {
165  return get_face_data (flags, mapping,
167 }
168 
169 
170 
171 
172 
173 //---------------------------------------------------------------------------
174 // Fill data of FEValues
175 //---------------------------------------------------------------------------
176 template <class POLY, int dim, int spacedim>
177 void
181  const Quadrature<dim> &,
185  CellSimilarity::Similarity &) const
186 {
187  // Do nothing, since we do not have
188  // values in the interior
189 }
190 
191 
192 
193 template <class POLY, int dim, int spacedim>
194 void
196  const Mapping<dim,spacedim> &,
198  const unsigned int face,
199  const Quadrature<dim-1>& quadrature,
202  FEValuesData<dim,spacedim> &data) const
203 {
204  // convert data object to internal
205  // data for this class. fails with
206  // an exception if that is not
207  // possible
208  Assert (dynamic_cast<InternalData *> (&fedata) != 0, ExcInternalError());
209  InternalData &fe_data = static_cast<InternalData &> (fedata);
210 
211  const UpdateFlags flags(fe_data.update_once | fe_data.update_each);
212 
213  const unsigned int foffset = fe_data.shape_values.size() * face;
214  if (flags & update_values)
215  for (unsigned int i=0; i<quadrature.size(); ++i)
216  {
217  for (unsigned int k=0; k<this->dofs_per_cell; ++k)
218  data.shape_values(k,i) = 0.;
219  for (unsigned int k=0; k<fe_data.shape_values.size(); ++k)
220  data.shape_values(foffset+k,i) = fe_data.shape_values[k][i];
221  }
222 }
223 
224 
225 template <class POLY, int dim, int spacedim>
226 void
228  const Mapping<dim,spacedim> &,
230  const unsigned int face,
231  const unsigned int subface,
232  const Quadrature<dim-1>& quadrature,
235  FEValuesData<dim,spacedim> &data) const
236 {
237  // convert data object to internal
238  // data for this class. fails with
239  // an exception if that is not
240  // possible
241  Assert (dynamic_cast<InternalData *> (&fedata) != 0, ExcInternalError());
242  InternalData &fe_data = static_cast<InternalData &> (fedata);
243 
244  const UpdateFlags flags(fe_data.update_once | fe_data.update_each);
245 
246  const unsigned int foffset = fe_data.shape_values.size() * face;
247  const unsigned int offset = subface*quadrature.size();
248 
249  if (flags & update_values)
250  for (unsigned int i=0; i<quadrature.size(); ++i)
251  {
252  for (unsigned int k=0; k<this->dofs_per_cell; ++k)
253  data.shape_values(k,i) = 0.;
254  for (unsigned int k=0; k<fe_data.shape_values.size(); ++k)
255  data.shape_values(foffset+k,i) = fe_data.shape_values[k][i+offset];
256  }
257 
258  Assert (!(flags & update_gradients), ExcNotImplemented());
259  Assert (!(flags & update_hessians), ExcNotImplemented());
260 }
261 
262 DEAL_II_NAMESPACE_CLOSE
Transformed quadrature weights.
Shape function values.
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:858
std::vector< std::vector< double > > shape_values
Definition: fe_poly_face.h:345
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
unsigned int get_degree() const
No update.
unsigned int size() const
#define Assert(cond, exc)
Definition: exceptions.h:299
UpdateFlags
virtual UpdateFlags update_each(const UpdateFlags flags) const
Mapping< dim, spacedim >::InternalDataBase * get_face_data(const UpdateFlags, const Mapping< dim, spacedim > &mapping, const Quadrature< dim-1 > &quadrature) const
ShapeVector shape_values
Definition: fe_values.h:1171
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
FE_PolyFace(const POLY &poly_space, const FiniteElementData< dim > &fe_data, const std::vector< bool > &restriction_is_additive_flags)
virtual Mapping< dim, spacedim >::InternalDataBase * get_data(const UpdateFlags, const Mapping< dim, spacedim > &mapping, const Quadrature< dim > &quadrature) const
Second derivatives of shape functions.
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 UpdateFlags update_once(const UpdateFlags flags) const
UpdateFlags update_each
Definition: mapping.h:365
Mapping< dim, spacedim >::InternalDataBase * get_subface_data(const UpdateFlags, const Mapping< dim, spacedim > &mapping, const Quadrature< dim-1 > &quadrature) const
UpdateFlags update_once
Definition: mapping.h:359
Shape function gradients.
::ExceptionBase & ExcNotImplemented()
UpdateFlags update_flags
Definition: mapping.h:353
::ExceptionBase & ExcInternalError()
Covariant transformation.
const Point< dim > & point(const unsigned int i) const