22#ifndef INCLUDE_EXADG_POSTPROCESSOR_WRITE_OUTPUT_H_
23#define INCLUDE_EXADG_POSTPROCESSOR_WRITE_OUTPUT_H_
29#include <deal.II/base/bounding_box.h>
30#include <deal.II/grid/grid_generator.h>
31#include <deal.II/grid/grid_out.h>
32#include <deal.II/numerics/data_out.h>
33#include <deal.II/numerics/data_out_faces.h>
34#include <deal.II/particles/data_out.h>
35#include <deal.II/particles/particle_handler.h>
41write_surface_mesh(dealii::Triangulation<dim>
const & triangulation,
42 dealii::Mapping<dim>
const & mapping,
43 unsigned int const n_subdivisions,
44 std::string
const & folder,
45 std::string
const & file,
46 unsigned int const counter,
47 MPI_Comm
const & mpi_comm)
50 dealii::DataOutFaces<dim> data_out_surface(
true );
51 data_out_surface.attach_triangulation(triangulation);
52 data_out_surface.build_patches(mapping, n_subdivisions);
53 data_out_surface.write_vtu_with_pvtu_record(folder, file +
"_surface", counter, mpi_comm, 4);
58write_boundary_IDs(dealii::Triangulation<dim>
const & triangulation,
59 std::string
const & folder,
60 std::string
const & file,
61 MPI_Comm
const & mpi_communicator)
63 unsigned int const rank = dealii::Utilities::MPI::this_mpi_process(mpi_communicator);
64 unsigned int const n_ranks = dealii::Utilities::MPI::n_mpi_processes(mpi_communicator);
66 unsigned int const n_digits =
static_cast<int>(std::ceil(std::log10(std::fabs(n_ranks))));
68 std::string filename = folder + file +
"_boundary_IDs" +
"." +
69 dealii::Utilities::int_to_string(rank, n_digits) +
".vtk";
70 std::ofstream output(filename);
72 dealii::GridOut grid_out;
73 dealii::GridOutFlags::Vtk flags;
74 flags.output_cells =
false;
75 flags.output_faces =
true;
76 flags.output_edges =
false;
77 flags.output_only_relevant =
false;
78 grid_out.set_flags(flags);
79 grid_out.write_vtk(triangulation, output);
84write_grid(dealii::Triangulation<dim>
const & triangulation,
85 dealii::Mapping<dim>
const & mapping,
86 unsigned int const n_subdivisions,
87 std::string
const & folder,
88 std::string
const & file,
89 unsigned int const & counter,
90 MPI_Comm
const & mpi_comm)
92 std::string filename = file +
"_grid";
94 dealii::DataOut<dim> data_out;
96 dealii::DataOutBase::VtkFlags flags;
97 flags.write_higher_order_cells = n_subdivisions > 1;
98 data_out.set_flags(flags);
100 data_out.attach_triangulation(triangulation);
101 data_out.build_patches(mapping, n_subdivisions, dealii::DataOut<dim>::curved_inner_cells);
102 data_out.write_vtu_with_pvtu_record(folder, filename, counter, mpi_comm, 4);
107write_points(dealii::Triangulation<dim>
const & triangulation,
108 dealii::Mapping<dim>
const & mapping,
109 std::vector<dealii::Point<dim>>
const & points,
110 std::string
const & folder,
111 std::string
const & file,
112 unsigned int const counter,
113 MPI_Comm
const & mpi_comm)
115 std::string filename = file +
"_points";
117 dealii::Particles::ParticleHandler<dim, dim> particle_handler(triangulation, mapping);
119 particle_handler.insert_particles(points);
121 dealii::Particles::DataOut<dim, dim> particle_output;
122 particle_output.build_patches(particle_handler);
123 particle_output.write_vtu_with_pvtu_record(folder, filename, counter, mpi_comm);
128write_points_in_dummy_triangulation(std::vector<dealii::Point<dim>>
const & points,
129 std::string
const & folder,
130 std::string
const & file,
131 unsigned int const counter,
132 MPI_Comm
const & mpi_comm)
134 dealii::BoundingBox<dim> bounding_box(points);
135 auto const boundary_points =
136 bounding_box.create_extended(1e-3 * std::pow(bounding_box.volume(), 1.0 / ((
double)dim)))
137 .get_boundary_points();
139 dealii::Triangulation<dim> particle_dummy_tria;
140 dealii::GridGenerator::hyper_rectangle(particle_dummy_tria,
141 boundary_points.first,
142 boundary_points.second);
144 dealii::MappingQGeneric<dim> particle_dummy_mapping(1 );
147 particle_dummy_tria, particle_dummy_mapping, points, folder, file, counter, mpi_comm);