Reference documentation for deal.II version 9.5.0
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
Loading...
Searching...
No Matches
mapping_q1_eulerian.cc
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2001 - 2021 by the deal.II authors
4//
5// This file is part of the deal.II library.
6//
7// The deal.II library is free software; you can use it, redistribute
8// it, and/or modify it under the terms of the GNU Lesser General
9// Public License as published by the Free Software Foundation; either
10// version 2.1 of the License, or (at your option) any later version.
11// The full text of the license can be found in the file LICENSE.md at
12// the top level directory of deal.II.
13//
14// ---------------------------------------------------------------------
15
16
19
20#include <deal.II/fe/fe.h>
22
24
33#include <deal.II/lac/vector.h>
34
35#include <array>
36#include <memory>
37
39
40
41template <int dim, class VectorType, int spacedim>
43 const DoFHandler<dim, spacedim> &shiftmap_dof_handler,
44 const VectorType & euler_transform_vectors)
45 : MappingQ<dim, spacedim>(1)
46 , euler_transform_vectors(&euler_transform_vectors)
47 , shiftmap_dof_handler(&shiftmap_dof_handler)
48{}
49
50
51
52template <int dim, class VectorType, int spacedim>
53boost::container::small_vector<Point<spacedim>,
56 const typename Triangulation<dim, spacedim>::cell_iterator &cell) const
57{
58 boost::container::small_vector<Point<spacedim>,
61 // The assertions can not be in the constructor, since this would
62 // require to call dof_handler.distribute_dofs(fe) *before* the mapping
63 // object is constructed, which is not necessarily what we want.
64
65 // TODO: Only one of these two assertions should be relevant
66 AssertDimension(spacedim, shiftmap_dof_handler->get_fe().n_dofs_per_vertex());
67 AssertDimension(shiftmap_dof_handler->get_fe(0).n_components(), spacedim);
68
69 AssertDimension(shiftmap_dof_handler->n_dofs(),
70 euler_transform_vectors->size());
71
72 // cast the Triangulation<dim>::cell_iterator into a
73 // DoFHandler<dim>::cell_iterator which is necessary for access to
74 // DoFCellAccessor::get_dof_values()
76 *cell, shiftmap_dof_handler);
77
78 // We require the cell to be active since we can only then get nodal
79 // values for the shifts
80 Assert(dof_cell->is_active() == true, ExcInactiveCell());
81
82 // now get the values of the shift vectors at the vertices
84 shiftmap_dof_handler->get_fe().n_dofs_per_cell());
85 dof_cell->get_dof_values(*euler_transform_vectors, mapping_values);
86
87 for (const unsigned int i : GeometryInfo<dim>::vertex_indices())
88 {
89 Point<spacedim> shift_vector;
90
91 // pick out the value of the shift vector at the present
92 // vertex. since vertex dofs are always numbered first, we can
93 // access them easily
94 for (unsigned int j = 0; j < spacedim; ++j)
95 shift_vector[j] = mapping_values(i * spacedim + j);
96
97 // compute new support point by old (reference) value and added
98 // shift
99 vertices[i] = cell->vertex(i) + shift_vector;
100 }
101 return vertices;
102}
103
104
105
106template <int dim, class VectorType, int spacedim>
107std::vector<Point<spacedim>>
109 const typename Triangulation<dim, spacedim>::cell_iterator &cell) const
110{
111 const auto vertices = this->get_vertices(cell);
112
113 std::vector<Point<spacedim>> a(GeometryInfo<dim>::vertices_per_cell);
114 for (const unsigned int i : GeometryInfo<dim>::vertex_indices())
115 a[i] = vertices[i];
116
117 return a;
118}
119
120
121
122template <int dim, class VectorType, int spacedim>
123std::unique_ptr<Mapping<dim, spacedim>>
125{
126 return std::make_unique<MappingQ1Eulerian<dim, VectorType, spacedim>>(*this);
127}
128
129
130
131template <int dim, class VectorType, int spacedim>
136 const Quadrature<dim> & quadrature,
137 const typename Mapping<dim, spacedim>::InternalDataBase &internal_data,
139 &output_data) const
140{
141 // call the function of the base class, but ignoring
142 // any potentially detected cell similarity between
143 // the current and the previous cell
146 quadrature,
147 internal_data,
148 output_data);
149 // also return the updated flag since any detected
150 // similarity wasn't based on the mapped field, but
151 // the original vertices which are meaningless
153}
154
155
156
157// explicit instantiations
158#include "mapping_q1_eulerian.inst"
159
160
virtual boost::container::small_vector< Point< spacedim >, GeometryInfo< dim >::vertices_per_cell > get_vertices(const typename Triangulation< dim, spacedim >::cell_iterator &cell) const override
virtual std::vector< Point< spacedim > > compute_mapping_support_points(const typename Triangulation< dim, spacedim >::cell_iterator &cell) const override
virtual std::unique_ptr< Mapping< dim, spacedim > > clone() const override
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
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
Definition mapping_q.cc:786
Definition point.h:112
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
Point< 3 > vertices[4]
#define Assert(cond, exc)
#define AssertDimension(dim1, dim2)
typename ActiveSelector::cell_iterator cell_iterator
static std_cxx20::ranges::iota_view< unsigned int, unsigned int > vertex_indices()