deal.II version GIT relicensing-2901-g19332422bd 2025-03-23 19:50:00+00:00
\(\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// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 2001 - 2023 by the deal.II authors
5//
6// This file is part of the deal.II library.
7//
8// Part of the source code is dual licensed under Apache-2.0 WITH
9// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
10// governing the source code and code contributions can be found in
11// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
12//
13// ------------------------------------------------------------------------
14
15
18
19#include <deal.II/fe/fe.h>
21
23
31#include <deal.II/lac/vector.h>
32
33#include <boost/container/small_vector.hpp>
34
35#include <array>
36#include <memory>
37
38
40
41
42template <int dim, typename VectorType, int spacedim>
44 const DoFHandler<dim, spacedim> &shiftmap_dof_handler,
45 const VectorType &euler_transform_vectors)
46 : MappingQ<dim, spacedim>(1)
47 , euler_transform_vectors(&euler_transform_vectors)
48 , shiftmap_dof_handler(&shiftmap_dof_handler)
49{}
50
51
52
53template <int dim, typename VectorType, int spacedim>
54boost::container::small_vector<Point<spacedim>,
55#ifndef _MSC_VER
56 ReferenceCells::max_n_vertices<dim>()
57#else
59#endif
60 >
62 const typename Triangulation<dim, spacedim>::cell_iterator &cell) const
63{
64 boost::container::small_vector<Point<spacedim>,
65#ifndef _MSC_VER
66 ReferenceCells::max_n_vertices<dim>()
67#else
69#endif
70 >
71 vertices(cell->n_vertices());
72 // The assertions can not be in the constructor, since this would
73 // require to call dof_handler.distribute_dofs(fe) *before* the mapping
74 // object is constructed, which is not necessarily what we want.
75
76 // TODO: Only one of these two assertions should be relevant
77 AssertDimension(spacedim, shiftmap_dof_handler->get_fe().n_dofs_per_vertex());
78 AssertDimension(shiftmap_dof_handler->get_fe(0).n_components(), spacedim);
79
80 AssertDimension(shiftmap_dof_handler->n_dofs(),
81 euler_transform_vectors->size());
82
83 // cast the Triangulation<dim>::cell_iterator into a
84 // DoFHandler<dim>::cell_iterator which is necessary for access to
85 // DoFCellAccessor::get_dof_values()
87 *cell, shiftmap_dof_handler);
88
89 // We require the cell to be active since we can only then get nodal
90 // values for the shifts
91 Assert(dof_cell->is_active() == true, ExcInactiveCell());
92
93 // now get the values of the shift vectors at the vertices
95 shiftmap_dof_handler->get_fe().n_dofs_per_cell());
96 dof_cell->get_dof_values(*euler_transform_vectors, mapping_values);
97
98 for (const unsigned int i : GeometryInfo<dim>::vertex_indices())
99 {
100 Point<spacedim> shift_vector;
101
102 // pick out the value of the shift vector at the present
103 // vertex. since vertex dofs are always numbered first, we can
104 // access them easily
105 for (unsigned int j = 0; j < spacedim; ++j)
106 shift_vector[j] = mapping_values(i * spacedim + j);
107
108 // compute new support point by old (reference) value and added
109 // shift
110 vertices[i] = cell->vertex(i) + shift_vector;
111 }
112 return vertices;
113}
114
115
116
117template <int dim, typename VectorType, int spacedim>
118std::vector<Point<spacedim>>
120 const typename Triangulation<dim, spacedim>::cell_iterator &cell) const
121{
122 const auto vertices = this->get_vertices(cell);
123
124 std::vector<Point<spacedim>> a(GeometryInfo<dim>::vertices_per_cell);
125 for (const unsigned int i : GeometryInfo<dim>::vertex_indices())
126 a[i] = vertices[i];
127
128 return a;
129}
130
131
132
133template <int dim, typename VectorType, int spacedim>
134std::unique_ptr<Mapping<dim, spacedim>>
136{
137 return std::make_unique<MappingQ1Eulerian<dim, VectorType, spacedim>>(*this);
138}
139
140
141
142template <int dim, typename VectorType, int spacedim>
147 const Quadrature<dim> &quadrature,
148 const typename Mapping<dim, spacedim>::InternalDataBase &internal_data,
150 &output_data) const
151{
152 // call the function of the base class, but ignoring
153 // any potentially detected cell similarity between
154 // the current and the previous cell
157 quadrature,
158 internal_data,
159 output_data);
160 // also return the updated flag since any detected
161 // similarity wasn't based on the mapped field, but
162 // the original vertices which are meaningless
164}
165
166
167
168// explicit instantiations
169#include "fe/mapping_q1_eulerian.inst"
170
171
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 boost::container::small_vector< Point< spacedim >, ReferenceCells::max_n_vertices< dim >() > get_vertices(const typename Triangulation< dim, spacedim >::cell_iterator &cell) 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
Definition mapping_q.cc:833
Definition point.h:113
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:40
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:41
#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()