23#include <deal.II/numerics/data_out_dof_data.templates.h>
32template <
int dim,
int patch_dim,
int spacedim>
42template <
int dim,
int patch_dim,
int spacedim>
46 const unsigned int n_subdivisions)
55 std::vector<std::pair<types::global_dof_index, Point<spacedim>>> points_all;
68 partitioner = std::make_shared<Utilities::MPI::Partitioner>(
76 cell->get_dof_indices(dof_indices);
78 for (
unsigned int i = 0; i < dof_indices.size(); ++i)
79 points_all.emplace_back(
partitioner->global_to_local(dof_indices[i]),
83 std::sort(points_all.begin(),
85 [](
const auto &a,
const auto &b) { return a.first < b.first; });
86 points_all.erase(std::unique(points_all.begin(),
88 [](
const auto &a,
const auto &b) {
89 return a.first == b.first;
93 std::vector<Point<spacedim>> points;
95 for (
const auto &i : points_all)
98 points.push_back(i.second);
106template <
int dim,
int patch_dim,
int spacedim>
110 const unsigned int n_subdivisions,
119template <
int dim,
int patch_dim,
int spacedim>
126 if (
rpe.is_ready() ==
false)
131 "Mapping is not valid anymore! Please register a new mapping via "
132 "update_mapping() or the other build_patches() function."));
136 std::vector<std::shared_ptr<LinearAlgebra::distributed::Vector<double>>>
141 unsigned int counter = 0;
143 for (
const auto &data : this->
dof_data)
145 const auto data_ptr =
dynamic_cast<
146 internal::DataOutImplementation::DataEntry<dim, spacedim, double> *
>(
151 const auto &dh = *data_ptr->dof_handler;
154 for (
const auto &fe : dh.get_fe_collection())
156 fe.n_base_elements() == 1,
158 "This class currently only supports scalar elements and elements "
159 "with a single base element."));
162 for (
unsigned int comp = 0; comp < dh.get_fe_collection().n_components();
168 vectors.emplace_back(
172 for (
unsigned int j = 0; j < values.size(); ++j)
176 vectors.back()->set_ghost_state(
true);
183 std::string(
"temp_" + std::to_string(counter)),
185 DataVectorType::type_dof_data);
198template <
int dim,
int patch_dim,
int spacedim>
199const std::vector<typename DataOutBase::Patch<patch_dim, spacedim>> &
207#include "data_out_resample.inst"
virtual const std::vector< typename DataOutBase::Patch< patch_dim, spacedim > > & get_patches() const override
std::shared_ptr< Utilities::MPI::Partitioner > partitioner
DataOutResample(const Triangulation< patch_dim, spacedim > &patch_tria, const Mapping< patch_dim, spacedim > &patch_mapping)
DoFHandler< patch_dim, spacedim > patch_dof_handler
std::vector< types::global_dof_index > point_to_local_vector_indices
SmartPointer< const Mapping< dim, spacedim > > 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)
Utilities::MPI::RemotePointEvaluation< dim, spacedim > rpe
void update_mapping(const Mapping< dim, spacedim > &mapping, const unsigned int n_subdivisions=0)
const SmartPointer< const Mapping< patch_dim, spacedim > > patch_mapping
DataOut< patch_dim, spacedim > patch_data_out
std::vector< std::shared_ptr< internal::DataOutImplementation::DataEntryBase< dim, spacedim > > > dof_data
friend class DataOut_DoFData
SmartPointer< const Triangulation< dim, spacedim > > triangulation
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.
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
static ::ExceptionBase & ExcNotImplemented()
#define Assert(cond, exc)
static ::ExceptionBase & ExcMessage(std::string arg1)
@ update_quadrature_points
Transformed quadrature points.
::VectorizedArray< Number, width > max(const ::VectorizedArray< Number, width > &, const ::VectorizedArray< Number, width > &)