Loading [MathJax]/extensions/TeX/newcommand.js
 Reference documentation for deal.II version 9.3.3
\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\}}
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages
utilities.h
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2020 - 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#ifndef dealii_particles_utilities
17#define dealii_particles_utilities
18
19#include <deal.II/base/config.h>
20
22#include <deal.II/base/point.h>
24
26
28#include <deal.II/fe/fe.h>
30
32
34
36
37
39
40
41
42namespace Particles
43{
48 namespace Utilities
49 {
94 template <int dim,
95 int spacedim,
96 typename SparsityType,
97 typename number = double>
98 void
100 const DoFHandler<dim, spacedim> & space_dh,
101 const Particles::ParticleHandler<dim, spacedim> &particle_handler,
102 SparsityType & sparsity,
103 const AffineConstraints<number> & constraints =
105 const ComponentMask &space_comps = ComponentMask());
106
150 template <int dim, int spacedim, typename MatrixType>
151 void
153 const DoFHandler<dim, spacedim> & space_dh,
154 const Particles::ParticleHandler<dim, spacedim> &particle_handler,
155 MatrixType & matrix,
158 const ComponentMask &space_comps = ComponentMask());
159
182 template <int dim,
183 int spacedim,
184 typename InputVectorType,
185 typename OutputVectorType>
186 void
188 const DoFHandler<dim, spacedim> & field_dh,
189 const Particles::ParticleHandler<dim, spacedim> &particle_handler,
190 const InputVectorType & field_vector,
191 OutputVectorType & interpolated_field,
192 const ComponentMask &field_comps = ComponentMask())
193 {
194 if (particle_handler.n_locally_owned_particles() == 0)
195 {
196 interpolated_field.compress(VectorOperation::add);
197 return; // nothing else to do here
198 }
199
200 const auto &tria = field_dh.get_triangulation();
201 const auto &fe = field_dh.get_fe();
202 auto particle = particle_handler.begin();
203
204 // Take care of components
205 const ComponentMask comps =
206 (field_comps.size() == 0 ? ComponentMask(fe.n_components(), true) :
207 field_comps);
208 AssertDimension(comps.size(), fe.n_components());
209 const auto n_comps = comps.n_selected_components();
210
211 AssertDimension(field_vector.size(), field_dh.n_dofs());
212 AssertDimension(interpolated_field.size(),
213 particle_handler.get_next_free_particle_index() *
214 n_comps);
215
216 // Global to local indices
217 std::vector<unsigned int> space_gtl(fe.n_components(),
219 for (unsigned int i = 0, j = 0; i < space_gtl.size(); ++i)
220 if (comps[i])
221 space_gtl[i] = j++;
222
223 std::vector<types::global_dof_index> dof_indices(fe.n_dofs_per_cell());
224
225 while (particle != particle_handler.end())
226 {
227 const auto &cell = particle->get_surrounding_cell(tria);
228 const auto &dh_cell =
229 typename DoFHandler<dim, spacedim>::cell_iterator(*cell, &field_dh);
230 dh_cell->get_dof_indices(dof_indices);
231 const auto pic = particle_handler.particles_in_cell(cell);
232
233 Assert(pic.begin() == particle, ExcInternalError());
234 for (unsigned int i = 0; particle != pic.end(); ++particle, ++i)
235 {
236 const Point<dim> reference_location =
237 particle->get_reference_location();
238
239 const auto id = particle->get_id();
240
241 for (unsigned int j = 0; j < fe.n_dofs_per_cell(); ++j)
242 {
243 const auto comp_j =
244 space_gtl[fe.system_to_component_index(j).first];
245 if (comp_j != numbers::invalid_unsigned_int)
246 interpolated_field[id * n_comps + comp_j] +=
247 fe.shape_value(j, reference_location) *
248 field_vector(dof_indices[j]);
249 }
250 }
251 }
252 interpolated_field.compress(VectorOperation::add);
253 }
254
255 } // namespace Utilities
256} // namespace Particles
258
259#endif
unsigned int size() const
unsigned int n_selected_components(const unsigned int overall_number_of_components=numbers::invalid_unsigned_int) const
const FiniteElement< dim, spacedim > & get_fe(const unsigned int index=0) const
const Triangulation< dim, spacedim > & get_triangulation() const
types::global_dof_index n_dofs() const
particle_iterator begin() const
particle_iterator end() const
types::particle_index n_locally_owned_particles() const
particle_iterator_range particles_in_cell(const typename Triangulation< dim, spacedim >::active_cell_iterator &cell)
types::particle_index get_next_free_particle_index() const
Definition: point.h:111
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:402
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:403
#define Assert(cond, exc)
Definition: exceptions.h:1465
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1622
static ::ExceptionBase & ExcInternalError()
typename ActiveSelector::cell_iterator cell_iterator
Definition: dof_handler.h:466
@ matrix
Contents is actually a matrix.
void interpolate_field_on_particles(const DoFHandler< dim, spacedim > &field_dh, const Particles::ParticleHandler< dim, spacedim > &particle_handler, const InputVectorType &field_vector, OutputVectorType &interpolated_field, const ComponentMask &field_comps=ComponentMask())
Definition: utilities.h:187
void create_interpolation_sparsity_pattern(const DoFHandler< dim, spacedim > &space_dh, const Particles::ParticleHandler< dim, spacedim > &particle_handler, SparsityType &sparsity, const AffineConstraints< number > &constraints=AffineConstraints< number >(), const ComponentMask &space_comps=ComponentMask())
Definition: utilities.cc:32
void create_interpolation_matrix(const DoFHandler< dim, spacedim > &space_dh, const Particles::ParticleHandler< dim, spacedim > &particle_handler, MatrixType &matrix, const AffineConstraints< typename MatrixType::value_type > &constraints=AffineConstraints< typename MatrixType::value_type >(), const ComponentMask &space_comps=ComponentMask())
Definition: utilities.cc:115
static const unsigned int invalid_unsigned_int
Definition: types.h:196