code improvements, clustering improvements

derohde · Jan 30, 2023 · 3da00c9 · 3da00c9
1 parent c723e61
commit 3da00c9
Show file tree

Hide file tree

Showing 8 changed files with 184 additions and 115 deletions.
diff --git a/Fred/stabbing.py b/Fred/stabbing.py
@@ -191,10 +191,10 @@ def is_stabbable(balls, old_samples, new_samples, start, end):
                 return False
     return True
 
-def stabbing_path(balls, epsilon=0.1, n_samples=None):
+def stabbing_path(balls, epsilon=0.5, n_samples=None):
     dim = len(balls[0][0])
     if n_samples is None:
-        n_samples = int(10 * 1/epsilon * np.log(len(balls)))
+        n_samples = int(100 * 1/epsilon * np.log(len(balls)))
     old_samples = list()
     new_samples = list()
     segments = list()

diff --git a/README.md b/README.md
@@ -63,6 +63,7 @@ All simplifications are vertex-restricted!
 - signature: `fred.discrete_klcenter(k, l, curves, distances, random_first_center, fast_simplification)` with parameters 
     - `k`: number of centers
     - `l`: maximum complexity of the centers
+    - `local_search`: number of iterations of local search to improve solution, defaults to `0`
     - `consecutive_call`: reuses distances and simplifications already computed in a previous call if `true`, defaults to `false`
     - `random_first_center`: determines if first center is chosen uniformly at random or first curve is used as first center, optional, defaults to true
     - `fast_simplification`: determines whether to use the minimum error simplification or the faster approximate minimum error simplification, defaults to `false`
@@ -241,7 +242,7 @@ curves.add(fred.Curve([[-10.90772872, -37.0858016 ],[-10.90769264, -37.0858147 ]
 curves.add(fred.Curve([[-10.92366833, -37.04470333],[-10.92388667, -37.04459   ],[-10.9243    , -37.04443667],[-10.92942   , -37.04500667],[-10.92949833, -37.04511667],[-10.92964167, -37.045355  ],[-10.92979667, -37.04559833],[-10.92994833, -37.04581333],[-10.93014333, -37.04608333],[-10.930215  , -37.04626   ],[-10.930225  , -37.04630833],[-10.93017167, -37.046335  ],[-10.93017667, -37.04637833],[-10.93021167, -37.046435  ],[-10.93022167, -37.046565  ],[-10.930385  , -37.046725  ],[-10.93042333, -37.04682667],[-10.93045333, -37.04683833],[-10.930475  , -37.04681833],[-10.93047167, -37.04683   ],[-10.93051333, -37.04682833],[-10.93051333, -37.04682833],[-10.930535  , -37.04677833],[-10.93051167, -37.04674667],[-10.930545  , -37.04672333],[-10.93058167, -37.046735  ],[-10.93061   , -37.04673167],[-10.930675  , -37.04672   ],[-10.93091167, -37.04697667],[-10.93120167, -37.04733833],[-10.93131   , -37.04729333],[-10.93144833, -37.047365  ],[-10.93161333, -37.04756167],[-10.93188833, -37.047885  ],[-10.93227333, -37.04835333],[-10.93259   , -37.0487    ],[-10.932965  , -37.04894333],[-10.93307667, -37.04904833],[-10.93308   , -37.04904833],[-10.93318667, -37.04907333],[-10.933335  , -37.04922667],[-10.93362   , -37.04955333],[-10.93400333, -37.04988667],[-10.93964   , -37.05037833],[-10.93996   , -37.05036   ],[-10.94004333, -37.05033167],[-10.94005667, -37.05033333],[-10.94004   , -37.05036   ],[-10.94003   , -37.05035667],[-10.94003   , -37.05035667],[-10.94009   , -37.050395  ],[-10.94013333, -37.05041167],[-10.940015  , -37.050435  ],[-10.93998167, -37.05052333],[-10.939915  , -37.05052167],[-10.94105   , -37.05060833],[-10.94120167, -37.05051667],[-10.941275  , -37.050325  ],[-10.94148667, -37.04996667],[-10.94171833, -37.049545  ],[-10.94205167, -37.049065  ],[-10.942365  , -37.04862167],[-10.94272167, -37.04815667],[-10.943265  , -37.04768167],[-10.94363333, -37.04743667],[-10.94393333, -37.04705833],[-10.94407833, -37.04694   ],[-10.944025  , -37.04700667],[-10.94415667, -37.04692   ],[-10.944685  , -37.04639667],[-10.94502   , -37.046285  ],[-10.945335  , -37.04617167],[-10.94551167, -37.04608667],[-10.94556167, -37.04607167],[-10.94560667, -37.04605833],[-10.94565333, -37.04604667],[-10.94564   , -37.04606833],[-10.945655  , -37.04607   ],[-10.94581333, -37.04604   ],[-10.94615833, -37.04593   ],[-10.94657333, -37.04581833],[-10.94716667, -37.04568833],[-10.94778   , -37.04560833],[-10.94858667, -37.04545833],[-10.949345  , -37.04543   ],[-10.95012   , -37.04549833],[-10.95087   , -37.045625  ],[-10.95162167, -37.045695  ],[-10.95236333, -37.045695  ],[-10.95323   , -37.04548   ],[-10.953905  , -37.045265  ],[-10.95445167, -37.04512333],[-10.95491167, -37.04514167],[-10.95536167, -37.04527333],[-10.95569833, -37.04549   ],[-10.956075  , -37.04541667],[-10.95645   , -37.04527833],[-10.9566    , -37.04522167],[-10.956585  , -37.04525   ],[-10.95676167, -37.04518667],[-10.95703667, -37.045075  ],[-10.95743167, -37.04488667],[-10.95790667, -37.04466333],[-10.95842167, -37.04445   ],[-10.959     , -37.04420167],[-10.95958167, -37.04394   ],[-10.96025667, -37.04360167],[-10.96096   , -37.04327667],[-10.96165333, -37.04298667],[-10.96249667, -37.04261833],[-10.96316   , -37.042405  ],[-10.96383667, -37.04225667],[-10.96469   , -37.04225667],[-10.96537333, -37.04236167],[-10.96593   , -37.0425    ],[-10.96629333, -37.04257167],[-10.96652167, -37.04256167],[-10.96662333, -37.04262833],[-10.96689333, -37.04273833],[-10.96724833, -37.04289833],[-10.967395  , -37.04297333],[-10.96738833, -37.04269   ],[-10.967495  , -37.04238833],[-10.96758   , -37.04203333],[-10.96761667, -37.04176333],[-10.96762   , -37.04172333],[-10.96762   , -37.04172   ],[-10.96760667, -37.04174   ],[-10.96763   , -37.04173833],[-10.96764   , -37.04174667],[-10.96767167, -37.04174333],[-10.96767333, -37.04172833],[-10.96771   , -37.04163167],[-10.96771167, -37.04162833],[-10.96771167, -37.04162167],[-10.96769   , -37.04167   ],[-10.96768833, -37.04166   ],[-10.96767667, -37.04166   ],[-10.96767167, -37.04166333],[-10.96770167, -37.04162167],[-10.96783667, -37.04138333],[-10.96813833, -37.04147667],[-10.96852833, -37.04157   ],[-10.96927   , -37.04167833],[-10.96989833, -37.04186   ],[-10.970665  , -37.04203333],[-10.97124167, -37.04213167],[-10.97162   , -37.042135  ],[-10.97151667, -37.04214667],[-10.97125167, -37.04228167],[-10.971275  , -37.04231833],[-10.97151167, -37.04235333],[-10.972     , -37.04243167],[-10.97252   , -37.04254333],[-10.97287   , -37.04260667],[-10.97314   , -37.04266667],[-10.97335667, -37.04272333],[-10.973515  , -37.04261833],[-10.97361167, -37.04260167],[-10.97369   , -37.0426    ],[-10.97369333, -37.04260167],[-10.97385833, -37.042605  ],[-10.97423167, -37.042715  ],[-10.97472333, -37.04288   ],[-10.97553   , -37.04310833],[-10.97617667, -37.04323667],[-10.97693833, -37.04342   ],[-10.97756167, -37.04356   ],[-10.97825667, -37.043745  ],[-10.97886333, -37.0439    ],[-10.97936   , -37.04401   ],[-10.97961333, -37.04409667],[-10.97974833, -37.044205  ],[-10.97982667, -37.044585  ],[-10.97963333, -37.04517   ],[-10.979385  , -37.04588667],[-10.97914833, -37.04651   ],[-10.97898   , -37.04679667],[-10.97893333, -37.04698667],[-10.97909667, -37.04712833],[-10.97912667, -37.047185  ],[-10.97908667, -37.047215  ],[-10.979085  , -37.04724333]]))
 curves.add(fred.Curve([[-10.89688234, -37.05339968],[-10.89689515, -37.05332973],[-10.89688042, -37.05333247],[-10.89688029, -37.05333484],[-10.89688058, -37.05333512],[-10.89688073, -37.05333528],[-10.8968805 , -37.05333564],[-10.89688045, -37.0533357 ],[-10.89688034, -37.05333584],[-10.89688022, -37.05333665],[-10.89687998, -37.05333665],[-10.8968798 , -37.05333672],[-10.8968771 , -37.05333475],[-10.89686571, -37.05335033],[-10.89686453, -37.05338245],[-10.89690625, -37.053412  ],[-10.89714799, -37.05345149],[-10.89731552, -37.05353814],[-10.89722363, -37.0537264 ],[-10.89712402, -37.05392281],[-10.89700272, -37.05391454],[-10.89669729, -37.05377559],[-10.89637283, -37.05373788],[-10.89628634, -37.05376364],[-10.89625468, -37.05380163],[-10.89626958, -37.05404301],[-10.89630975, -37.05439376],[-10.8963886 , -37.05509832],[-10.89644765, -37.05538955],[-10.8965611 , -37.05575682],[-10.89665346, -37.05613108],[-10.89677345, -37.05672517],[-10.89689319, -37.05728734],[-10.89692781, -37.05739169],[-10.89695023, -37.05753143],[-10.89702285, -37.05769891],[-10.8971346 , -37.057792  ],[-10.89734674, -37.05797843],[-10.89746565, -37.05809327],[-10.89745873, -37.05809087],[-10.89746601, -37.05809045],[-10.89746603, -37.05809043],[-10.89746598, -37.0580906 ],[-10.89746595, -37.05809065],[-10.89750503, -37.05812029],[-10.89776568, -37.05832533],[-10.89815091, -37.05824039],[-10.89867417, -37.05785587],[-10.89905788, -37.05755751],[-10.89965258, -37.05710931],[-10.90031391, -37.05666742],[-10.90072693, -37.05639989],[-10.90117312, -37.0560867 ],[-10.90167986, -37.05570044],[-10.90223155, -37.05532095],[-10.90267435, -37.05501594],[-10.90281205, -37.05492839],[-10.90281977, -37.05492466],[-10.9028111 , -37.05493056],[-10.90281184, -37.05493098],[-10.90281179, -37.05493125],[-10.90281184, -37.05493146],[-10.9028119 , -37.05493155],[-10.90281195, -37.05493173],[-10.9028117 , -37.05493193],[-10.90281944, -37.05492637],[-10.90318115, -37.05466438],[-10.90361849, -37.05433795],[-10.90395364, -37.05420787],[-10.90395573, -37.05453654],[-10.90394111, -37.05505263],[-10.90394185, -37.05587736],[-10.90392791, -37.05664783],[-10.90389758, -37.05750308],[-10.90388737, -37.05786664],[-10.90403   , -37.05799884],[-10.90431907, -37.05799646],[-10.90463524, -37.05800058],[-10.90492265, -37.0580025 ],[-10.90509223, -37.0580149 ],[-10.90513417, -37.05805725],[-10.90519558, -37.05826165],[-10.90527994, -37.05877507],[-10.90535515, -37.05906472],[-10.90619649, -37.05905899],[-10.90641173, -37.06018562],[-10.90644041, -37.06067098],[-10.90642897, -37.06086107],[-10.90641925, -37.06087587],[-10.90642255, -37.06088143],[-10.90642329, -37.06088235],[-10.90642264, -37.06088388],[-10.90642211, -37.06090121],[-10.90642028, -37.06112732],[-10.90649051, -37.06164138],[-10.90728268, -37.06183192],[-10.90753491, -37.06201037],[-10.90761894, -37.06208487],[-10.90763191, -37.06209038],[-10.90764293, -37.06208956],[-10.90764576, -37.06209575],[-10.90764926, -37.06209673],[-10.90764498, -37.06210186],[-10.90764444, -37.0621018 ],[-10.9076441 , -37.06210176],[-10.90764407, -37.06210173],[-10.90791968, -37.06228879],[-10.90840258, -37.06261833],[-10.90885022, -37.06293506],[-10.90941396, -37.06335066],[-10.91002845, -37.06379331],[-10.91021656, -37.06385412],[-10.91021107, -37.06382905],[-10.91020892, -37.0638263 ],[-10.91020903, -37.06382629],[-10.91020921, -37.06382624],[-10.91020941, -37.06382628],[-10.91020941, -37.06382636],[-10.91020942, -37.0638264 ],[-10.91020943, -37.06382641],[-10.91021067, -37.06382699],[-10.91051354, -37.06396407],[-10.91086507, -37.06407459],[-10.91112765, -37.06418145],[-10.91112402, -37.06418525],[-10.91112739, -37.064179  ],[-10.91113069, -37.06417743],[-10.91113068, -37.06417746],[-10.91113069, -37.06417745],[-10.91113068, -37.06417741],[-10.91113062, -37.06417743],[-10.91113027, -37.06417782],[-10.91113022, -37.06417796],[-10.91113014, -37.06417815],[-10.91113018, -37.06417815],[-10.9111301 , -37.06417819],[-10.91113003, -37.06417818],[-10.91113018, -37.06417821],[-10.91121402, -37.06421624],[-10.91152921, -37.06426605],[-10.91195342, -37.06430617],[-10.91243864, -37.06431631],[-10.91301115, -37.06431624],[-10.91373203, -37.06430726],[-10.9143199 , -37.06424356],[-10.91526178, -37.06407409],[-10.91568718, -37.06396253],[-10.91601655, -37.06391468],[-10.91614016, -37.06400304],[-10.91613907, -37.06426894],[-10.91611837, -37.0645705 ],[-10.91620218, -37.06472138],[-10.91657107, -37.0648288 ],[-10.9170171 , -37.06493571],[-10.91728016, -37.06502682],[-10.91734995, -37.06509569],[-10.91735968, -37.06534346],[-10.91731253, -37.06557038],[-10.91730704, -37.06556575],[-10.91730599, -37.06555167]]))
 
-clustering = fred.discrete_klcenter(6, 8, curves, fast_simplification=True)
+clustering = fred.discrete_klcenter(6, 8, curves, local_search = 10, fast_simplification=True)
 clustering.optimize_centers(curves)
 
 fred.plot_clustering(clustering, curves)

diff --git a/include/clustering.hpp b/include/clustering.hpp
@@ -104,9 +104,20 @@ inline distance_t _center_cost_sum(const Curves &in, const Curves &simplified_in
     return cost;
 }
 
-Clustering_Result kl_cluster(const curve_number_t, const curve_size_t, const Curves &, const bool, const bool, const bool, const bool);
+inline distance_t _center_cost_max(const Curves &in, const Curves &simplified_in, const Curve_Numbers &centers, Distance_Matrix &distances) {
+    distance_t cost = 0;
+
+    // for all curves
+    for (curve_number_t i = 0; i < in.size(); ++i) {
+        const auto min_cost_elem = _curve_cost(i, in, simplified_in, centers, distances);
+        cost = std::max(cost, min_cost_elem);
+    }
+    return cost;
+}
+
+Clustering_Result kl_cluster(const curve_number_t, const curve_size_t, const Curves &, unsigned int, const bool, const bool, const bool, const bool);
 
-Clustering_Result kl_center(const curve_number_t, const curve_size_t, const Curves &, const bool = false, const bool = true, const bool = false);
+Clustering_Result kl_center(const curve_number_t, const curve_size_t, const Curves &, unsigned int, const bool = false, const bool = true, const bool = false);
 
 Clustering_Result kl_median(const curve_number_t, const curve_size_t, const Curves &, const bool = false, const bool = false); 
 

diff --git a/include/simplification.hpp b/include/simplification.hpp
@@ -35,94 +35,9 @@ class Subcurve_Shortcut_Graph {
 
 public:
 
-    Subcurve_Shortcut_Graph(Curve &curve) : curve{curve}, edges{std::vector<Distances>(curve.complexity(), Distances(curve.complexity(), std::numeric_limits<distance_t>::infinity()))} {
-        if (Config::verbosity > 1) py::print("SIMPL: computing shortcut graph");
-        const curve_size_t complexity = curve.complexity();
-        Curve segment(2, curve.front().dimensions());
-        auto distance = Frechet::Continuous::Distance();
-
-        for (curve_size_t i = 0; i < complexity - 1; ++i) {
-            for (curve_size_t j = i + 1; j < complexity; ++j) {
-
-                if (Config::verbosity > 1) py::print("SIMPL: computing shortcut distance from vertex ", i, " to vertex ", j);
-
-                curve.set_subcurve(i, j);
-
-                segment[0] = curve.front();
-                segment[1] = curve.back();
-
-                distance = Frechet::Continuous::distance(curve, segment);
-                edges[i][j] = distance.value;
-
-                curve.reset_subcurve();
-            }
-
-        }
-    }
-
-    Curve minimum_error_simplification(const curve_size_t ll) const {
-        if (Config::verbosity > 1) py::print("SIMPL: computing exact minimum error simplification using shortcut graph");
-        if (ll >= curve.complexity()) return curve;
-
-        const curve_size_t l = ll - 1;
-
-        Curve result(curve.dimensions());
-
-        if (ll <= 2) {
-            result.push_back(curve.front());
-            result.push_back(curve.back());
-            return result;
-        }
-
-        std::vector<Distances> distances(curve.complexity(), Distances(l, std::numeric_limits<distance_t>::infinity()));
-        std::vector<std::vector<curve_size_t>> predecessors(curve.complexity(), std::vector<curve_size_t>(l));
-
-        Distances others;
-        curve_size_t best = 0;
-        for (curve_size_t i = 0; i < l; ++i) {
-
-            if (i == 0) {
-                if (Config::verbosity > 1) py::print("SIMPL: initializing arrays");
-                #pragma omp parallel for
-                for (curve_size_t j = 1; j < curve.complexity(); ++j) {
-                    distances[j][0] = edges[0][j];
-                    predecessors[j][0] = 0;
-                }
-            } else {
-                for (curve_size_t j = 1; j < curve.complexity(); ++j) {
-                    if (Config::verbosity > 1) py::print("SIMPL: computing shortcut using ", i, " jumps");
-                    others.resize(j);
-                    #pragma omp parallel for
-                    for (curve_size_t k = 0; k < j; ++k) {
-                        others[k] = std::max(distances[k][i - 1], edges[k][j]);
-                    }
-
-                    best = std::distance(others.begin(), std::min_element(others.begin(), others.end()));
-
-                    distances[j][i] = others[best];
-                    predecessors[j][i] = best;
-                }
-            }
-        }
-
-        if (Config::verbosity > 1) py::print("SIMPL: backwards constructing simplification");
-
-        curve_size_t ell = l - 1;
-
-        result.push_back(curve.back());
-        curve_size_t predecessor = predecessors[curve.complexity() - 1][ell];
-
-        for (curve_size_t i = 0; i < l - 1; ++i) {
-            result.push_back(curve[predecessor]);
-            predecessor = predecessors[predecessor][--ell];
-        }
-
-        result.push_back(curve.front());
-
-        std::reverse(result.begin(), result.end());
-        return result;
-    }
+    Subcurve_Shortcut_Graph(const Curve&);
 
+    Curve minimum_error_simplification(const curve_size_t) const;
 };
 
 Curve approximate_minimum_link_simplification(const Curve&, const distance_t);

diff --git a/setup.py b/setup.py
@@ -80,7 +80,7 @@ def build_extension(self, ext):
 
 setup(
     name='Fred-Frechet',
-    version='1.10.9',
+    version='1.10.10',
     author='Dennis Rohde',
     author_email='[email protected]',
     description='A fast, scalable and light-weight C++ Fréchet distance library, exposed to python and focused on (k,l)-clustering of polygonal curves.',