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
data_out_resample.cc
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2021 - 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
17
19#include <deal.II/fe/fe_tools.h>
20
23
24#include <deal.II/numerics/data_out_dof_data.templates.h>
27
28#include <sstream>
29
31
32
33template <int dim, int patch_dim, int spacedim>
35 const Triangulation<patch_dim, spacedim> &patch_tria,
36 const Mapping<patch_dim, spacedim> & patch_mapping)
37 : patch_dof_handler(patch_tria)
38 , patch_mapping(&patch_mapping)
39{}
40
41
42
43template <int dim, int patch_dim, int spacedim>
44void
46 const Mapping<dim, spacedim> &mapping,
47 const unsigned int n_subdivisions)
48{
49 this->mapping = &mapping;
50 this->point_to_local_vector_indices.clear();
51
53 std::max<unsigned int>(1, n_subdivisions));
54 patch_dof_handler.distribute_dofs(fe);
55
56 std::vector<std::pair<types::global_dof_index, Point<spacedim>>> points_all;
57
59
60 FEValues<patch_dim, spacedim> fe_values(*patch_mapping,
61 fe,
62 quadrature,
64
65 std::vector<types::global_dof_index> dof_indices(fe.n_dofs_per_cell());
66
67 IndexSet active_dofs;
68 DoFTools::extract_locally_active_dofs(patch_dof_handler, active_dofs);
69 partitioner = std::make_shared<Utilities::MPI::Partitioner>(
70 patch_dof_handler.locally_owned_dofs(), active_dofs, MPI_COMM_WORLD);
71
72 for (const auto &cell : patch_dof_handler.active_cell_iterators() |
74 {
75 fe_values.reinit(cell);
76
77 cell->get_dof_indices(dof_indices);
78
79 for (unsigned int i = 0; i < dof_indices.size(); ++i)
80 points_all.emplace_back(partitioner->global_to_local(dof_indices[i]),
81 fe_values.quadrature_point(i));
82 }
83
84 std::sort(points_all.begin(),
85 points_all.end(),
86 [](const auto &a, const auto &b) { return a.first < b.first; });
87 points_all.erase(std::unique(points_all.begin(),
88 points_all.end(),
89 [](const auto &a, const auto &b) {
90 return a.first == b.first;
91 }),
92 points_all.end());
93
94 std::vector<Point<spacedim>> points;
95
96 for (const auto &i : points_all)
97 {
98 point_to_local_vector_indices.push_back(i.first);
99 points.push_back(i.second);
100 }
101
102 rpe.reinit(points, *this->triangulation, *this->mapping);
103}
104
105
106
107template <int dim, int patch_dim, int spacedim>
108void
110 const Mapping<dim, spacedim> & mapping,
111 const unsigned int n_subdivisions,
112 const typename DataOut<patch_dim, spacedim>::CurvedCellRegion curved_region)
113{
114 this->update_mapping(mapping, n_subdivisions);
115 this->build_patches(curved_region);
116}
117
118
119
120template <int dim, int patch_dim, int spacedim>
121void
123 const typename DataOut<patch_dim, spacedim>::CurvedCellRegion curved_region)
124{
125 patch_data_out.clear();
126
127 if (rpe.is_ready() == false)
128 {
129 Assert(
130 this->mapping,
132 "Mapping is not valid anymore! Please register a new mapping via "
133 "update_mapping() or the other build_patches() function."));
134 update_mapping(*this->mapping, patch_dof_handler.get_fe().degree);
135 }
136
137 std::vector<std::shared_ptr<LinearAlgebra::distributed::Vector<double>>>
138 vectors;
139
140 patch_data_out.attach_dof_handler(patch_dof_handler);
141
142 unsigned int counter = 0;
143
144 for (const auto &data : this->dof_data)
145 {
146 const auto data_ptr = dynamic_cast<
147 internal::DataOutImplementation::DataEntry<dim, spacedim, double> *>(
148 data.get());
149
150 Assert(data_ptr, ExcNotImplemented());
151
152 const auto &dh = *data_ptr->dof_handler;
153
154#ifdef DEBUG
155 for (const auto &fe : dh.get_fe_collection())
156 Assert(
157 fe.n_base_elements() == 1,
159 "This class currently only supports scalar elements and elements "
160 "with a single base element."));
161#endif
162
163 for (unsigned int comp = 0; comp < dh.get_fe_collection().n_components();
164 ++comp)
165 {
166 const auto values = VectorTools::point_values<1>(
167 rpe, dh, data_ptr->vector, VectorTools::EvaluationFlags::avg, comp);
168
169 vectors.emplace_back(
171 partitioner));
172
173 for (unsigned int j = 0; j < values.size(); ++j)
174 vectors.back()->local_element(point_to_local_vector_indices[j]) =
175 values[j];
176
177 vectors.back()->set_ghost_state(true);
178
179 // we can give the vectors arbitrary names ("temp_*") here, since
180 // these are only used internally (by patch_data_out) but not later on
181 // during the actual output to file
182 patch_data_out.add_data_vector(
183 *vectors.back(),
184 std::string("temp_" + std::to_string(counter)),
186 DataVectorType::type_dof_data);
187
188 counter++;
189 }
190 }
191
192 patch_data_out.build_patches(*patch_mapping,
193 patch_dof_handler.get_fe().degree,
194 curved_region);
195}
196
197
198
199template <int dim, int patch_dim, int spacedim>
200const std::vector<typename DataOutBase::Patch<patch_dim, spacedim>> &
202{
203 return patch_data_out.get_patches();
204}
205
206
207// explicit instantiations
208#include "data_out_resample.inst"
209
virtual const std::vector< typename DataOutBase::Patch< patch_dim, spacedim > > & get_patches() const override
DataOutResample(const Triangulation< patch_dim, spacedim > &patch_tria, const Mapping< patch_dim, spacedim > &patch_mapping)
void build_patches(const Mapping< dim, spacedim > &mapping, const unsigned int n_subdivisions=0, const typename DataOut< patch_dim, spacedim >::CurvedCellRegion curved_region=DataOut< patch_dim, spacedim >::CurvedCellRegion::curved_boundary)
void update_mapping(const Mapping< dim, spacedim > &mapping, const unsigned int n_subdivisions=0)
CurvedCellRegion
Definition data_out.h:186
const Point< spacedim > & quadrature_point(const unsigned int q_point) const
void reinit(const TriaIterator< DoFCellAccessor< dim, spacedim, level_dof_access > > &cell)
unsigned int n_dofs_per_cell() const
const std::vector< Point< dim > > & get_unit_support_points() const
Abstract base class for mapping classes.
Definition mapping.h:317
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
static ::ExceptionBase & ExcNotImplemented()
#define Assert(cond, exc)
static ::ExceptionBase & ExcMessage(std::string arg1)
@ update_mapping
@ update_quadrature_points
Transformed quadrature points.
IndexSet extract_locally_active_dofs(const DoFHandler< dim, spacedim > &dof_handler)
const ::parallel::distributed::Triangulation< dim, spacedim > * triangulation