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
utilities.cc
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2020 - 2022 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#include <deal.II/base/config.h>
17
19
21
22
24
25
26namespace Particles
27{
28 namespace Utilities
29 {
30 template <int dim, int spacedim, typename number>
31 void
33 const DoFHandler<dim, spacedim> & space_dh,
34 const Particles::ParticleHandler<dim, spacedim> &particle_handler,
35 SparsityPatternBase & sparsity,
36 const AffineConstraints<number> & constraints,
37 const ComponentMask & space_comps)
38 {
39 if (particle_handler.n_locally_owned_particles() == 0)
40 return; // nothing to do here
41
42 const auto &fe = space_dh.get_fe();
43 const auto max_particles_per_cell =
44 particle_handler.n_global_max_particles_per_cell();
45
46 // Take care of components
47 const ComponentMask comps =
48 (space_comps.size() == 0 ? ComponentMask(fe.n_components(), true) :
49 space_comps);
50 AssertDimension(comps.size(), fe.n_components());
51
52 const auto n_comps = comps.n_selected_components();
53
54 // Global to local indices
55 std::vector<unsigned int> space_gtl(fe.n_components(),
57 for (unsigned int i = 0, j = 0; i < space_gtl.size(); ++i)
58 if (comps[i])
59 space_gtl[i] = j++;
60
61 // [TODO]: when the add_entries_local_to_global below will implement
62 // the version with the dof_mask, this should be uncommented.
63 // // Construct a dof_mask, used to distribute entries to the sparsity
64 // Table<2, bool> dof_mask(max_particles_per_cell * n_comps,
65 // fe.n_dofs_per_cell());
66 // dof_mask.fill(false);
67 // for (unsigned int i = 0; i < space_fe.n_dofs_per_cell(); ++i)
68 // {
69 // const auto comp_i = space_fe.system_to_component_index(i).first;
70 // if (space_gtl[comp_i] != numbers::invalid_unsigned_int)
71 // for (unsigned int j = 0; j < max_particles_per_cell; ++j)
72 // dof_mask(i, j * n_comps + space_gtl[comp_i]) = true;
73 // }
74
75 std::vector<types::global_dof_index> dof_indices(fe.n_dofs_per_cell());
76 std::vector<types::particle_index> particle_indices(
77 max_particles_per_cell * n_comps);
78
79 auto particle = particle_handler.begin();
80 while (particle != particle_handler.end())
81 {
82 const auto &cell = particle->get_surrounding_cell();
83 const auto dh_cell =
84 typename DoFHandler<dim, spacedim>::cell_iterator(*cell, &space_dh);
85 dh_cell->get_dof_indices(dof_indices);
86 const auto pic = particle_handler.particles_in_cell(cell);
87 const auto n_particles = particle_handler.n_particles_in_cell(cell);
88 particle_indices.resize(n_particles * n_comps);
89 Assert(pic.begin() == particle, ExcInternalError());
90 for (unsigned int i = 0; particle != pic.end(); ++particle, ++i)
91 {
92 const auto p_id = particle->get_id();
93 for (unsigned int j = 0; j < fe.n_dofs_per_cell(); ++j)
94 {
95 const auto comp_j =
96 space_gtl[fe.system_to_component_index(j).first];
99 {p_id * n_comps + comp_j}, {dof_indices[j]}, sparsity);
100 }
101 }
102 // [TODO]: when this works, use this:
103 // constraints.add_entries_local_to_global(particle_indices,
104 // dof_indices,
105 // sparsity,
106 // dof_mask);
107 }
108 }
109
110
111
112 template <int dim, int spacedim, typename MatrixType>
113 void
115 const DoFHandler<dim, spacedim> & space_dh,
116 const Particles::ParticleHandler<dim, spacedim> &particle_handler,
117 MatrixType & matrix,
119 const ComponentMask & space_comps)
120 {
121 if (particle_handler.n_locally_owned_particles() == 0)
122 {
123 matrix.compress(VectorOperation::add);
124 return; // nothing else to do here
125 }
126
127 AssertDimension(matrix.n(), space_dh.n_dofs());
128
129 const auto &fe = space_dh.get_fe();
130 const auto max_particles_per_cell =
131 particle_handler.n_global_max_particles_per_cell();
132
133 // Take care of components
134 const ComponentMask comps =
135 (space_comps.size() == 0 ? ComponentMask(fe.n_components(), true) :
136 space_comps);
137 AssertDimension(comps.size(), fe.n_components());
138 const auto n_comps = comps.n_selected_components();
139
140 AssertDimension(matrix.m(),
141 particle_handler.n_global_particles() * n_comps);
142
143
144 // Global to local indices
145 std::vector<unsigned int> space_gtl(fe.n_components(),
147 for (unsigned int i = 0, j = 0; i < space_gtl.size(); ++i)
148 if (comps[i])
149 space_gtl[i] = j++;
150
151 // [TODO]: when the add_entries_local_to_global below will implement
152 // the version with the dof_mask, this should be uncommented.
153 // // Construct a dof_mask, used to distribute entries to the sparsity
154 // Table<2, bool> dof_mask(max_particles_per_cell * n_comps,
155 // fe.n_dofs_per_cell());
156 // dof_mask.fill(false);
157 // for (unsigned int i = 0; i < space_fe.n_dofs_per_cell(); ++i)
158 // {
159 // const auto comp_i = space_fe.system_to_component_index(i).first;
160 // if (space_gtl[comp_i] != numbers::invalid_unsigned_int)
161 // for (unsigned int j = 0; j < max_particles_per_cell; ++j)
162 // dof_mask(i, j * n_comps + space_gtl[comp_i]) = true;
163 // }
164
165 std::vector<types::global_dof_index> dof_indices(fe.n_dofs_per_cell());
166 std::vector<types::particle_index> particle_indices(
167 max_particles_per_cell * n_comps);
168
170 max_particles_per_cell * n_comps, fe.n_dofs_per_cell());
171
172 auto particle = particle_handler.begin();
173 while (particle != particle_handler.end())
174 {
175 const auto &cell = particle->get_surrounding_cell();
176 const auto &dh_cell =
177 typename DoFHandler<dim, spacedim>::cell_iterator(*cell, &space_dh);
178 dh_cell->get_dof_indices(dof_indices);
179 const auto pic = particle_handler.particles_in_cell(cell);
180 const auto n_particles = particle_handler.n_particles_in_cell(cell);
181 particle_indices.resize(n_particles * n_comps);
182 local_matrix.reinit({n_particles * n_comps, fe.n_dofs_per_cell()});
183 Assert(pic.begin() == particle, ExcInternalError());
184 for (unsigned int i = 0; particle != pic.end(); ++particle, ++i)
185 {
186 const auto &reference_location =
187 particle->get_reference_location();
188
189 for (unsigned int d = 0; d < n_comps; ++d)
190 particle_indices[i * n_comps + d] =
191 particle->get_id() * n_comps + d;
192
193 for (unsigned int j = 0; j < fe.n_dofs_per_cell(); ++j)
194 {
195 const auto comp_j =
196 space_gtl[fe.system_to_component_index(j).first];
197 if (comp_j != numbers::invalid_unsigned_int)
198 local_matrix(i * n_comps + comp_j, j) =
199 fe.shape_value(j, reference_location);
200 }
201 }
202 constraints.distribute_local_to_global(local_matrix,
203 particle_indices,
204 dof_indices,
205 matrix);
206 }
207 matrix.compress(VectorOperation::add);
208 }
209
210#include "utilities.inst"
211
212 } // namespace Utilities
213} // namespace Particles
void distribute_local_to_global(const InVector &local_vector, const std::vector< size_type > &local_dof_indices, OutVector &global_vector) const
void add_entries_local_to_global(const std::vector< size_type > &local_dof_indices, SparsityPatternBase &sparsity_pattern, const bool keep_constrained_entries=true, const Table< 2, bool > &dof_mask=Table< 2, bool >()) const
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 types::fe_index index=0) const
types::global_dof_index n_dofs() const
types::particle_index n_global_particles() const
particle_iterator begin() const
particle_iterator end() const
types::particle_index n_locally_owned_particles() const
types::particle_index n_particles_in_cell(const typename Triangulation< dim, spacedim >::active_cell_iterator &cell) const
particle_iterator_range particles_in_cell(const typename Triangulation< dim, spacedim >::active_cell_iterator &cell)
types::particle_index n_global_max_particles_per_cell() const
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
#define Assert(cond, exc)
#define AssertDimension(dim1, dim2)
static ::ExceptionBase & ExcInternalError()
typename ActiveSelector::cell_iterator cell_iterator
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:114
void create_interpolation_sparsity_pattern(const DoFHandler< dim, spacedim > &space_dh, const Particles::ParticleHandler< dim, spacedim > &particle_handler, SparsityPatternBase &sparsity, const AffineConstraints< number > &constraints=AffineConstraints< number >(), const ComponentMask &space_comps=ComponentMask())
Definition utilities.cc:32
static const unsigned int invalid_unsigned_int
Definition types.h:213