Reference documentation for deal.II version 8.1.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
mapping_q_eulerian.h
1 // ---------------------------------------------------------------------
2 // @f$Id: mapping_q_eulerian.h 30450 2013-08-23 15:48:29Z kronbichler @f$
3 //
4 // Copyright (C) 2001 - 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 #ifndef __deal2__mapping_q_eulerian_h
19 #define __deal2__mapping_q_eulerian_h
20 
21 #include <deal.II/base/smartpointer.h>
22 #include <deal.II/base/thread_management.h>
23 #include <deal.II/grid/tria_iterator.h>
24 #include <deal.II/dofs/dof_handler.h>
25 #include <deal.II/dofs/dof_accessor.h>
26 #include <deal.II/fe/fe.h>
27 #include <deal.II/fe/fe_values.h>
28 #include <deal.II/fe/mapping_q.h>
29 
30 
31 DEAL_II_NAMESPACE_OPEN
32 
33 
36 
93 template <int dim, class VECTOR = Vector<double>, int spacedim=dim >
94 class MappingQEulerian : public MappingQ<dim, spacedim>
95 {
96 public:
109  MappingQEulerian (const unsigned int degree,
110  const VECTOR &euler_vector,
112 
117  virtual
118  Mapping<dim,spacedim> *clone () const;
119 
125  bool preserves_vertex_locations () const;
126 
130  DeclException0 (ExcInactiveCell);
131 
132 protected:
137  virtual void
139  const Quadrature<dim> &quadrature,
140  typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
141  typename std::vector<Point<spacedim> > &quadrature_points,
142  std::vector<double> &JxW_values,
143  std::vector<DerivativeForm<1,dim,spacedim> > &jacobians,
144  std::vector<DerivativeForm<2,dim,spacedim> > &jacobian_grads,
145  std::vector<DerivativeForm<1,spacedim,dim> > &inverse_jacobians,
146  std::vector<Point<spacedim> > &cell_normal_vectors,
147  CellSimilarity::Similarity &cell_similarity) const;
148 
154 
159 
160 
161 private:
162 
168  class SupportQuadrature : public Quadrature<dim>
169  {
170  public:
175  SupportQuadrature (const unsigned int map_degree);
176 
177  };
178 
183 
192 
197 
201  virtual void compute_mapping_support_points(
202  const typename Triangulation<dim,spacedim>::cell_iterator &cell,
203  std::vector<Point<spacedim> > &a) const;
204 
205 };
206 
210 /*----------------------------------------------------------------------*/
211 
212 #ifndef DOXYGEN
213 
214 template <int dim, class VECTOR, int spacedim>
215 inline
216 bool
218 {
219  return false;
220 }
221 
222 #endif // DOXYGEN
223 
224 
225 DEAL_II_NAMESPACE_CLOSE
226 
227 
228 #endif // __deal2__mapping_q_eulerian_h
229 
DeclException0(ExcInactiveCell)
virtual void compute_mapping_support_points(const typename Triangulation< dim, spacedim >::cell_iterator &cell, std::vector< Point< spacedim > > &a) const
virtual void fill_fe_values(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Quadrature< dim > &quadrature, typename Mapping< dim, spacedim >::InternalDataBase &mapping_data, typename std::vector< Point< spacedim > > &quadrature_points, std::vector< double > &JxW_values, std::vector< DerivativeForm< 1, dim, spacedim > > &jacobians, std::vector< DerivativeForm< 2, dim, spacedim > > &jacobian_grads, std::vector< DerivativeForm< 1, spacedim, dim > > &inverse_jacobians, std::vector< Point< spacedim > > &cell_normal_vectors, CellSimilarity::Similarity &cell_similarity) const
FEValues< dim, spacedim > fe_values
const SupportQuadrature support_quadrature
bool preserves_vertex_locations() const
MappingQEulerian(const unsigned int degree, const VECTOR &euler_vector, const DoFHandler< dim, spacedim > &euler_dof_handler)
const unsigned int degree
Definition: mapping_q.h:393
Threads::Mutex fe_values_mutex
SupportQuadrature(const unsigned int map_degree)
SmartPointer< const DoFHandler< dim, spacedim >, MappingQEulerian< dim, VECTOR, spacedim > > euler_dof_handler
SmartPointer< const VECTOR, MappingQEulerian< dim, VECTOR, spacedim > > euler_vector
virtual Mapping< dim, spacedim > * clone() const