[overlay] fix cluster by adapting tolerance

Fixes #1081

doc/release_notes.qbk

@@ -26,6 +26,7 @@
 [*Solved issues]
 
 * [@https://github.com/boostorg/geometry/issues/1076 1076] Union: in rare cases it might miss one polygon
+* [@https://github.com/boostorg/geometry/issues/1081 1081] Union: due to precision it might miss interior rings
 
 [/=================]
 [heading Boost 1.80]

include/boost/geometry/algorithms/detail/overlay/get_clusters.hpp

@@ -35,20 +35,28 @@ namespace detail { namespace overlay
 template <typename Tag = no_rescale_policy_tag, bool Integral = false>
 struct sweep_equal_policy
 {
+private:
+    template <typename T>
+    static inline T threshold()
+    {
+        // Tuned by cases of #1081
+        return T(3);
+    }
+public:
     template <typename P>
     static inline bool equals(P const& p1, P const& p2)
     {
-        // Points within a kilo epsilon are considered as equal
+        // Points within a some epsilons are considered as equal
         using coor_t = typename coordinate_type<P>::type;
-        return approximately_equals(p1, p2, coor_t(1000));
+        return approximately_equals(p1, p2, threshold<coor_t>());
     }
 
     template <typename T>
     static inline bool exceeds(T value)
     {
         // This threshold is an arbitrary value
-        // as long as it is than the used kilo-epsilon
-        T const limit = T(1) / T(1000);
+        // as long as it is bigger than the used value above
+        T const limit = T(1) / threshold<T>();
         return value > limit;
     }
 };

test/algorithms/overlay/get_clusters.cpp

@@ -68,7 +68,7 @@ void do_test(std::string const& case_id,
 }
 
 template <typename Point>
-void test_get_clusters(typename bg::coordinate_type<Point>::type eps)
+void test_get_clusters()
 {
     do_test<Point>("no", {{1.0, 1.0}, {1.0, 2.0}}, 0);
     do_test<Point>("simplex", {{1.0, 1.0}, {1.0, 1.0}}, 1);
@@ -80,8 +80,11 @@ void test_get_clusters(typename bg::coordinate_type<Point>::type eps)
                    6);
     do_test<Point>("buffer3", {{6.41421356237, 5},{6.41421356236, 5},{6.70710678119, 5.29289321881},{6.41421356237, 5},{6, 5},{6.41421356238, 5},{7, 5},{8, 10},{8.41421356237, 10},{8, 9.58578643763},{8.41421356237, 10},{7.41421356237, 9},{7.41421356237, 9},{7, 5.58578643763},{7, 5.58578643763},{6, 5},{6, 5},{6, 5},{6, 5},{6, 5},{6, 6},{4, 6},{4, 6},{3.41421356237, 3},{3, 5},{6, 5},{5, 3},{4, 6},{4, 6},{4, 7},{4, 8},{10.9142135624, 5.5},{8, 5},{10.4142135624, 5},{8, 5},{8, 3.58578643763},{8, 5},{9.41421356237, 7},{9.41421356237, 7},{8.91421356237, 7.5},{10, 7},{8, 9},{7.41421356237, 9},{11, 7}},
                    8);
+}
 
-    // Border cases
+template <typename Point>
+void test_get_clusters_border_cases(typename bg::coordinate_type<Point>::type eps)
+{
     do_test<Point>("borderx_no",  {{1, 1}, {1, 2}, {1 + eps * 10, 1}}, 0);
     do_test<Point>("borderx_yes", {{1, 1}, {1, 2}, {1 + eps, 1}}, 1);
     do_test<Point>("bordery_no",  {{1, 1}, {2, 1}, {1 + eps * 10, 1}}, 0);
@@ -93,8 +96,18 @@ int test_main(int, char* [])
     using fp = bg::model::point<float, 2, bg::cs::cartesian>;
     using dp = bg::model::point<double, 2, bg::cs::cartesian>;
     using ep = bg::model::point<long double, 2, bg::cs::cartesian>;
-    test_get_clusters<fp>(1.0e-4);
-    test_get_clusters<dp>(1.0e-13);
-    test_get_clusters<ep>(1.0e-16);
+
+    test_get_clusters<fp>();
+    test_get_clusters<dp>();
+    test_get_clusters<ep>();
+
+    // This constant relates to the threshold in get_clusters,
+    // which is now return T(3) (earlier it was 1000)
+    double const multiplier = 1.0 / 333.0;
+
+    test_get_clusters_border_cases<fp>(1.0e-4 * multiplier);
+    test_get_clusters_border_cases<dp>(1.0e-13 * multiplier);
+    test_get_clusters_border_cases<ep>(1.0e-16 * multiplier);
+
     return 0;
 }

test/algorithms/overlay/overlay_cases.hpp

@@ -1083,6 +1083,18 @@ static std::string issue_1076[2] =
     "POLYGON((927.152135631899114 0, 981.792858339935151 98, 995 98, 927.152135631899114 0))"
 };
 
+static std::string issue_1081a[2] =
+{
+    "POLYGON((40 100,50 200,75 160,60 160,75 140,40 100),(46.317137830014417 152.14683289288962,59.9312673826029823 167.24190656486104,59.9312673826029823 173.678978034969646,60 175,46.317137830014417 152.14683289288962))",
+    "POLYGON((65 170,67 140,65 163,56.7651907424660749 163.715020193369355,59.9312673826029823 167.241906564861154,65 170))"
+};
+
+static std::string issue_1081b[2] =
+{
+    "POLYGON((21907 0,21917 100,21942 60,21927 60,21942 40,21907 0),(21913.317137830014 52.1468328928896199,21926.9312673826025 67.2419065648610399,21926.9312673826025 73.6789780349696457,21927 75,21913.317137830014 52.1468328928896199))",
+    "POLYGON((21932 70,21934 40,21932 63,21923.7651907424661 63.7150201933693552,21926.9312673826025 67.2419065648611536,21932 70))"
+};
+
 static std::string ggl_list_20120229_volker[3] =
     {
     "POLYGON((1716 1554,2076 2250,2436 2352,2796 1248,3156 2484,3516 2688,3516 2688,3156 2484,2796 1248,2436 2352,2076 2250, 1716 1554))",

test/algorithms/set_operations/difference/difference.cpp

@@ -67,8 +67,9 @@ void test_all()
         1, 5, 8.0);
 
     {
+        // Sym difference works, but expectations are different for rescaling
         ut_settings settings;
-        settings.validity_false_negative_sym = true;
+        settings.sym_difference = false;
         test_one<polygon, polygon, polygon>("star_comb_15",
             star_comb_15[0], star_comb_15[1],
             30, -1, 227.658275102812,

test/algorithms/set_operations/union/union.cpp

@@ -447,6 +447,12 @@ void test_areal()
     TEST_UNION(issue_1076, 1, 0, -1, 1225.0);
     TEST_UNION_REV(issue_1076, 1, 0, -1, 1225.0);
 
+    TEST_UNION(issue_1081a, 1, 2, -1, 1600.56);
+    TEST_UNION_REV(issue_1081a, 1, 2, -1, 1600.56);
+
+    TEST_UNION(issue_1081b, 1, 2, -1, 1600.56);
+    TEST_UNION_REV(issue_1081b, 1, 2, -1, 1600.56);
+
     {
         // Rescaling produces an invalid result
         ut_settings settings;