16 #include <deal.II/base/std_cxx14/memory.h> 18 #include <deal.II/dofs/dof_accessor.h> 19 #include <deal.II/dofs/dof_handler.h> 21 #include <deal.II/fe/fe.h> 22 #include <deal.II/fe/mapping_q1_eulerian.h> 24 #include <deal.II/grid/tria_iterator.h> 26 #include <deal.II/lac/block_vector.h> 27 #include <deal.II/lac/la_parallel_block_vector.h> 28 #include <deal.II/lac/la_parallel_vector.h> 29 #include <deal.II/lac/la_vector.h> 30 #include <deal.II/lac/petsc_block_vector.h> 31 #include <deal.II/lac/petsc_vector.h> 32 #include <deal.II/lac/trilinos_parallel_block_vector.h> 33 #include <deal.II/lac/trilinos_vector.h> 34 #include <deal.II/lac/vector.h> 38 DEAL_II_NAMESPACE_OPEN
41 template <
int dim,
class VectorType,
int spacedim>
44 const VectorType & euler_transform_vectors)
46 , euler_transform_vectors(&euler_transform_vectors)
47 , shiftmap_dof_handler(&shiftmap_dof_handler)
52 template <
int dim,
class VectorType,
int spacedim>
83 for (
unsigned int i = 0; i < GeometryInfo<dim>::vertices_per_cell; ++i)
90 for (
unsigned int j = 0; j < spacedim; ++j)
91 shift_vector[j] = mapping_values(i * spacedim + j);
95 vertices[i] = cell->vertex(i) + shift_vector;
102 template <
int dim,
class VectorType,
int spacedim>
103 std::vector<Point<spacedim>>
111 for (
unsigned int i = 0; i < GeometryInfo<dim>::vertices_per_cell; ++i)
119 template <
int dim,
class VectorType,
int spacedim>
120 std::unique_ptr<Mapping<dim, spacedim>>
123 return std_cxx14::make_unique<MappingQ1Eulerian<dim, VectorType, spacedim>>(
129 template <
int dim,
class VectorType,
int spacedim>
157 #include "mapping_q1_eulerian.inst" 160 DEAL_II_NAMESPACE_CLOSE
MappingQ1Eulerian(const DoFHandler< dim, spacedim > &euler_dof_handler, const VectorType &euler_vector)
virtual CellSimilarity::Similarity fill_fe_values(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const CellSimilarity::Similarity cell_similarity, const Quadrature< dim > &quadrature, const typename Mapping< dim, spacedim >::InternalDataBase &internal_data, internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &output_data) const override
#define AssertDimension(dim1, dim2)
static::ExceptionBase & ExcInactiveCell()
virtual std::vector< Point< spacedim > > compute_mapping_support_points(const typename Triangulation< dim, spacedim >::cell_iterator &cell) const override
#define Assert(cond, exc)
SmartPointer< const DoFHandler< dim, spacedim >, MappingQ1Eulerian< dim, VectorType, spacedim > > shiftmap_dof_handler
SmartPointer< const VectorType, MappingQ1Eulerian< dim, VectorType, spacedim > > euler_transform_vectors
virtual std::unique_ptr< Mapping< dim, spacedim > > clone() const override
typename ActiveSelector::cell_iterator cell_iterator
virtual CellSimilarity::Similarity fill_fe_values(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const CellSimilarity::Similarity cell_similarity, const Quadrature< dim > &quadrature, const typename Mapping< dim, spacedim >::InternalDataBase &internal_data,::internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &output_data) const override
virtual std::array< Point< spacedim >, GeometryInfo< dim >::vertices_per_cell > get_vertices(const typename Triangulation< dim, spacedim >::cell_iterator &cell) const override