SemiDiscreteOT 1.0
Semi-Discrete Optimal Transport Library
Loading...
Searching...
No Matches
PowerDiagram.h
Go to the documentation of this file.
1#ifndef POWER_DIAGRAM_H
2#define POWER_DIAGRAM_H
3
4#include <deal.II/base/point.h>
5#include <deal.II/grid/tria.h>
6#include <deal.II/grid/tria_accessor.h>
7#include <deal.II/grid/tria_iterator.h>
8#include <deal.II/base/geometry_info.h>
9#include <deal.II/numerics/data_out.h>
10#include <vector>
11#include <string>
12#include <deal.II/lac/vector.h>
13#include <memory>
14#include <geogram/basic/smart_pointer.h>
15#include <geogram/delaunay/delaunay.h>
16#include <geogram/voronoi/RVD.h>
17
19
20// Forward declarations for Geogram types
21namespace GEO {
22 class Mesh;
23 class MeshCellsAABB;
24}
25
27
28using namespace dealii;
29
39template <int dim, int spacedim=dim>
41public:
42 virtual ~PowerDiagramBase() = default;
43
49 virtual void set_generators(const std::vector<Point<spacedim>> &points,
50 const Vector<double> &potentials) = 0;
51
55 virtual void compute_power_diagram() = 0;
56
61 virtual void output_vtu(const std::string& filename) const = 0;
62
66 virtual void compute_cell_centroids() = 0;
71 virtual void save_centroids_to_file(const std::string& filename) const = 0;
75 virtual const std::vector<Point<spacedim>>& get_cell_centroids() const = 0;
76
77protected:
78 std::vector<Point<spacedim>> generator_points;
79 std::vector<double> generator_potentials;
80 std::vector<Point<spacedim>> cell_centroids;
81};
82
89template <int dim, int spacedim=dim>
90class DealIIPowerDiagram : public PowerDiagramBase<dim, spacedim> {
91public:
96 DealIIPowerDiagram(const Triangulation<dim, spacedim> &source_mesh);
97
103 void set_generators(const std::vector<Point<spacedim>> &points,
104 const Vector<double> &potentials) override;
105
111 const std::function<double(const Point<spacedim>&, const Point<spacedim>&)>& dist)
112 {
113 distance_function = dist;
114 }
115
122 double power_distance(const Point<spacedim> &point,
123 const unsigned int generator_idx) const;
124
128 void compute_power_diagram() override;
129
134 void output_vtu(const std::string& filename) const override;
135
141 unsigned int get_cell_assignment(const unsigned int cell_index) const;
145 const std::vector<unsigned int>& get_cell_assignments() const;
146
150 void compute_cell_centroids() override;
155 void save_centroids_to_file(const std::string& filename) const override;
159 const std::vector<Point<spacedim>>& get_cell_centroids() const override;
160
161private:
162 const Triangulation<dim, spacedim>* source_triangulation;
163 std::vector<unsigned int> cell_assignments;
164 std::function<double(const Point<spacedim>&, const Point<spacedim>&)> distance_function;
165};
166
167// TODO enable if dim=3 ?
174template <int dim, int spacedim = dim>
175class GeogramPowerDiagram : public PowerDiagramBase<dim, spacedim> {
176public:
181 GeogramPowerDiagram(const std::string& source_mesh_file);
183
189 void set_generators(const std::vector<Point<spacedim>> &points,
190 const Vector<double> &potentials) override;
191
195 void compute_power_diagram() override;
196
201 void output_vtu(const std::string& filename) const override;
202
206 void compute_cell_centroids() override;
211 void save_centroids_to_file(const std::string& filename) const override;
215 const std::vector<Point<spacedim>>& get_cell_centroids() const override;
216
217private:
218 void init_power_diagram();
219 bool load_volume_mesh(const std::string& filename, GEO::Mesh& mesh);
220
221 std::unique_ptr<GEO::Mesh> source_mesh;
223 std::unique_ptr<GEO::Mesh> RVD_mesh;
224 GEO::Delaunay_var delaunay;
225 GEO::RestrictedVoronoiDiagram_var RVD;
226 std::vector<double> points_mesh_target;
227 std::vector<double> points_mesh_target_lifted;
228 bool save_RVD = false;
229 bool save_morph = false;
230};
231
242template <int dim, int spacedim>
243std::unique_ptr<PowerDiagramBase<dim, spacedim>> create_power_diagram(
244 const std::string& implementation_type,
245 const Triangulation<dim, spacedim>* dealii_mesh = nullptr,
246 const std::string& geogram_mesh_file = "") {
247
248 if (implementation_type == "dealii" && dealii_mesh != nullptr) {
249 return std::make_unique<DealIIPowerDiagram<dim, spacedim>>(*dealii_mesh);
250 } else if (implementation_type == "geogram" && !geogram_mesh_file.empty()) {
251 if constexpr (dim==3 && dim==spacedim)
252 return std::make_unique<GeogramPowerDiagram<dim, spacedim>>(geogram_mesh_file);
253 else
254 throw std::runtime_error("Geogram power diagram is only available for dim==spacedim");
255 }
256 throw std::runtime_error("Invalid power diagram implementation type or missing required parameters");
257}
258
259} // namespace PowerDiagramSpace
260
261#endif
A class for computing power diagrams using deal.II.
void output_vtu(const std::string &filename) const override
Outputs the power diagram to a VTU file.
unsigned int get_cell_assignment(const unsigned int cell_index) const
Returns the assignment of a cell to a generator.
const Triangulation< dim, spacedim > * source_triangulation
The source triangulation.
void compute_cell_centroids() override
Computes the centroids of the power cells.
void set_generators(const std::vector< Point< spacedim > > &points, const Vector< double > &potentials) override
Sets the generator points and their potentials.
void set_distance_function(const std::function< double(const Point< spacedim > &, const Point< spacedim > &)> &dist)
Sets the distance function to be used.
std::vector< unsigned int > cell_assignments
The assignments of cells to generators.
void save_centroids_to_file(const std::string &filename) const override
Saves the centroids of the power cells to a file.
std::function< double(const Point< spacedim > &, const Point< spacedim > &)> distance_function
The distance function.
double power_distance(const Point< spacedim > &point, const unsigned int generator_idx) const
Computes the power distance between a point and a generator.
const std::vector< unsigned int > & get_cell_assignments() const
Returns the assignments of all cells to generators.
void compute_power_diagram() override
Computes the power diagram.
const std::vector< Point< spacedim > > & get_cell_centroids() const override
Returns the centroids of the power cells.
A class for computing power diagrams using Geogram.
bool save_morph
A flag to indicate whether to save the morphed mesh.
void output_vtu(const std::string &filename) const override
Outputs the power diagram to a VTU file.
const std::vector< Point< spacedim > > & get_cell_centroids() const override
Returns the centroids of the power cells.
std::vector< double > points_mesh_target_lifted
The lifted target points.
void compute_cell_centroids() override
Computes the centroids of the power cells.
void compute_power_diagram() override
Computes the power diagram.
void set_generators(const std::vector< Point< spacedim > > &points, const Vector< double > &potentials) override
Sets the generator points and their potentials.
std::vector< double > points_mesh_target
The target points.
void save_centroids_to_file(const std::string &filename) const override
Saves the centroids of the power cells to a file.
bool save_RVD
A flag to indicate whether to save the restricted Voronoi diagram.
std::unique_ptr< GEO::Mesh > source_mesh
The source mesh.
GEO::RestrictedVoronoiDiagram_var RVD
The restricted Voronoi diagram.
bool load_volume_mesh(const std::string &filename, GEO::Mesh &mesh)
int dimension_voronoi
The dimension of the Voronoi diagram.
std::unique_ptr< GEO::Mesh > RVD_mesh
The restricted Voronoi diagram mesh.
GEO::Delaunay_var delaunay
The Delaunay triangulation.
A base class for power diagrams.
std::vector< double > generator_potentials
The potentials of the generator points.
virtual void compute_power_diagram()=0
Computes the power diagram.
virtual const std::vector< Point< spacedim > > & get_cell_centroids() const =0
Returns the centroids of the power cells.
virtual void output_vtu(const std::string &filename) const =0
Outputs the power diagram to a VTU file.
std::vector< Point< spacedim > > generator_points
The generator points.
virtual void set_generators(const std::vector< Point< spacedim > > &points, const Vector< double > &potentials)=0
Sets the generator points and their potentials.
virtual void compute_cell_centroids()=0
Computes the centroids of the power cells.
virtual void save_centroids_to_file(const std::string &filename) const =0
Saves the centroids of the power cells to a file.
std::vector< Point< spacedim > > cell_centroids
The centroids of the power cells.
std::unique_ptr< PowerDiagramBase< dim, spacedim > > create_power_diagram(const std::string &implementation_type, const Triangulation< dim, spacedim > *dealii_mesh=nullptr, const std::string &geogram_mesh_file="")
A factory function to create a power diagram.