Parabola length calculation.

src/libslic3r/SLA/SupportIslands/Parabola.hpp

@@ -0,0 +1,25 @@
+#ifndef slic3r_SLA_SuppotstIslands_Parabola_hpp_
+#define slic3r_SLA_SuppotstIslands_Parabola_hpp_
+
+#include <libslic3r/Line.hpp>
+#include <libslic3r/Point.hpp>
+
+namespace Slic3r::sla {
+
+/// <summary>
+/// DTO store prabola params
+/// A parabola can be defined geometrically as a set of points (locus of points) in the Euclidean plane:
+/// Where distance from focus point is same as distance from line(directrix).
+/// </summary>
+struct Parabola
+{
+    Line  directrix;
+    Point focus;
+
+    Parabola(Line directrix, Point focus)
+        : directrix(std::move(directrix)), focus(std::move(focus))
+    {}
+};
+
+} // namespace Slic3r::sla
+#endif // slic3r_SLA_SuppotstIslands_Parabola_hpp_

src/libslic3r/SLA/SupportIslands/ParabolaUtils.cpp

@@ -0,0 +1,76 @@
+#include "ParabolaUtils.hpp"
+
+using namespace Slic3r::sla;
+
+double ParabolaUtils::calculate_length_of_parabola(
+    const Parabola &parabola, const Point &from, const Point &to)
+{
+    const Point &point = parabola.focus;
+    const Line & line  = parabola.directrix;
+    Line norm_line(point, point + line.normal());
+
+    // sign of distance is resolved by dot product in function is_over_zero()
+    double scaled_x1 = norm_line.perp_distance_to(from);
+    double scaled_x2 = norm_line.perp_distance_to(to);
+
+    double parabola_scale = 1. / (4. * focal_length(parabola));
+
+    double x1 = scaled_x1 * parabola_scale;
+    double x2 = scaled_x2 * parabola_scale;
+
+    double length_x1 = parabola_arc_length(x1) / parabola_scale;
+    double length_x2 = parabola_arc_length(x2) / parabola_scale;
+
+    return (is_over_zero(parabola, from, to)) ?
+               (length_x1 + length_x2) :    // interval is over zero
+               fabs(length_x1 - length_x2); // interval is on same side of parabola
+}
+
+
+#include <Libslic3r/Geometry.hpp>
+#include <Libslic3r/VoronoiOffset.hpp>
+#include <Libslic3r/VoronoiVisualUtils.hpp>
+double ParabolaUtils::calculate_length_of_parabola_by_sampling(
+    const Parabola &parabola,
+    const Point &   from,
+    const Point &   to,
+    double          discretization_step)
+{
+    using VD = Slic3r::Geometry::VoronoiDiagram;
+    std::vector<Voronoi::Internal::point_type> parabola_samples({from, to});
+
+    VD::point_type   source_point = parabola.focus;
+    VD::segment_type source_segment(parabola.directrix.a, parabola.directrix.b);
+    ::boost::polygon::voronoi_visual_utils<double>::discretize(
+        source_point, source_segment, discretization_step, &parabola_samples);
+
+    double sumLength = 0;
+    for (size_t index = 1; index < parabola_samples.size(); ++index) {
+        double diffX = parabola_samples[index - 1].x() -
+                       parabola_samples[index].x();
+        double diffY = parabola_samples[index - 1].y() -
+                       parabola_samples[index].y();
+        double length = sqrt(diffX * diffX + diffY * diffY);
+        sumLength += length;
+    }
+    return sumLength;
+}
+
+double ParabolaUtils::focal_length(const Parabola &parabola)
+{
+    // https://en.wikipedia.org/wiki/Parabola
+    // p = 2f; y = 1/(4f) * x^2; y = 1/(2p) * x^2
+    double p = parabola.directrix.perp_distance_to(parabola.focus);
+    double f = p / 2.;
+    return f;
+}
+
+bool ParabolaUtils::is_over_zero(const Parabola &parabola,
+                                 const Point &   from,
+                                 const Point &   to)
+{
+    Point line_direction = parabola.directrix.b - parabola.directrix.a;
+    bool  is_positive_x1 = line_direction.dot(parabola.focus - from) > 0.;
+    bool  is_positive_x2 = line_direction.dot(parabola.focus - to) > 0.;
+    return is_positive_x1 != is_positive_x2;
+}

src/libslic3r/SLA/SupportIslands/ParabolaUtils.hpp

@@ -0,0 +1,74 @@
+#ifndef slic3r_SLA_SuppotstIslands_ParabolaUtils_hpp_
+#define slic3r_SLA_SuppotstIslands_ParabolaUtils_hpp_
+
+#include "Parabola.hpp"
+
+namespace Slic3r::sla {
+
+/// <summary>
+/// Class which contain collection of static function
+/// for work with Parabola.
+/// </summary>
+class ParabolaUtils
+{
+public:
+    ParabolaUtils() = delete;
+
+    /// <summary>
+    /// Integrate length over interval defined by points from and to
+    /// </summary>
+    /// <param name="parabola">Input parabola</param>
+    /// <param name="from">Input point lay on parabola</param>
+    /// <param name="to">Input point lay on parabola</param>
+    /// <returns>Length of parabola arc</returns>
+    static double calculate_length_of_parabola(const Parabola &parabola,
+                                               const Point &   from,
+                                               const Point &   to);
+
+    /// <summary>
+    /// Sample parabola between points from and to by step.
+    /// </summary>
+    /// <param name="parabola">Input parabola</param>
+    /// <param name="from">Input point lay on parabola</param>
+    /// <param name="to">Input point lay on parabola</param>
+    /// <param name="discretization_step">Define sampling</param>
+    /// <returns>Length of parabola arc</returns>
+    static double calculate_length_of_parabola_by_sampling(
+        const Parabola &parabola,
+        const Point &   from,
+        const Point &   to,
+        double          discretization_step = 0.0002 * 1e6);
+
+    /// <summary>
+    /// calculate focal length of parabola
+    /// </summary>
+    /// <param name="parabola">input parabola</param>
+    /// <returns>Focal length</returns>
+    static double focal_length(const Parabola &parabola);
+
+    /// <summary>
+    /// Check if parabola interval (from, to) contains top of parabola
+    /// </summary>
+    /// <param name="parabola">input parabola</param>
+    /// <param name="from">Start of interval, point lay on parabola</param>
+    /// <param name="to">End of interval, point lay on parabola</param>
+    /// <returns>True when interval contain top of parabola otherwise False</returns>
+    static bool is_over_zero(const Parabola &parabola,
+                             const Point &   from,
+                             const Point &   to);
+
+private:
+    /// <summary>
+    /// Integral of parabola: y = x^2 from zero to point x
+    /// https://ocw.mit.edu/courses/mathematics/18-01sc-single-variable-calculus-fall-2010/unit-4-techniques-of-integration/part-b-partial-fractions-integration-by-parts-arc-length-and-surface-area/session-78-computing-the-length-of-a-curve/MIT18_01SCF10_Ses78d.pdf
+    /// </summary>
+    /// <param name="x">x coordinate of parabola, Positive number</param>
+    /// <returns>Length of parabola from zero to x</returns> 
+    static double parabola_arc_length(double x) {
+        double sqrtRes = sqrt(1 + 4 * x * x);
+        return 1 / 4. * log(2 * x + sqrtRes) + 1 / 2. * x * sqrtRes;
+    };
+};
+
+} // namespace Slic3r::sla
+#endif // slic3r_SLA_SuppotstIslands_ParabolaUtils_hpp_

src/libslic3r/SLA/SupportIslands/SampleConfig.hpp

@@ -17,11 +17,10 @@ struct SampleConfig
 
     // maximal length of longest path in voronoi diagram to be island
     // supported only by one single support point this point will be in center
-    // of path suggestion: smaller than 2* SampleConfig.start_distance
+    // of path.
+    // suggestion 1: Smaller than 2 * SampleConfig.start_distance
+    // suggestion 2: Bigger than 2 * Head diameter
     double max_length_for_one_support_point = 1.;
-
-    // each curve is sampled by this value to test distance to edge of island
-    double curve_sample = 1.; // must be bigger than zero
 };
 
 } // namespace Slic3r::sla

src/libslic3r/SLA/SupportIslands/VoronoiGraph.hpp

@@ -52,13 +52,18 @@ struct VoronoiGraph::Node
 struct VoronoiGraph::Node::Neighbor
 {
     const VD::edge_type *edge;
+    // pointer on graph node structure
+    const Node *node;
+
     // length edge between vertices
     double edge_length;
-    // pointer on graph node structure
-    Node *node;
+
+    // maximal widht of sland(distance to outline)
+    double max_width;
+
 public:
-    Neighbor(const VD::edge_type *edge, double edge_length, Node *node)
-        : edge(edge), edge_length(edge_length), node(node)
+    Neighbor(const VD::edge_type *edge, const Node *node, double edge_length)
+        : edge(edge), node(node), edge_length(edge_length)
     {}
 };
 

src/libslic3r/SLA/SupportIslands/VoronoiGraphUtils.cpp

@@ -1,17 +1,19 @@
 #include "VoronoiGraphUtils.hpp"
 
+#include <cmath>
 #include <libslic3r/VoronoiOffset.hpp>
 #include "IStackFunction.hpp"
 #include "EvaluateNeighbor.hpp"
+#include "ParabolaUtils.hpp"
 
 #include <libslic3r/VoronoiVisualUtils.hpp>
 
 using namespace Slic3r::sla;
 
 VoronoiGraph::Node *VoronoiGraphUtils::getNode(VoronoiGraph &         graph,
-                            const VD::vertex_type *vertex,
-                            const VD::edge_type *  edge,
-                            const Lines &          lines)
+                                               const VD::vertex_type *vertex,
+                                               const VD::edge_type *  edge,
+                                               const Lines &          lines)
 {
     std::map<const VD::vertex_type *, VoronoiGraph::Node> &data = graph.data;
     auto &mapItem = data.find(vertex);
@@ -34,31 +36,104 @@ VoronoiGraph::Node *VoronoiGraphUtils::getNode(VoronoiGraph &         graph,
     return &iterator->second;
 }
 
-double VoronoiGraphUtils::calculate_length_of_parabola(
-    const VD::edge_type &                             edge,
-    const std::vector<Voronoi::Internal::segment_type> &segments)
+Slic3r::Point VoronoiGraphUtils::retrieve_point(const Lines &        lines,
+                                                const VD::cell_type &cell)
 {
-    // TODO: len by param not sampling of parabola
+    assert(cell.source_category() ==
+               boost::polygon::SOURCE_CATEGORY_SEGMENT_START_POINT ||
+           cell.source_category() ==
+               boost::polygon::SOURCE_CATEGORY_SEGMENT_END_POINT);
+    return (cell.source_category() ==
+            boost::polygon::SOURCE_CATEGORY_SEGMENT_START_POINT) ?
+               lines[cell.source_index()].a :
+               lines[cell.source_index()].b;
+}
 
-    double discretization_step = 0.0002 * 1e7;
-    Points points; // voronoi created by line segments only
+Slic3r::Point VoronoiGraphUtils::get_parabola_point(
+    const VD::edge_type &parabola, const Slic3r::Lines &lines)
+{
+    using namespace boost::polygon;
+    assert(parabola.is_curved());
+    const VD::cell_type& cell = (parabola.cell()->contains_point())?
+                                *parabola.cell() : *parabola.twin()->cell();
+    assert(cell.contains_point());
+    assert(cell.source_category() == SOURCE_CATEGORY_SEGMENT_START_POINT ||
+           cell.source_category() == SOURCE_CATEGORY_SEGMENT_END_POINT);
+    return (cell.source_category() == SOURCE_CATEGORY_SEGMENT_START_POINT) ?
+           lines[cell.source_index()].a :
+           lines[cell.source_index()].b;
+}
 
-    std::vector<Voronoi::Internal::point_type> samples;
-    samples.push_back(Voronoi::Internal::point_type(edge.vertex0()->x(),
-                                                    edge.vertex0()->y()));
-    samples.push_back(Voronoi::Internal::point_type(edge.vertex1()->x(),
-                                                    edge.vertex1()->y()));
-    Voronoi::Internal::sample_curved_edge(points, segments, edge, samples,
-                                          discretization_step); 
-    
-    double sumLength = 0;
-    for (size_t index = 1; index < samples.size(); ++index) {
-        double diffX = samples[index - 1].x() - samples[index].x();
-        double diffY = samples[index - 1].y() - samples[index].y();
-        double length = sqrt(diffX * diffX + diffY * diffY);
-        sumLength += length;
+Slic3r::Line VoronoiGraphUtils::get_parabola_line(
+    const VD::edge_type &parabola, const Slic3r::Lines &lines)
+{
+    assert(parabola.is_curved());
+    const VD::cell_type& cell = (parabola.cell()->contains_segment())?
+                                *parabola.cell() : *parabola.twin()->cell();
+    assert(cell.contains_segment());
+    return lines[cell.source_index()];
+}
+
+Parabola VoronoiGraphUtils::get_parabola(
+    const VD::edge_type &edge, const Lines &lines)
+{
+    Point    point = get_parabola_point(edge, lines);
+    Line     line  = get_parabola_line(edge, lines);
+    return Parabola(line, point);
+}
+
+double VoronoiGraphUtils::calculate_length_of_parabola(
+    const VD::edge_type &                               edge,
+    const Lines &                                       lines)
+{
+    Point v0{edge.vertex0()->x(), edge.vertex0()->y()};
+    Point v1{edge.vertex1()->x(), edge.vertex1()->y()};
+    Parabola parabola = get_parabola(edge, lines);
+    return ParabolaUtils::calculate_length_of_parabola(parabola, v0, v1);
+}
+
+double VoronoiGraphUtils::calculate_length(
+    const VD::edge_type &edge, const Lines &lines)
+{
+    if (edge.is_linear()) {
+        const VD::vertex_type* v0 = edge.vertex0();
+        const VD::vertex_type* v1 = edge.vertex1();
+        double diffX = v0->x() - v1->x();
+        double diffY = v0->y() - v1->y();
+        return sqrt(diffX * diffX + diffY * diffY);
+    }    
+    assert(edge.is_curved());
+    return calculate_length_of_parabola(edge, lines);
+}
+
+double VoronoiGraphUtils::calculate_max_width(
+    const VD::edge_type &edge, const Lines &lines)
+{
+    Point v0{edge.vertex0()->x(), edge.vertex0()->y()};
+    Point v1{edge.vertex1()->x(), edge.vertex1()->y()};
+
+    if (edge.is_linear()) {
+        // line is initialized by 2 line segments only
+        assert(!edge.cell()->contains_point());
+        assert(edge.cell()->contains_segment());
+        assert(!edge.twin()->cell()->contains_point());
+        assert(edge.twin()->cell()->contains_segment());
+
+        const Line &line = lines[edge.cell()->source_index()];
+
+        double distance0 = line.perp_distance_to(v0);
+        double distance1 = line.perp_distance_to(v1);
+        return 2 * std::max(distance0, distance1);
     }
-    return sumLength;
+    assert(edge.is_curved());
+    Parabola parabola = get_parabola(edge, lines);
+    // distance to point and line is same
+    // vector from edge vertex to parabola focus point
+    Point  vec0 = parabola.focus - v0;
+    Point  vec1 = parabola.focus - v1;
+    double distance0 = sqrt(vec0.x() * vec0.x() + vec0.y() * vec0.y());
+    double distance1 = sqrt(vec0.x() * vec0.x() + vec0.y() * vec0.y());
+    return 2 * std::max(distance0, distance1);
 }
 
 VoronoiGraph VoronoiGraphUtils::getSkeleton(const VD &vd, const Lines &lines)
@@ -66,13 +141,6 @@ VoronoiGraph VoronoiGraphUtils::getSkeleton(const VD &vd, const Lines &lines)
     // vd should be annotated.
     // assert(Voronoi::debug::verify_inside_outside_annotations(vd));
 
-        std::vector<Voronoi::Internal::segment_type> segments;
-    for (Lines::const_iterator it = lines.begin(); it != lines.end(); ++it)
-        segments.push_back(Voronoi::Internal::segment_type(
-            Voronoi::Internal::point_type(double(it->a(0)), double(it->a(1))),
-            Voronoi::Internal::point_type(double(it->b(0)),
-                                          double(it->b(1)))));
-
     VoronoiGraph         skeleton;
     const VD::edge_type *first_edge = &vd.edges().front();
     for (const VD::edge_type &edge : vd.edges()) {
@@ -102,31 +170,22 @@ VoronoiGraph VoronoiGraphUtils::getSkeleton(const VD &vd, const Lines &lines)
             category1 == Voronoi::VertexCategory::Unknown)
             return {}; // vd must be annotated
 
-        double length = 0;
-        if (edge.is_linear()) {
-            double diffX = v0->x() - v1->x();
-            double diffY = v0->y() - v1->y();
-            length       = sqrt(diffX * diffX + diffY * diffY);
-        } else { // if (edge.is_curved())
-            assert(edge.is_curved());
-            length = calculate_length_of_parabola(edge, segments);
-        }
-
+        double length = calculate_length(edge, lines);
+        double max_width = calculate_max_width(edge, lines);
         VoronoiGraph::Node *node0 = getNode(skeleton, v0, &edge, lines);
         VoronoiGraph::Node *node1 = getNode(skeleton, v1, &edge, lines);
 
         // add extended Edge to graph, both side
-        VoronoiGraph::Node::Neighbor neighbor0(&edge, length, node1);
+        VoronoiGraph::Node::Neighbor neighbor0(&edge, node1, length);
         node0->neighbors.push_back(neighbor0);
-        VoronoiGraph::Node::Neighbor neighbor1(edge.twin(), length, node0);
+        VoronoiGraph::Node::Neighbor neighbor1(edge.twin(), node0, length);
         node1->neighbors.push_back(neighbor1);
     }
     return skeleton;
 }
 
 Slic3r::Point VoronoiGraphUtils::get_offseted_point(
-    const VoronoiGraph::Node &node,
-                                            double padding)
+    const VoronoiGraph::Node &node, double padding)
 {
     assert(node.neighbors.size() == 1);
     const VoronoiGraph::Node::Neighbor &neighbor = node.neighbors.front();
@@ -153,7 +212,7 @@ const VoronoiGraph::Node::Neighbor *VoronoiGraphUtils::get_neighbor(
 }
 
 double VoronoiGraphUtils::get_neighbor_distance(const VoronoiGraph::Node *from,
-                             const VoronoiGraph::Node *to)
+                                                const VoronoiGraph::Node *to)
 {
     const VoronoiGraph::Node::Neighbor *neighbor = get_neighbor(from, to);
     assert(neighbor != nullptr);
@@ -220,7 +279,8 @@ VoronoiGraph::Path VoronoiGraphUtils::find_longest_path_on_circle(
                                               circle_iterator + 1);
     }
     // append longest side branch
-    circle_path.insert(circle_path.end(), longest_circle_branch->nodes.begin(),
+    circle_path.insert(circle_path.end(),
+                       longest_circle_branch->nodes.begin(),
                        longest_circle_branch->nodes.end());
     return {circle_path, longest_branch_length};
 }
@@ -290,7 +350,7 @@ VoronoiGraph::Path VoronoiGraphUtils::find_longest_path_on_circles(
             double length = longest_branch.length + neighbor_path.length;
             if (longest_path.length < length) {
                 longest_path.length = length;
-                longest_path.nodes   = neighbor_path.nodes; // copy path
+                longest_path.nodes  = neighbor_path.nodes; // copy path
             }
         }
     }
@@ -306,8 +366,8 @@ VoronoiGraph::Path VoronoiGraphUtils::find_longest_path_on_circles(
     }
     const VoronoiGraph::Path &longest_branch = branches_item->second.top();
     longest_path.nodes.insert(longest_path.nodes.end(),
-                             longest_branch.nodes.begin(),
-                             longest_branch.nodes.end());
+                              longest_branch.nodes.begin(),
+                              longest_branch.nodes.end());
     return longest_path;
 }
 
@@ -334,7 +394,7 @@ std::optional<VoronoiGraph::Circle> VoronoiGraphUtils::create_circle(
 void VoronoiGraphUtils::merge_connected_circle(
     VoronoiGraph::ExPath::ConnectedCircles &dst,
     VoronoiGraph::ExPath::ConnectedCircles &src,
-    size_t dst_circle_count)
+    size_t                                  dst_circle_count)
 {
     std::set<size_t> done;
     for (const auto &item : src) {
@@ -362,8 +422,8 @@ void VoronoiGraphUtils::merge_connected_circle(
     }
 }
 
-void VoronoiGraphUtils::append_neighbor_branch(
-    VoronoiGraph::ExPath &dst, VoronoiGraph::ExPath &src)
+void VoronoiGraphUtils::append_neighbor_branch(VoronoiGraph::ExPath &dst,
+                                               VoronoiGraph::ExPath &src)
 {
     // move side branches
     if (!src.side_branches.empty())
@@ -415,14 +475,15 @@ void VoronoiGraphUtils::reshape_longest_path(VoronoiGraph::ExPath &path)
         std::reverse(side_branch.nodes.begin(), side_branch.nodes.end());
         VoronoiGraph::Path new_main_branch(std::move(branches.top()));
         branches.pop();
-        std::reverse(new_main_branch.nodes.begin(), new_main_branch.nodes.end());
+        std::reverse(new_main_branch.nodes.begin(),
+                     new_main_branch.nodes.end());
         // add old main path store into side branches - may be it is not neccessary
         branches.push(std::move(side_branch));
 
         // swap side branch with main branch
         path.nodes.erase(path.nodes.begin(), end_path);
         path.nodes.insert(path.nodes.begin(), new_main_branch.nodes.begin(),
-                         new_main_branch.nodes.end());
+                          new_main_branch.nodes.end());
 
         path.length += new_main_branch.length;
         path.length -= actual_length;
@@ -434,8 +495,8 @@ void VoronoiGraphUtils::reshape_longest_path(VoronoiGraph::ExPath &path)
 VoronoiGraph::ExPath VoronoiGraphUtils::create_longest_path(
     const VoronoiGraph::Node *start_node)
 {
-    VoronoiGraph::ExPath      longest_path;
-    CallStack call_stack;
+    VoronoiGraph::ExPath longest_path;
+    CallStack            call_stack;
     call_stack.emplace(
         std::make_unique<EvaluateNeighbor>(longest_path, start_node));
 
@@ -453,8 +514,7 @@ VoronoiGraph::ExPath VoronoiGraphUtils::create_longest_path(
     return longest_path;
 }
 
-Slic3r::Point VoronoiGraphUtils::get_edge_point(const VD::edge_type *edge,
-                                                double               ratio)
+Slic3r::Point VoronoiGraphUtils::get_edge_point(const VD::edge_type *edge, double ratio)
 {
     const VD::vertex_type *v0 = edge->vertex0();
     const VD::vertex_type *v1 = edge->vertex1();
@@ -464,7 +524,9 @@ Slic3r::Point VoronoiGraphUtils::get_edge_point(const VD::edge_type *edge,
         return Point(v1->x(), v1->y());
 
     if (edge->is_linear()) {
-        Point dir(v1->x() - v0->x(), v1->y() - v0->y());
+        Point dir(
+            v1->x() - v0->x(), 
+            v1->y() - v0->y());
         // normalize
         dir *= ratio;
         return Point(v0->x() + dir.x(), v0->y() + dir.y());
@@ -472,42 +534,62 @@ Slic3r::Point VoronoiGraphUtils::get_edge_point(const VD::edge_type *edge,
 
     assert(edge->is_curved());
     // TODO: distance on curve
-    return Point(v0->x(), v0->y());
+
+    // approx by line
+    Point dir(v1->x() - v0->x(), v1->y() - v0->y());
+    dir *= ratio;
+    return Point(v0->x() + dir.x(), v0->y() + dir.y());
 }
 
-Slic3r::Point VoronoiGraphUtils::get_center_of_path(
-    const VoronoiGraph::Nodes &path,
-                                            double path_length)
+Slic3r::Point VoronoiGraphUtils::get_point_on_path(const VoronoiGraph::Nodes &path, double distance)
 {
     const VoronoiGraph::Node *prev_node        = nullptr;
-    double                    half_path_length = path_length / 2.;
-    double                    distance         = 0.;
+    double                    actual_distance         = 0.;
     for (const VoronoiGraph::Node *node : path) {
         if (prev_node == nullptr) { // first call
             prev_node = node;
             continue;
         }
-        const VoronoiGraph::Node::Neighbor *neighbor = get_neighbor(prev_node,
-                                                                    node);
-        distance += neighbor->edge_length;
-        if (distance >= half_path_length) {
+        const VoronoiGraph::Node::Neighbor *neighbor = get_neighbor(prev_node, node);
+        actual_distance += neighbor->edge_length;
+        if (actual_distance >= distance) {
             // over half point is on
-            double ratio = 1. - (distance - half_path_length) /
-                                    neighbor->edge_length;
+            double previous_distance = actual_distance - distance;
+            double over_ratio = previous_distance / neighbor->edge_length;
+            double ratio = 1. - over_ratio;
             return get_edge_point(neighbor->edge, ratio);
         }
         prev_node = node;
     }
-    // half_path_length must be inside path
+    // distance must be inside path
     // this means bad input params
     assert(false);
     return Point(0, 0);
 }
 
+std::vector<Slic3r::Point> VoronoiGraphUtils::sample_longest_path(
+    const VoronoiGraph::ExPath &longest_path, const SampleConfig &config)
+{
+    // 1) One support point
+    if (longest_path.length <
+        config.max_length_for_one_support_point) { // create only one
+                                                   // point in center
+        // sample in center of voronoi
+        return {get_center_of_path(longest_path.nodes, longest_path.length)};
+    }
+    // 2) Two support points
+    //if (longest_path.length < config.max_distance) {}
+
+    std::vector<Point> points;
+    points.push_back(get_offseted_point(*longest_path.nodes.front(), config.start_distance));
+
+    return points;
+}
+
 std::vector<Slic3r::Point> VoronoiGraphUtils::sample_voronoi_graph(
     const VoronoiGraph &  graph,
-                                        const SampleConfig &  config,
-                                        VoronoiGraph::ExPath &longest_path)
+    const SampleConfig &  config,
+    VoronoiGraph::ExPath &longest_path)
 {
     // first vertex on contour:
     const VoronoiGraph::Node *start_node = nullptr;
@@ -521,23 +603,14 @@ std::vector<Slic3r::Point> VoronoiGraphUtils::sample_voronoi_graph(
     }
     // every island has to have a point on contour
     assert(start_node != nullptr);
-
     longest_path = create_longest_path(start_node);
     // longest_path = create_longest_path_recursive(start_node);
-    if (longest_path.length <
-        config.max_length_for_one_support_point) { // create only one
-                                                   // point in center
-        // sample in center of voronoi
-        return {get_center_of_path(longest_path.nodes, longest_path.length)};
-    }
-
-    std::vector<Point> points;
-    points.push_back(get_offseted_point(*start_node, config.start_distance));
-
-    return points;
+    return sample_longest_path(longest_path, config);
 }
 
-void VoronoiGraphUtils::draw(SVG &svg, const VoronoiGraph &graph, coord_t width)
+void VoronoiGraphUtils::draw(SVG &               svg,
+                             const VoronoiGraph &graph,
+                             coord_t             width)
 {
     for (const auto &[key, value] : graph.data) {
         svg.draw(Point(key->x(), key->y()), "lightgray", width);
@@ -618,3 +691,4 @@ void VoronoiGraphUtils::draw(SVG &                       svg,
 
     draw(svg, path.nodes, width, mainPathColor);
 }
+

src/libslic3r/SLA/SupportIslands/VoronoiGraphUtils.hpp

@@ -6,11 +6,11 @@
 #include <set>
 #include <libslic3r/Geometry.hpp>
 #include <libslic3r/Point.hpp>
-#include "VoronoiGraph.hpp"
-#include "SampleConfig.hpp"
+#include <libslic3r/SVG.hpp>
 
-// for debug draw purpose
-#include "SVG.hpp"
+#include "VoronoiGraph.hpp"
+#include "Parabola.hpp"
+#include "SampleConfig.hpp"
 
 namespace Slic3r::sla {
 
@@ -18,33 +18,63 @@ namespace Slic3r::sla {
 /// Class which contain collection of static function
 /// for work with Voronoi Graph.
 /// </summary>
-class VoronoiGraphUtils {
+class VoronoiGraphUtils
+{
     using VD = Slic3r::Geometry::VoronoiDiagram;
+
 public:
     VoronoiGraphUtils() = delete;
 
     // return node from graph by vertex, when no exists create one
-    static VoronoiGraph::Node *getNode(
-        VoronoiGraph &         graph,
-        const VD::vertex_type *vertex,
-        const VD::edge_type *  edge,
-        const Lines &          lines
-    );
-      
+    static VoronoiGraph::Node *getNode(VoronoiGraph &         graph,
+                                       const VD::vertex_type *vertex,
+                                       const VD::edge_type *  edge,
+                                       const Lines &          lines);
+
+    /// <summary>
+    /// extract parabola focus point from lines belongs to cell
+    /// </summary>
+    /// <param name="lines">source of Voronoi diagram</param>
+    /// <param name="cell">cell inside of Voronoi diagram</param>
+    /// <returns>Point from lines which create parabola.</returns>
+    static Point retrieve_point(const Lines &lines, const VD::cell_type &cell);
+    static Slic3r::Point get_parabola_point(const VD::edge_type &parabola, const Slic3r::Lines &lines);
+    static Slic3r::Line get_parabola_line(const VD::edge_type &parabola, const Lines &lines);
+    static Parabola get_parabola(const VD::edge_type &edge, const Lines &lines);
+
     /// <summary>
     /// Calculate length
     /// </summary>
     /// <param name="edge">curved edge</param>
     /// <returns>edge length</returns>
-    static double calculate_length_of_parabola(const VD::edge_type &edge,
-        const std::vector<VD::segment_type> &segments);
+    static double calculate_length_of_parabola(
+        const VD::edge_type &                edge,
+        const Lines &                        lines);
+
+    /// <summary>
+    /// Calculate length of edge line segment or curve - parabola.
+    /// </summary>
+    /// <param name="edge">Input edge to calcuate length</param>
+    /// <param name="lines">Source for Voronoi diagram. It contains parabola parameters</param>
+    /// <returns>The length of edge</returns>
+    static double calculate_length(const VD::edge_type &edge,
+                                   const Lines &        lines);
+
+    /// <summary>
+    /// Calculate maximal distance to outline and multiply by two(must be similar on both side)
+    /// </summary>
+    /// <param name="edge">Input edge.</param>
+    /// <param name="lines">Source for Voronoi diagram. It contains parabola parameters</param>
+    /// <returns>Maximal island width along edge</returns>
+    static double calculate_max_width(const VD::edge_type &edge,
+                                      const Lines &        lines);
 
     /// <summary>
     /// calculate distances to border of island and length on skeleton
     /// </summary>
     /// <param name="voronoi_diagram">Input anotated voronoi diagram
     /// (use function Slic3r::Voronoi::annotate_inside_outside)</param>
-    /// <param name="lines">Source lines for voronoi diagram</param> 
+    /// <param name="lines">Source lines for voronoi diagram</param>
     /// <returns>Extended voronoi graph by distances and length</returns>
     static VoronoiGraph getSkeleton(const VD &vd, const Lines &lines);
 
@@ -55,8 +85,8 @@ public:
     /// <param name="padding">Distance from outline</param>
     /// <returns></returns>
     static Slic3r::Point get_offseted_point(const VoronoiGraph::Node &node,
-                                          double                    padding);
-    
+                                            double padding);
+
     /// <summary>
     /// find neighbor and return distance between nodes
     /// </summary>
@@ -75,27 +105,29 @@ public:
     /// <param name="to">destination Node</param>
     /// <returns>distance between Nodes or Assert when not neighbor</returns>
     static double get_neighbor_distance(const VoronoiGraph::Node *from,
-                                 const VoronoiGraph::Node *to);
+                                        const VoronoiGraph::Node *to);
 
     /// <summary>
     /// Create longest node path over circle together with side branches
     /// </summary>
-    /// <param name="circle">Source circle, can't be connected with another circle</param>
-    /// <param name="side_branches">Circle side branches from nodes of circle</param>
-    /// <param name="start_path">Path before circle - defince input point to circle</param>
-    /// <returns>Longest nodes path and its length</returns>
+    /// <param name="circle">Source circle, can't be connected with another
+    /// circle</param> <param name="side_branches">Circle side branches from
+    /// nodes of circle</param> <param name="start_path">Path before circle -
+    /// defince input point to circle</param> <returns>Longest nodes path and
+    /// its length</returns>
     static VoronoiGraph::Path find_longest_path_on_circle(
         const VoronoiGraph::Circle &                 circle,
         const VoronoiGraph::ExPath::SideBranchesMap &side_branches);
-        
+
     /// <summary>
     /// Serach longest path from input_node throw Nodes in connected circles,
     /// when circle is alone call find_longest_path_on_circle.
     /// </summary>
     /// <param name="input_node">Node on circle</param>
-    /// <param name="finished_circle_index">index of circle with input node</param> 
-    /// <param name="ex_path">Hold Circles, connection of circles and Side branches</param> 
-    /// <returns>Longest path from input node</returns>
+    /// <param name="finished_circle_index">index of circle with input
+    /// node</param> <param name="ex_path">Hold Circles, connection of circles
+    /// and Side branches</param> <returns>Longest path from input
+    /// node</returns>
     static VoronoiGraph::Path find_longest_path_on_circles(
         const VoronoiGraph::Node &  input_node,
         size_t                      finished_circle_index,
@@ -115,21 +147,22 @@ public:
     /// Move source connected circles into destination
     /// </summary>
     /// <param name="dst">In/Out param</param>
-    /// <param name="src">Input possible modified, do not use it after this function</param>
-    /// <param name="dst_circle_count">Count of destination circles before
-    /// merge Source circle are append afted destination, therfore all src
-    /// indexes must be increased by destination circle count</param>
+    /// <param name="src">Input possible modified, do not use it after this
+    /// function</param> <param name="dst_circle_count">Count of destination
+    /// circles before merge Source circle are append afted destination, therfore
+    /// all src indexes must be increased by destination circle count</param>
     static void merge_connected_circle(
         VoronoiGraph::ExPath::ConnectedCircles &dst,
         VoronoiGraph::ExPath::ConnectedCircles &src,
-        size_t dst_circle_count);
+        size_t                                  dst_circle_count);
 
     /// <summary>
     /// move data from source to destination
     /// side_branches + circles + connected_circle
     /// </summary>
-    /// <param name="dst">destination extended path - append data from source</param>
-    /// <param name="src">source extended path - data will be moved to dst</param>
+    /// <param name="dst">destination extended path - append data from
+    /// source</param> <param name="src">source extended path - data will be
+    /// moved to dst</param>
     static void append_neighbor_branch(VoronoiGraph::ExPath &dst,
                                        VoronoiGraph::ExPath &src);
 
@@ -150,23 +183,42 @@ public:
     /// <param name="start_node">Random point from outline.</param>
     static VoronoiGraph::ExPath create_longest_path(
         const VoronoiGraph::Node *start_node);
-        
+
     /// <summary>
     /// Create point on edge defined by neighbor
     /// in distance defined by edge length ratio
     /// </summary>
     static Point get_edge_point(const VD::edge_type *edge, double ratio);
 
+    /// <summary>
+    /// Find point lay on path with distance from first point on path
+    /// </summary>
+    /// <param name="path">Neighbor connected Nodes</param>
+    /// <param name="distance">Distance to final point</param>
+    /// <returns>Points with distance to first node</returns>
+    static Point get_point_on_path(const VoronoiGraph::Nodes &path, double distance);
+
     /// <summary>
     /// Find point lay in center of path
-    /// Distance from this point to front of path 
+    /// Distance from this point to front of path
     /// is same as distance to back of path
     /// </summary>
     /// <param name="path">Queue of neighbor nodes.(must be neighbor)</param>
     /// <param name="path_length">length of path</param>
     /// <returns>Point laying on voronoi diagram</returns>
-    static Point get_center_of_path(const VoronoiGraph::Nodes &path,
-                             double                     path_length);
+    static Point get_center_of_path(const VoronoiGraph::Nodes &path, double path_length)
+    { return get_point_on_path(path, path_length / 2); }
+
+    /// <summary>
+    /// decide how to sample longest path
+    /// </summary>
+    /// <param name="longest_path">Path inside voronoi diagram with all side branches and circles</param>
+    /// <param name="config">Definition how to sample</param>
+    /// <returns>Support points for island</returns>
+    static std::vector<Point> sample_longest_path(
+        const VoronoiGraph::ExPath &longest_path, 
+        const SampleConfig &config
+    );
 
     /// <summary>
     /// Sample voronoi skeleton
@@ -176,18 +228,21 @@ public:
     /// <param name="longest_path">OUTPUT: longest path in graph</param>
     /// <returns>Vector of sampled points or Empty when distance from edge is
     /// bigger than max_distance</returns>
-    static std::vector<Point> sample_voronoi_graph(const VoronoiGraph &  graph,
-                                            const SampleConfig &  config,
-                                            VoronoiGraph::ExPath &longest_path);
+    static std::vector<Point> sample_voronoi_graph(
+        const VoronoiGraph &  graph,
+        const SampleConfig &  config,
+        VoronoiGraph::ExPath &longest_path);
 
-public: //draw function for debug
+public: // draw function for debug
     static void draw(SVG &svg, const VoronoiGraph &graph, coord_t width);
     static void draw(SVG &                      svg,
                      const VoronoiGraph::Nodes &path,
                      coord_t                    width,
                      const char *               color,
                      bool                       finish = false);
-    static void draw(SVG &svg, const VoronoiGraph::ExPath &path, coord_t width);
+    static void draw(SVG &                       svg,
+                     const VoronoiGraph::ExPath &path,
+                     coord_t                     width);
 };
 
 } // namespace Slic3r::sla

tests/sla_print/CMakeLists.txt

@@ -4,6 +4,7 @@ add_executable(${_TEST_NAME}_tests ${_TEST_NAME}_tests_main.cpp
     sla_test_utils.hpp sla_test_utils.cpp
     sla_supptgen_tests.cpp
     sla_raycast_tests.cpp
+    sla_parabola_tests.cpp
     sla_supptreeutils_tests.cpp
     sla_archive_readwrite_tests.cpp
     sla_zcorrection_tests.cpp)
@@ -18,3 +19,11 @@ endif()
 
 # catch_discover_tests(${_TEST_NAME}_tests TEST_PREFIX "${_TEST_NAME}: ")
 add_test(${_TEST_NAME}_tests ${_TEST_NAME}_tests ${CATCH_EXTRA_ARGS})
+
+if (WIN32)
+    # Adds a post-build copy of libgmp-10.dll
+    add_custom_command(TARGET ${_TEST_NAME}_tests POST_BUILD
+          COMMAND ${CMAKE_COMMAND} -E copy_if_different
+          "${CMAKE_PREFIX_PATH}/bin/libgmp-10.dll"
+          $<TARGET_FILE_DIR:${_TEST_NAME}_tests>)
+endif()

tests/sla_print/sla_parabola_tests.cpp

@@ -0,0 +1,54 @@
+#include "sla_test_utils.hpp"
+
+#include <libslic3r/SLA/SupportIslands/ParabolaUtils.hpp>
+
+using namespace Slic3r;
+using namespace Slic3r::sla;
+
+void parabola_check(const Parabola &parabola,
+                    const Point &        from,
+                    const Point &        to)
+{
+    auto   diffPoint = to - from;
+    double min       = sqrt(diffPoint.x() * diffPoint.x() +
+                      diffPoint.y() * diffPoint.y());
+    double max       = static_cast<double>(diffPoint.x()) + diffPoint.y();
+
+    double len = ParabolaUtils::calculate_length_of_parabola(parabola, from, to);
+    double len2 = ParabolaUtils::calculate_length_of_parabola_by_sampling(
+            parabola, from, to, 1.);
+
+    CHECK(fabs(len2 - len) < 1.);
+    CHECK(len >= min);
+    CHECK(len <= max);
+}
+
+// after generalization put to ParabolaUtils
+double getParabolaY(const Parabola &parabola, double x)
+{
+    double f    = ParabolaUtils::focal_length(parabola);
+    Vec2d  perp = parabola.directrix.normal().cast<double>();
+    // work only for test cases
+    if (perp.y() > 0.) perp *= -1.;
+    perp.normalize();
+    Vec2d v = parabola.focus.cast<double>() + perp * f;
+    return 1 / (4 * f) * (x - v.x()) * (x - v.x()) + v.y();
+}
+
+TEST_CASE("Parabola length", "[SupGen][Voronoi][Parabola]")
+{
+    using namespace Slic3r::sla;
+    double scale = 1e6;
+    // U shape parabola
+    Parabola parabola_x2(Line({-1. * scale, -.25 * scale},
+                              {1. * scale, -.25 * scale}),
+                         Point(0. * scale, .25 * scale));
+
+    double from_x = 1 * scale;
+    double to_x   = 3 * scale;
+    Point  from(from_x, getParabolaY(parabola_x2, from_x));
+    Point  to(to_x, getParabolaY(parabola_x2, to_x));
+    parabola_check(parabola_x2, from, to);
+} 
+
+

tests/sla_print/sla_supptgen_tests.cpp

@@ -6,6 +6,7 @@
 #include <libslic3r/SLA/SpatIndex.hpp>
 
 #include <libslic3r/SLA/SupportIslands/SampleConfig.hpp>
+#include <libslic3r/SLA/SupportIslands/VoronoiGraphUtils.hpp>
 
 #include "sla_test_utils.hpp"
 
@@ -250,7 +251,6 @@ TEST_CASE("Sampling speed test on FrogLegs", "[VoronoiSkeleton]")
     cfg.max_distance   = size + 0.1;
     cfg.sample_size    = size / 5;
     cfg.start_distance = 0.2 * size; // radius of support head
-    cfg.curve_sample   = 0.1 * size;
     cfg.max_length_for_one_support_point = 3 * size;
 
     for (int i = 0; i < 100; ++i) {
@@ -259,7 +259,6 @@ TEST_CASE("Sampling speed test on FrogLegs", "[VoronoiSkeleton]")
     }
 }
 
-
 TEST_CASE("Small islands should be supported in center", "[SupGen][VoronoiSkeleton]")
 {
     double size = 3e7;
@@ -267,7 +266,6 @@ TEST_CASE("Small islands should be supported in center", "[SupGen][VoronoiSkelet
     cfg.max_distance   = size + 0.1;
     cfg.sample_size    = size / 5;
     cfg.start_distance = 0.2 * size; // radius of support head
-    cfg.curve_sample   = 0.1 * size;
     cfg.max_length_for_one_support_point = 3 * size;
 
     ExPolygons islands = createTestIslands(size);