Fixed conflicts after merge with master

src/libslic3r/AppConfig.cpp

@@ -101,6 +101,9 @@ void AppConfig::set_defaults()
     if (get("use_inches").empty())
         set("use_inches", "0");
 
+    if (get("show_splash_screen").empty())
+        set("show_splash_screen", "1");
+
     // Remove legacy window positions/sizes
     erase("", "main_frame_maximized");
     erase("", "main_frame_pos");

src/libslic3r/CMakeLists.txt

@@ -215,7 +215,9 @@ add_library(libslic3r STATIC
     SimplifyMeshImpl.hpp
     SimplifyMesh.cpp
     MarchingSquares.hpp
-    Optimizer.hpp
+    Optimize/Optimizer.hpp
+    Optimize/NLoptOptimizer.hpp
+    Optimize/BruteforceOptimizer.hpp
     ${OpenVDBUtils_SOURCES}
     SLA/Pad.hpp
     SLA/Pad.cpp

src/libslic3r/Fill/Fill.cpp

@@ -318,7 +318,7 @@ void export_group_fills_to_svg(const char *path, const std::vector<SurfaceFill>
 #endif
 
 // friend to Layer
-void Layer::make_fills(FillAdaptive_Internal::Octree* adaptive_fill_octree)
+void Layer::make_fills(FillAdaptive_Internal::Octree* adaptive_fill_octree, FillAdaptive_Internal::Octree* support_fill_octree)
 {
 	for (LayerRegion *layerm : m_regions)
 		layerm->fills.clear();
@@ -346,6 +346,7 @@ void Layer::make_fills(FillAdaptive_Internal::Octree* adaptive_fill_octree)
         f->z 		= this->print_z;
         f->angle 	= surface_fill.params.angle;
         f->adapt_fill_octree = adaptive_fill_octree;
+        f->support_fill_octree = support_fill_octree;
 
         // calculate flow spacing for infill pattern generation
         bool using_internal_flow = ! surface_fill.surface.is_solid() && ! surface_fill.params.flow.bridge;

src/libslic3r/Fill/FillAdaptive.cpp

@@ -35,7 +35,7 @@ std::pair<double, double> adaptive_fill_line_spacing(const PrintObject &print_ob
         const PrintRegionConfig &config   = region->config();
         bool                     nonempty = config.fill_density > 0;
         bool                     has_adaptive_infill = nonempty && config.fill_pattern == ipAdaptiveCubic;
-        bool                     has_support_infill  = nonempty && false; // config.fill_pattern == icSupportCubic;
+        bool                     has_support_infill  = nonempty && config.fill_pattern == ipSupportCubic;
         region_fill_data.push_back(RegionFillData({
             has_adaptive_infill ? Tristate::Maybe : Tristate::No,
             has_support_infill ? Tristate::Maybe : Tristate::No,
@@ -90,19 +90,32 @@ std::pair<double, double> adaptive_fill_line_spacing(const PrintObject &print_ob
     return std::make_pair(adaptive_line_spacing, support_line_spacing);
 }
 
-void FillAdaptive::_fill_surface_single(
-    const FillParams                &params, 
-    unsigned int                     thickness_layers,
-    const std::pair<float, Point>   &direction, 
-    ExPolygon                       &expolygon, 
-    Polylines                       &polylines_out)
+void FillAdaptive::_fill_surface_single(const FillParams &             params,
+                                        unsigned int                   thickness_layers,
+                                        const std::pair<float, Point> &direction,
+                                        ExPolygon &                    expolygon,
+                                        Polylines &                    polylines_out)
 {
+    if(this->adapt_fill_octree != nullptr)
+        this->generate_infill(params, thickness_layers, direction, expolygon, polylines_out, this->adapt_fill_octree);
+}
+
+void FillAdaptive::generate_infill(const FillParams &             params,
+                                        unsigned int                   thickness_layers,
+                                        const std::pair<float, Point> &direction,
+                                        ExPolygon &                    expolygon,
+                                        Polylines &                    polylines_out,
+                                        FillAdaptive_Internal::Octree *octree)
+{
+    Vec3d rotation = Vec3d((5.0 * M_PI) / 4.0, Geometry::deg2rad(215.264), M_PI / 6.0);
+    Transform3d rotation_matrix = Geometry::assemble_transform(Vec3d::Zero(), rotation, Vec3d::Ones(), Vec3d::Ones());
+
     // Store grouped lines by its direction (multiple of 120°)
     std::vector<Lines> infill_lines_dir(3);
-    this->generate_infill_lines(this->adapt_fill_octree->root_cube.get(),
-                                this->z, this->adapt_fill_octree->origin,infill_lines_dir,
-                                this->adapt_fill_octree->cubes_properties,
-                                int(this->adapt_fill_octree->cubes_properties.size()) - 1);
+    this->generate_infill_lines(octree->root_cube.get(),
+                                this->z, octree->origin, rotation_matrix,
+                                infill_lines_dir, octree->cubes_properties,
+                                int(octree->cubes_properties.size()) - 1);
 
     Polylines all_polylines;
     all_polylines.reserve(infill_lines_dir[0].size() * 3);
@@ -186,6 +199,7 @@ void FillAdaptive::generate_infill_lines(
         FillAdaptive_Internal::Cube *cube,
         double z_position,
         const Vec3d &origin,
+        const Transform3d &rotation_matrix,
         std::vector<Lines> &dir_lines_out,
         const std::vector<FillAdaptive_Internal::CubeProperties> &cubes_properties,
         int depth)
@@ -197,7 +211,8 @@ void FillAdaptive::generate_infill_lines(
         return;
     }
 
-    double z_diff = std::abs(z_position - cube->center.z());
+    Vec3d cube_center_tranformed = rotation_matrix * cube->center;
+    double z_diff = std::abs(z_position - cube_center_tranformed.z());
 
     if (z_diff > cubes_properties[depth].height / 2)
     {
@@ -208,14 +223,14 @@ void FillAdaptive::generate_infill_lines(
     {
         Point from(
                 scale_((cubes_properties[depth].diagonal_length / 2) * (cubes_properties[depth].line_z_distance - z_diff) / cubes_properties[depth].line_z_distance),
-                scale_(cubes_properties[depth].line_xy_distance - ((z_position - (cube->center.z() - cubes_properties[depth].line_z_distance)) / sqrt(2))));
+                scale_(cubes_properties[depth].line_xy_distance - ((z_position - (cube_center_tranformed.z() - cubes_properties[depth].line_z_distance)) / sqrt(2))));
         Point to(-from.x(), from.y());
         // Relative to cube center
 
         double rotation_angle = (2.0 * M_PI) / 3.0;
         for (Lines &lines : dir_lines_out)
         {
-            Vec3d offset = cube->center - origin;
+            Vec3d offset = cube_center_tranformed - (rotation_matrix * origin);
             Point from_abs(from), to_abs(to);
 
             from_abs.x() += int(scale_(offset.x()));
@@ -235,7 +250,7 @@ void FillAdaptive::generate_infill_lines(
     {
         if(child != nullptr)
         {
-            generate_infill_lines(child.get(), z_position, origin, dir_lines_out, cubes_properties, depth - 1);
+            generate_infill_lines(child.get(), z_position, origin, rotation_matrix, dir_lines_out, cubes_properties, depth - 1);
         }
     }
 }
@@ -301,14 +316,11 @@ std::unique_ptr<FillAdaptive_Internal::Octree> FillAdaptive::build_octree(
         triangle_mesh.require_shared_vertices();
     }
 
-    Vec3d rotation = Vec3d((5.0 * M_PI) / 4.0, Geometry::deg2rad(215.264), M_PI / 6.0);
-    Transform3d rotation_matrix = Geometry::assemble_transform(Vec3d::Zero(), rotation, Vec3d::Ones(), Vec3d::Ones());
-
     AABBTreeIndirect::Tree3f aabbTree = AABBTreeIndirect::build_aabb_tree_over_indexed_triangle_set(
             triangle_mesh.its.vertices, triangle_mesh.its.indices);
     auto octree = std::make_unique<Octree>(std::make_unique<Cube>(cube_center), cube_center, cubes_properties);
 
-    FillAdaptive::expand_cube(octree->root_cube.get(), cubes_properties, rotation_matrix, aabbTree, triangle_mesh, int(cubes_properties.size()) - 1);
+    FillAdaptive::expand_cube(octree->root_cube.get(), cubes_properties, aabbTree, triangle_mesh, int(cubes_properties.size()) - 1);
 
     return octree;
 }
@@ -316,7 +328,6 @@ std::unique_ptr<FillAdaptive_Internal::Octree> FillAdaptive::build_octree(
 void FillAdaptive::expand_cube(
     FillAdaptive_Internal::Cube *cube,
     const std::vector<FillAdaptive_Internal::CubeProperties> &cubes_properties,
-    const Transform3d &rotation_matrix,
     const AABBTreeIndirect::Tree3f &distance_tree,
     const TriangleMesh &triangle_mesh, int depth)
 {
@@ -328,8 +339,8 @@ void FillAdaptive::expand_cube(
     }
 
     std::vector<Vec3d> child_centers = {
-            Vec3d(-1, -1, -1), Vec3d( 1, -1, -1), Vec3d(-1,  1, -1), Vec3d(-1, -1,  1),
-            Vec3d( 1,  1,  1), Vec3d(-1,  1,  1), Vec3d( 1, -1,  1), Vec3d( 1,  1, -1)
+        Vec3d(-1, -1, -1), Vec3d( 1, -1, -1), Vec3d(-1,  1, -1), Vec3d( 1,  1, -1),
+        Vec3d(-1, -1,  1), Vec3d( 1, -1,  1), Vec3d(-1,  1,  1), Vec3d( 1,  1,  1)
     };
 
     double cube_radius_squared = (cubes_properties[depth].height * cubes_properties[depth].height) / 16;
@@ -337,15 +348,173 @@ void FillAdaptive::expand_cube(
     for (size_t i = 0; i < 8; ++i)
     {
         const Vec3d &child_center = child_centers[i];
-        Vec3d child_center_transformed = cube->center + rotation_matrix * (child_center * (cubes_properties[depth].edge_length / 4));
+        Vec3d child_center_transformed = cube->center + (child_center * (cubes_properties[depth].edge_length / 4));
 
         if(AABBTreeIndirect::is_any_triangle_in_radius(triangle_mesh.its.vertices, triangle_mesh.its.indices,
             distance_tree, child_center_transformed, cube_radius_squared))
         {
             cube->children[i] = std::make_unique<Cube>(child_center_transformed);
-            FillAdaptive::expand_cube(cube->children[i].get(), cubes_properties, rotation_matrix, distance_tree, triangle_mesh, depth - 1);
+            FillAdaptive::expand_cube(cube->children[i].get(), cubes_properties, distance_tree, triangle_mesh, depth - 1);
         }
     }
 }
 
+void FillAdaptive_Internal::Octree::propagate_point(
+    Vec3d                                                     point,
+    FillAdaptive_Internal::Cube *                             current,
+    int                                                       depth,
+    const std::vector<FillAdaptive_Internal::CubeProperties> &cubes_properties)
+{
+    using namespace FillAdaptive_Internal;
+
+    if(depth <= 0)
+    {
+        return;
+    }
+
+    size_t octant_idx = Octree::find_octant(point, current->center);
+    Cube * child = current->children[octant_idx].get();
+
+    // Octant not exists, then create it
+    if(child == nullptr) {
+        std::vector<Vec3d> child_centers = {
+            Vec3d(-1, -1, -1), Vec3d( 1, -1, -1), Vec3d(-1,  1, -1), Vec3d( 1,  1, -1),
+            Vec3d(-1, -1,  1), Vec3d( 1, -1,  1), Vec3d(-1,  1,  1), Vec3d( 1,  1,  1)
+        };
+
+        const Vec3d &child_center = child_centers[octant_idx];
+        Vec3d child_center_transformed = current->center + (child_center * (cubes_properties[depth].edge_length / 4));
+
+        current->children[octant_idx] = std::make_unique<Cube>(child_center_transformed);
+        child = current->children[octant_idx].get();
+    }
+
+    Octree::propagate_point(point, child, (depth - 1), cubes_properties);
+}
+
+std::unique_ptr<FillAdaptive_Internal::Octree> FillSupportCubic::build_octree(
+    TriangleMesh &     triangle_mesh,
+    coordf_t           line_spacing,
+    const Vec3d &      cube_center,
+    const Transform3d &rotation_matrix)
+{
+    using namespace FillAdaptive_Internal;
+
+    if(line_spacing <= 0 || std::isnan(line_spacing))
+    {
+        return nullptr;
+    }
+
+    Vec3d bb_size = triangle_mesh.bounding_box().size();
+    // The furthest point from the center of the bottom of the mesh bounding box.
+    double furthest_point = std::sqrt(((bb_size.x() * bb_size.x()) / 4.0) +
+                                      ((bb_size.y() * bb_size.y()) / 4.0) +
+                                      (bb_size.z() * bb_size.z()));
+    double max_cube_edge_length = furthest_point * 2;
+
+    std::vector<CubeProperties> cubes_properties;
+    for (double edge_length = (line_spacing * 2); edge_length < (max_cube_edge_length * 2); edge_length *= 2)
+    {
+        CubeProperties props{};
+        props.edge_length = edge_length;
+        props.height = edge_length * sqrt(3);
+        props.diagonal_length = edge_length * sqrt(2);
+        props.line_z_distance = edge_length / sqrt(3);
+        props.line_xy_distance = edge_length / sqrt(6);
+        cubes_properties.push_back(props);
+    }
+
+    if (triangle_mesh.its.vertices.empty())
+    {
+        triangle_mesh.require_shared_vertices();
+    }
+
+    AABBTreeIndirect::Tree3f aabbTree = AABBTreeIndirect::build_aabb_tree_over_indexed_triangle_set(
+        triangle_mesh.its.vertices, triangle_mesh.its.indices);
+
+    auto octree = std::make_unique<Octree>(std::make_unique<Cube>(cube_center), cube_center, cubes_properties);
+
+    double cube_edge_length = line_spacing / 2.0;
+    int max_depth = int(octree->cubes_properties.size()) - 1;
+    BoundingBoxf3 mesh_bb = triangle_mesh.bounding_box();
+    Vec3f vertical(0, 0, 1);
+
+    for (size_t facet_idx = 0; facet_idx < triangle_mesh.stl.facet_start.size(); ++facet_idx)
+    {
+        if(triangle_mesh.stl.facet_start[facet_idx].normal.dot(vertical) <= 0.707)
+        {
+            // The angle is smaller than PI/4, than infill don't to be there
+            continue;
+        }
+
+        stl_vertex v_1 = triangle_mesh.stl.facet_start[facet_idx].vertex[0];
+        stl_vertex v_2 = triangle_mesh.stl.facet_start[facet_idx].vertex[1];
+        stl_vertex v_3 = triangle_mesh.stl.facet_start[facet_idx].vertex[2];
+
+        std::vector<Vec3d> triangle_vertices =
+            {Vec3d(v_1.x(), v_1.y(), v_1.z()),
+             Vec3d(v_2.x(), v_2.y(), v_2.z()),
+             Vec3d(v_3.x(), v_3.y(), v_3.z())};
+
+        BoundingBoxf3 triangle_bb(triangle_vertices);
+
+        Vec3d triangle_start_relative = triangle_bb.min - mesh_bb.min;
+        Vec3d triangle_end_relative   = triangle_bb.max - mesh_bb.min;
+
+        Vec3crd triangle_start_idx = Vec3crd(
+            int(std::floor(triangle_start_relative.x() / cube_edge_length)),
+            int(std::floor(triangle_start_relative.y() / cube_edge_length)),
+            int(std::floor(triangle_start_relative.z() / cube_edge_length)));
+        Vec3crd triangle_end_idx = Vec3crd(
+            int(std::floor(triangle_end_relative.x() / cube_edge_length)),
+            int(std::floor(triangle_end_relative.y() / cube_edge_length)),
+            int(std::floor(triangle_end_relative.z() / cube_edge_length)));
+
+        for (int z = triangle_start_idx.z(); z <= triangle_end_idx.z(); ++z)
+        {
+            for (int y = triangle_start_idx.y(); y <= triangle_end_idx.y(); ++y)
+            {
+                for (int x = triangle_start_idx.x(); x <= triangle_end_idx.x(); ++x)
+                {
+                    Vec3d cube_center_relative(x * cube_edge_length + (cube_edge_length / 2.0), y * cube_edge_length + (cube_edge_length / 2.0), z * cube_edge_length);
+                    Vec3d cube_center_absolute = cube_center_relative + mesh_bb.min;
+
+                    double cube_center_absolute_arr[3] = {cube_center_absolute.x(), cube_center_absolute.y(), cube_center_absolute.z()};
+                    double distance = 0, cord_u = 0, cord_v = 0;
+
+                    double dir[3] = {0.0, 0.0, 1.0};
+
+                    double vert_0[3] = {triangle_vertices[0].x(),
+                                        triangle_vertices[0].y(),
+                                        triangle_vertices[0].z()};
+                    double vert_1[3] = {triangle_vertices[1].x(),
+                                        triangle_vertices[1].y(),
+                                        triangle_vertices[1].z()};
+                    double vert_2[3] = {triangle_vertices[2].x(),
+                                        triangle_vertices[2].y(),
+                                        triangle_vertices[2].z()};
+
+                    if(intersect_triangle(cube_center_absolute_arr, dir, vert_0, vert_1, vert_2, &distance, &cord_u, &cord_v) && distance > 0 && distance <= cube_edge_length)
+                    {
+                        Vec3d cube_center_transformed(cube_center_absolute.x(), cube_center_absolute.y(), cube_center_absolute.z() + (cube_edge_length / 2.0));
+                        Octree::propagate_point(rotation_matrix * cube_center_transformed, octree->root_cube.get(), max_depth, octree->cubes_properties);
+                    }
+                }
+            }
+        }
+    }
+
+    return octree;
+}
+
+void FillSupportCubic::_fill_surface_single(const FillParams &             params,
+                                            unsigned int                   thickness_layers,
+                                            const std::pair<float, Point> &direction,
+                                            ExPolygon &                    expolygon,
+                                            Polylines &                    polylines_out)
+{
+    if (this->support_fill_octree != nullptr)
+        this->generate_infill(params, thickness_layers, direction, expolygon, polylines_out, this->support_fill_octree);
+}
+
 } // namespace Slic3r

src/libslic3r/Fill/FillAdaptive.hpp

@@ -35,6 +35,17 @@ namespace FillAdaptive_Internal
 
         Octree(std::unique_ptr<Cube> rootCube, const Vec3d &origin, const std::vector<CubeProperties> &cubes_properties)
             : root_cube(std::move(rootCube)), origin(origin), cubes_properties(cubes_properties) {}
+
+        inline static int find_octant(const Vec3d &i_cube, const Vec3d &current)
+        {
+            return (i_cube.z() > current.z()) * 4 + (i_cube.y() > current.y()) * 2 + (i_cube.x() > current.x());
+        }
+
+        static void propagate_point(
+            Vec3d                        point,
+            FillAdaptive_Internal::Cube *current_cube,
+            int                          depth,
+            const std::vector<FillAdaptive_Internal::CubeProperties> &cubes_properties);
     };
 }; // namespace FillAdaptive_Internal
 
@@ -63,12 +74,20 @@ protected:
         FillAdaptive_Internal::Cube *cube,
         double                       z_position,
         const Vec3d &                origin,
+        const Transform3d &          rotation_matrix,
         std::vector<Lines> &         dir_lines_out,
         const std::vector<FillAdaptive_Internal::CubeProperties> &cubes_properties,
         int  depth);
 
     static void connect_lines(Lines &lines, Line new_line);
 
+    void generate_infill(const FillParams &             params,
+                         unsigned int                   thickness_layers,
+                         const std::pair<float, Point> &direction,
+                         ExPolygon &                    expolygon,
+                         Polylines &                    polylines_out,
+                         FillAdaptive_Internal::Octree *octree);
+
 public:
     static std::unique_ptr<FillAdaptive_Internal::Octree> build_octree(
         TriangleMesh &triangle_mesh,
@@ -78,12 +97,36 @@ public:
     static void expand_cube(
         FillAdaptive_Internal::Cube *cube,
         const std::vector<FillAdaptive_Internal::CubeProperties> &cubes_properties,
-        const Transform3d &             rotation_matrix,
         const AABBTreeIndirect::Tree3f &distance_tree,
         const TriangleMesh &            triangle_mesh,
         int                             depth);
 };
 
+class FillSupportCubic : public FillAdaptive
+{
+public:
+    virtual ~FillSupportCubic() = default;
+
+protected:
+    virtual Fill* clone() const { return new FillSupportCubic(*this); };
+
+    virtual bool no_sort() const { return true; }
+
+    virtual void _fill_surface_single(
+        const FillParams                &params,
+        unsigned int                     thickness_layers,
+        const std::pair<float, Point>   &direction,
+        ExPolygon                       &expolygon,
+        Polylines                       &polylines_out);
+
+public:
+    static std::unique_ptr<FillAdaptive_Internal::Octree> build_octree(
+        TriangleMesh &     triangle_mesh,
+        coordf_t           line_spacing,
+        const Vec3d &      cube_center,
+        const Transform3d &rotation_matrix);
+};
+
 // Calculate line spacing for
 // 1) adaptive cubic infill
 // 2) adaptive internal support cubic infill

src/libslic3r/Fill/FillBase.cpp

@@ -39,6 +39,7 @@ Fill* Fill::new_from_type(const InfillPattern type)
     case ipHilbertCurve:        return new FillHilbertCurve();
     case ipOctagramSpiral:      return new FillOctagramSpiral();
     case ipAdaptiveCubic:       return new FillAdaptive();
+    case ipSupportCubic:        return new FillSupportCubic();
     default: throw std::invalid_argument("unknown type");
     }
 }

src/libslic3r/Fill/FillBase.hpp

@@ -73,7 +73,10 @@ public:
     // In scaled coordinates. Bounding box of the 2D projection of the object.
     BoundingBox bounding_box;
 
+    // Octree builds on mesh for usage in the adaptive cubic infill
     FillAdaptive_Internal::Octree* adapt_fill_octree = nullptr;
+    // Octree builds on mesh for usage in the support cubic infill
+    FillAdaptive_Internal::Octree* support_fill_octree = nullptr;
 
 public:
     virtual ~Fill() {}

src/libslic3r/Layer.hpp

@@ -138,8 +138,8 @@ public:
         return false;
     }
     void                    make_perimeters();
-    void                    make_fills() { this->make_fills(nullptr); };
-    void                    make_fills(FillAdaptive_Internal::Octree* adaptive_fill_octree);
+    void                    make_fills() { this->make_fills(nullptr, nullptr); };
+    void                    make_fills(FillAdaptive_Internal::Octree* adaptive_fill_octree, FillAdaptive_Internal::Octree* support_fill_octree);
     void 					make_ironing();
 
     void                    export_region_slices_to_svg(const char *path) const;

src/libslic3r/Optimize/BruteforceOptimizer.hpp

@@ -0,0 +1,140 @@
+#ifndef BRUTEFORCEOPTIMIZER_HPP
+#define BRUTEFORCEOPTIMIZER_HPP
+
+#include <libslic3r/Optimize/Optimizer.hpp>
+
+namespace Slic3r { namespace opt {
+
+namespace detail {
+// Implementing a bruteforce optimizer
+
+// Return the number of iterations needed to reach a specific grid position (idx)
+template<size_t N>
+long num_iter(const std::array<size_t, N> &idx, size_t gridsz)
+{
+    long ret = 0;
+    for (size_t i = 0; i < N; ++i) ret += idx[i] * std::pow(gridsz, i);
+    return ret;
+}
+
+// Implementation of a grid search where the search interval is sampled in
+// equidistant points for each dimension. Grid size determines the number of
+// samples for one dimension so the number of function calls is gridsize ^ dimension.
+struct AlgBurteForce {
+    bool to_min;
+    StopCriteria stc;
+    size_t gridsz;
+
+    AlgBurteForce(const StopCriteria &cr, size_t gs): stc{cr}, gridsz{gs} {}
+
+    // This function is called recursively for each dimension and generates
+    // the grid values for the particular dimension. If D is less than zero,
+    // the object function input values are generated for each dimension and it
+    // can be evaluated. The current best score is compared with the newly
+    // returned score and changed appropriately.
+    template<int D, size_t N, class Fn, class Cmp>
+    bool run(std::array<size_t, N> &idx,
+             Result<N> &result,
+             const Bounds<N> &bounds,
+             Fn &&fn,
+             Cmp &&cmp)
+    {
+        if (stc.stop_condition()) return false;
+
+        if constexpr (D < 0) { // Let's evaluate fn
+            Input<N> inp;
+
+            auto max_iter = stc.max_iterations();
+            if (max_iter && num_iter(idx, gridsz) >= max_iter)
+                return false;
+
+            for (size_t d = 0; d < N; ++d) {
+                const Bound &b = bounds[d];
+                double step = (b.max() - b.min()) / (gridsz - 1);
+                inp[d] = b.min() + idx[d] * step;
+            }
+
+            auto score = fn(inp);
+            if (cmp(score, result.score)) { // Change current score to the new
+                double absdiff = std::abs(score - result.score);
+
+                result.score = score;
+                result.optimum = inp;
+
+                // Check if the required precision is reached.
+                if (absdiff < stc.abs_score_diff() ||
+                    absdiff < stc.rel_score_diff() * std::abs(score))
+                    return false;
+            }
+
+        } else {
+            for (size_t i = 0; i < gridsz; ++i) {
+                idx[D] = i; // Mark the current grid position and dig down
+                if (!run<D - 1>(idx, result, bounds, std::forward<Fn>(fn),
+                                std::forward<Cmp>(cmp)))
+                    return false;
+            }
+        }
+
+        return true;
+    }
+
+    template<class Fn, size_t N>
+    Result<N> optimize(Fn&& fn,
+                       const Input<N> &/*initvals*/,
+                       const Bounds<N>& bounds)
+    {
+        std::array<size_t, N> idx = {};
+        Result<N> result;
+
+        if (to_min) {
+            result.score = std::numeric_limits<double>::max();
+            run<int(N) - 1>(idx, result, bounds, std::forward<Fn>(fn),
+                            std::less<double>{});
+        }
+        else {
+            result.score = std::numeric_limits<double>::lowest();
+            run<int(N) - 1>(idx, result, bounds, std::forward<Fn>(fn),
+                            std::greater<double>{});
+        }
+
+        return result;
+    }
+};
+
+} // namespace detail
+
+using AlgBruteForce = detail::AlgBurteForce;
+
+template<>
+class Optimizer<AlgBruteForce> {
+    AlgBruteForce m_alg;
+
+public:
+
+    Optimizer(const StopCriteria &cr = {}, size_t gridsz = 100)
+        : m_alg{cr, gridsz}
+    {}
+
+    Optimizer& to_max() { m_alg.to_min = false; return *this; }
+    Optimizer& to_min() { m_alg.to_min = true;  return *this; }
+
+    template<class Func, size_t N>
+    Result<N> optimize(Func&& func,
+                       const Input<N> &initvals,
+                       const Bounds<N>& bounds)
+    {
+        return m_alg.optimize(std::forward<Func>(func), initvals, bounds);
+    }
+
+    Optimizer &set_criteria(const StopCriteria &cr)
+    {
+        m_alg.stc = cr; return *this;
+    }
+
+    const StopCriteria &get_criteria() const { return m_alg.stc; }
+};
+
+}} // namespace Slic3r::opt
+
+#endif // BRUTEFORCEOPTIMIZER_HPP

src/libslic3r/Optimizer.hpp → src/libslic3r/Optimize/NLoptOptimizer.hpp

@@ -12,133 +12,10 @@
 #endif
 
 #include <utility>
-#include <tuple>
-#include <array>
-#include <cmath>
-#include <functional>
-#include <limits>
-#include <cassert>
 
-namespace Slic3r { namespace opt {
-
-// A type to hold the complete result of the optimization.
-template<size_t N> struct Result {
-    int resultcode;
-    std::array<double, N> optimum;
-    double score;
-};
-
-// An interval of possible input values for optimization
-class Bound {
-    double m_min, m_max;
-
-public:
-    Bound(double min = std::numeric_limits<double>::min(),
-          double max = std::numeric_limits<double>::max())
-        : m_min(min), m_max(max)
-    {}
-
-    double min() const noexcept { return m_min; }
-    double max() const noexcept { return m_max; }
-};
+#include <libslic3r/Optimize/Optimizer.hpp>
 
-// Helper types for optimization function input and bounds
-template<size_t N> using Input = std::array<double, N>;
-template<size_t N> using Bounds = std::array<Bound, N>;
-
-// A type for specifying the stop criteria. Setter methods can be concatenated
-class StopCriteria {
-
-    // If the absolute value difference between two scores.
-    double m_abs_score_diff = std::nan("");
-
-    // If the relative value difference between two scores.
-    double m_rel_score_diff = std::nan("");
-
-    // Stop if this value or better is found.
-    double m_stop_score = std::nan("");
-
-    // A predicate that if evaluates to true, the optimization should terminate
-    // and the best result found prior to termination should be returned.
-    std::function<bool()> m_stop_condition = [] { return false; };
-
-    // The max allowed number of iterations.
-    unsigned m_max_iterations = 0;
-
-public:
-
-    StopCriteria & abs_score_diff(double val)
-    {
-        m_abs_score_diff = val; return *this;
-    }
-
-    double abs_score_diff() const { return m_abs_score_diff; }
-
-    StopCriteria & rel_score_diff(double val)
-    {
-        m_rel_score_diff = val; return *this;
-    }
-
-    double rel_score_diff() const { return m_rel_score_diff; }
-
-    StopCriteria & stop_score(double val)
-    {
-        m_stop_score = val; return *this;
-    }
-
-    double stop_score() const { return m_stop_score; }
-
-    StopCriteria & max_iterations(double val)
-    {
-        m_max_iterations = val; return *this;
-    }
-
-    double max_iterations() const { return m_max_iterations; }
-
-    template<class Fn> StopCriteria & stop_condition(Fn &&cond)
-    {
-        m_stop_condition = cond; return *this;
-    }
-
-    bool stop_condition() { return m_stop_condition(); }
-};
-
-// Helper class to use optimization methods involving gradient.
-template<size_t N> struct ScoreGradient {
-    double score;
-    std::optional<std::array<double, N>> gradient;
-
-    ScoreGradient(double s, const std::array<double, N> &grad)
-        : score{s}, gradient{grad}
-    {}
-};
-
-// Helper to be used in static_assert.
-template<class T> struct always_false { enum { value = false }; };
-
-// Basic interface to optimizer object
-template<class Method, class Enable = void> class Optimizer {
-public:
-
-    Optimizer(const StopCriteria &)
-    {
-        static_assert (always_false<Method>::value,
-                       "Optimizer unimplemented for given method!");
-    }
-
-    Optimizer<Method> &to_min() { return *this; }
-    Optimizer<Method> &to_max() { return *this; }
-    Optimizer<Method> &set_criteria(const StopCriteria &) { return *this; }
-    StopCriteria get_criteria() const { return {}; };
-
-    template<class Func, size_t N>
-    Result<N> optimize(Func&& func,
-                       const Input<N> &initvals,
-                       const Bounds<N>& bounds) { return {}; }
-
-    // optional for randomized methods:
-    void seed(long /*s*/) {}
-};
+namespace Slic3r { namespace opt {
 
 namespace detail {
 
@@ -166,19 +43,6 @@ struct IsNLoptAlg<NLoptAlgComb<a1, a2>> {
 template<class M, class T = void>
 using NLoptOnly = std::enable_if_t<IsNLoptAlg<M>::value, T>;
 
-// Helper to convert C style array to std::array. The copy should be optimized
-// away with modern compilers.
-template<size_t N, class T> auto to_arr(const T *a)
-{
-    std::array<T, N> r;
-    std::copy(a, a + N, std::begin(r));
-    return r;
-}
-
-template<size_t N, class T> auto to_arr(const T (&a) [N])
-{
-    return to_arr<N>(static_cast<const T *>(a));
-}
 
 enum class OptDir { MIN, MAX }; // Where to optimize
 
@@ -357,23 +221,12 @@ public:
     void seed(long s) { m_opt.seed(s); }
 };
 
-template<size_t N> Bounds<N> bounds(const Bound (&b) [N]) { return detail::to_arr(b); }
-template<size_t N> Input<N> initvals(const double (&a) [N]) { return detail::to_arr(a); }
-template<size_t N> auto score_gradient(double s, const double (&grad)[N])
-{
-    return ScoreGradient<N>(s, detail::to_arr(grad));
-}
-
-// Predefinded NLopt algorithms that are used in the codebase
+// Predefinded NLopt algorithms
 using AlgNLoptGenetic = detail::NLoptAlgComb<NLOPT_GN_ESCH>;
 using AlgNLoptSubplex = detail::NLoptAlg<NLOPT_LN_SBPLX>;
 using AlgNLoptSimplex = detail::NLoptAlg<NLOPT_LN_NELDERMEAD>;
-
-// TODO: define others if needed...
-
-// Helper defs for pre-crafted global and local optimizers that work well.
-using DefaultGlobalOptimizer = Optimizer<AlgNLoptGenetic>;
-using DefaultLocalOptimizer  = Optimizer<AlgNLoptSubplex>;
+using AlgNLoptDIRECT  = detail::NLoptAlg<NLOPT_GN_DIRECT>;
+using AlgNLoptMLSL    = detail::NLoptAlg<NLOPT_GN_MLSL>;
 
 }} // namespace Slic3r::opt