summaryrefslogtreecommitdiff
path: root/test/geometry
diff options
context:
space:
mode:
Diffstat (limited to 'test/geometry')
-rw-r--r--test/geometry/antipodalPoints.cpp70
-rw-r--r--test/geometry/circle.cpp116
-rw-r--r--test/geometry/closestPair.cpp69
-rw-r--r--test/geometry/closestPair.double.cpp66
-rw-r--r--test/geometry/convexHull.cpp79
-rw-r--r--test/geometry/delaunay.cpp144
-rw-r--r--test/geometry/formulas.cpp127
-rw-r--r--test/geometry/linesAndSegments.cpp240
-rw-r--r--test/geometry/polygon.cpp296
-rw-r--r--test/geometry/segmentIntersection.cpp88
-rw-r--r--test/geometry/sortAround.cpp83
-rw-r--r--test/geometry/triangle.cpp146
12 files changed, 1524 insertions, 0 deletions
diff --git a/test/geometry/antipodalPoints.cpp b/test/geometry/antipodalPoints.cpp
new file mode 100644
index 0000000..d20dfb6
--- /dev/null
+++ b/test/geometry/antipodalPoints.cpp
@@ -0,0 +1,70 @@
+#include "../util.h"
+constexpr ll EPS = 0;
+#define double ll
+#define polar polar<ll>
+#include <geometry/formulas.cpp>
+#undef polar
+#undef double
+#include <geometry/antipodalPoints.cpp>
+#include "../geometry.h"
+
+vector<pair<int, int>> naive(vector<pt> ps) {
+ ll n = sz(ps);
+ auto test = [&](int i, int j){
+ if (dot(ps[j] - ps[i], ps[i - 1] - ps[i]) <= 0) return false;
+ if (dot(ps[j] - ps[i], ps[i + 1] - ps[i]) <= 0) return false;
+ return true;
+ };
+ ps.push_back(ps[0]);
+ ps.push_back(ps[1]);
+ vector<pair<int, int>> res;
+ for (ll i = 1; i <= n; i++) {
+ for (ll j = 1; j < i; j++) {
+ if (test(i, j) && test(j, i)) res.emplace_back(i % n, j % n);
+ }
+ }
+ return res;
+}
+
+void stress_test(ll range) {
+ ll queries = 0;
+ for (int tries = 0; tries < 100'000; tries++) {
+ int n = Random::integer<int>(3, 30);
+ auto ps = Random::convex(n, range);
+
+ auto got = antipodalPoints(ps);
+ for (auto& [a, b] : got) if (a > b) swap(a, b);
+ sort(all(got));
+
+ auto expected = naive(ps);
+ for (auto& [a, b] : expected) if (a > b) swap(a, b);
+
+ for (auto x : expected) {
+ auto it = lower_bound(all(got), x);
+ if (it == got.end() || *it != x) cerr << "error" << FAIL;
+ }
+ queries += n;
+ }
+ cerr << "tested random queries: " << queries << endl;
+}
+
+constexpr int N = 99'000;
+void performance_test() {
+ timer t;
+
+ auto ps = Random::convex(N, 1'000'000'000);
+
+ t.start();
+ auto got = antipodalPoints(ps);
+ t.stop();
+
+ hash_t hash = sz(got);
+ if (t.time > 50) cerr << "too slow: " << t.time << FAIL;
+ cerr << "tested performance: " << t.time << "ms (hash: " << hash << ")" << endl;
+}
+
+int main() {
+ stress_test(100);
+ stress_test(1'000'000'000);
+ performance_test();
+}
diff --git a/test/geometry/circle.cpp b/test/geometry/circle.cpp
new file mode 100644
index 0000000..3d3d27d
--- /dev/null
+++ b/test/geometry/circle.cpp
@@ -0,0 +1,116 @@
+#include "../util.h"
+constexpr double EPS = 1e-6;
+#define ll double
+double gcd(double x, double /**/) {return x;} //hacky
+#include <geometry/formulas.cpp>
+#undef ll
+#include <geometry/circle.cpp>
+
+// Entfernung von Punkt p zur Geraden durch a-b. 2d und 3d
+double distToLine(pt a, pt b, pt p) {
+ return abs(cross(p - a, b - a)) / abs(b - a);
+}
+
+pt randomIntegerPT(ll range) {
+ return pt(Random::integer<ll>(-range, range), Random::integer<ll>(-range, range));
+}
+
+ll sq(ll x) {
+ return x*x;
+}
+
+int expectedCount(ll x1, ll y1, ll r1, ll x2, ll y2, ll r2) {
+ if (x1 == x2 && y1 == y2){
+ return r1 == r2 ? -1 : 0;
+ } else {
+ ll d = sq(x1 - x2) + sq(y1 - y2);
+
+ if (d > sq(r1 + r2) || d < sq(r1 - r2)) {
+ return 0;
+ } else if (d == sq(r1 + r2) || d == sq(r1 - r2)) {
+ return 1;
+ } else{
+ return 2;
+ }
+ }
+}
+
+void test_circleIntersection(ll range) {
+ int queries = 0;
+ for (int tries = 0; tries < 1'000'000; tries++) {
+ auto c1 = randomIntegerPT(range);
+ auto c2 = c1;
+ while (c1 == c2) c2 = randomIntegerPT(range);
+ double r1 = Random::integer<ll>(1, range);
+ double r2 = Random::integer<ll>(1, range);
+
+ auto got = circleIntersection(c1, r1, c2, r2);
+
+ if (sz(got) != expectedCount(real(c1), imag(c1), r1, real(c2), imag(c2), r2)) cerr << "error: wrong count" << FAIL;
+
+ for (int i = 0; i < sz(got); i++) {
+ for (int j = 0; j < i; j++) {
+ if (abs(got[i] - got[j]) < 1e-6) cerr << "error: identical" << FAIL;
+ }
+ }
+
+ for (auto p : got) {
+ if (float_error(abs(c1 - p), r1) > 1e-6) cerr << "error: 1" << FAIL;
+ if (float_error(abs(c2 - p), r2) > 1e-6) cerr << "error: 2" << FAIL;
+ }
+ queries += sz(got);
+ }
+ cerr << "tested circleIntersection: " << queries << endl;
+}
+
+void test_circleRayIntersection(ll range) {
+ int queries = 0;
+ for (int tries = 0; tries < 1'000'000; tries++) {
+ auto c = randomIntegerPT(range);
+ double r = Random::integer<ll>(1, range);
+
+ pt orig = randomIntegerPT(range);
+ pt dir = 0;
+ while (abs(dir) < 0.5) dir = randomIntegerPT(range);
+
+ auto got = circleRayIntersection(c, r, orig, dir);
+
+ double dist = distToLine(orig, orig + dir, c);
+ int lineIntersections = 0;
+ if (dist <= r) lineIntersections = 2;
+ if (abs(dist - r) < 1e-9) lineIntersections = 1;
+
+ int expected = 0;
+ if (abs(orig - c) < r) expected = 1; //starts inside
+ if (abs(orig - c) > r) { //starts outside
+ if (dot(dir, c - orig) >= 0) expected = lineIntersections;
+ else expected = 0;
+ }
+ if (abs(abs(orig - c) - r) < 1e-9) { //starts on circle
+ if (dot(dir, c - orig) >= 0) expected = lineIntersections;
+ else expected = 1;
+ }
+
+ if (sz(got) != expected) cerr << "error: wrong count" << FAIL;
+
+ for (int i = 0; i < sz(got); i++) {
+ for (int j = 0; j < i; j++) {
+ if (abs(got[i] - got[j]) < 1e-6) cerr << "error: identical" << FAIL;
+ }
+ }
+
+ for (auto p : got) {
+ if (float_error(abs(c - p), r) > 1e-6) cerr << "error: 1" << FAIL;
+ if (distToLine(orig, orig + dir, p) > 1e-6) cerr << "error: 2" << FAIL;
+ }
+ queries += sz(got);
+ }
+ cerr << "tested circleIntersection: " << queries << endl;
+}
+
+int main() {
+ test_circleIntersection(10);
+ test_circleIntersection(100);
+ test_circleRayIntersection(10);
+ test_circleRayIntersection(100);
+}
diff --git a/test/geometry/closestPair.cpp b/test/geometry/closestPair.cpp
new file mode 100644
index 0000000..5959b21
--- /dev/null
+++ b/test/geometry/closestPair.cpp
@@ -0,0 +1,69 @@
+#include "../util.h"
+constexpr ll EPS = 0;
+#define double ll
+#define polar polar<ll>
+#include <geometry/formulas.cpp>
+#undef polar
+#undef double
+constexpr ll INF = LL::INF;
+ll sq(ll x) {return x*x;}
+ll isqrt(ll x) {return (ll)sqrtl(x);}
+#include <geometry/closestPair.cpp>
+
+//strict convex hull
+ll naive(const vector<pt>& ps) {
+ ll opt = LL::INF;
+ for (ll i = 0; i < sz(ps); i++) {
+ for (ll j = 0; j < i; j++) {
+ opt = min(opt, norm(ps[i] - ps[j]));
+ }
+ }
+ return opt;
+}
+
+void stress_test(ll range) {
+ ll queries = 0;
+ for (int tries = 0; tries < 100'000; tries++) {
+ int n = Random::integer<int>(2, 100);
+ auto ps = Random::points<ll>(n, -range, range);
+ auto got = shortestDist(ps);
+ auto expected = naive(ps);
+ if (got != expected) cerr << "got: " << got << ", expected: " << expected << FAIL;
+ queries += n;
+ }
+ cerr << "tested random queries: " << queries << endl;
+}
+
+constexpr int N = 1'000'000;
+void performance_test() {
+ timer t;
+ hash_t hash = 0;
+ double maxTime = 0;
+
+ vector<pt> ps;
+ for (int i = 0; i*i <= N; i++) {
+ for (int j = 0; j*j <= N; j++) {
+ ps.emplace_back(i, j);
+ }
+ }
+ t.start();
+ hash = shortestDist(ps);
+ t.stop();
+ maxTime = max(maxTime, t.time);
+
+ ps = Random::points<ll>(N, -1'000'000'000, 1'000'000'000);
+ t.reset();
+ t.start();
+ hash += shortestDist(ps);
+ t.stop();
+ maxTime = max(maxTime, t.time);
+
+ if (maxTime > 500) cerr << "too slow: " << maxTime << FAIL;
+ cerr << "tested performance: " << maxTime << "ms (hash: " << hash << ")" << endl;
+}
+
+int main() {
+ stress_test(100);
+ stress_test(1'000'000'000);
+ performance_test();
+}
diff --git a/test/geometry/closestPair.double.cpp b/test/geometry/closestPair.double.cpp
new file mode 100644
index 0000000..2f8a1ab
--- /dev/null
+++ b/test/geometry/closestPair.double.cpp
@@ -0,0 +1,66 @@
+#include "../util.h"
+constexpr double EPS = 1e-9;
+#define ll double
+double gcd(double x, double /**/) {return x;} //hacky
+#include <geometry/formulas.cpp>
+constexpr ll INF = LL::INF;
+#include <geometry/closestPair.cpp>
+#undef ll
+
+//strict convex hull
+double naive(const vector<pt>& ps) {
+ double opt = LL::INF;
+ for (ll i = 0; i < sz(ps); i++) {
+ for (ll j = 0; j < i; j++) {
+ opt = min(opt, norm(ps[i] - ps[j]));
+ }
+ }
+ return opt;
+}
+
+void stress_test(ll range) {
+ ll queries = 0;
+ for (int tries = 0; tries < 100'000; tries++) {
+ int n = Random::integer<int>(2, 100);
+ auto ps = Random::points<double>(n, -range, range);
+ auto got = shortestDist(ps);
+ auto expected = naive(ps);
+ if (got != expected) cerr << "got: " << got << ", expected: " << expected << FAIL;
+ queries += n;
+ }
+ cerr << "tested random queries: " << queries << endl;
+}
+
+constexpr int N = 1'000'000;
+void performance_test() {
+ timer t;
+ hash_t hash = 0;
+ double maxTime = 0;
+
+ vector<pt> ps;
+ for (int i = 0; i*i <= N; i++) {
+ for (int j = 0; j*j <= N; j++) {
+ ps.emplace_back(i, j);
+ }
+ }
+ t.start();
+ hash = shortestDist(ps);
+ t.stop();
+ maxTime = max(maxTime, t.time);
+
+ ps = Random::points<double>(N, -1'000'000'000, 1'000'000'000);
+ t.reset();
+ t.start();
+ hash += shortestDist(ps);
+ t.stop();
+ maxTime = max(maxTime, t.time);
+
+ if (maxTime > 500) cerr << "too slow: " << maxTime << FAIL;
+ cerr << "tested performance: " << maxTime << "ms (hash: " << hash << ")" << endl;
+}
+
+int main() {
+ stress_test(100);
+ stress_test(1'000'000'000);
+ performance_test();
+}
diff --git a/test/geometry/convexHull.cpp b/test/geometry/convexHull.cpp
new file mode 100644
index 0000000..788a634
--- /dev/null
+++ b/test/geometry/convexHull.cpp
@@ -0,0 +1,79 @@
+#include "../util.h"
+constexpr ll EPS = 0;
+#define double ll
+#define polar polar<ll>
+#include <geometry/formulas.cpp>
+#undef polar
+#undef double
+#include <geometry/convexHull.cpp>
+
+//strict convex hull
+ll isConvexHull(const vector<pt>& ps, const vector<pt>& hull) {
+ ll n = sz(hull) - 1;
+ if (n == 0) {
+ for (pt p : ps) if (p != hull[0]) return 1;
+ return 0;
+ } else {
+ if (hull[0] != hull[n]) return 2;
+ //hull has no duplicates
+ for (ll i = 0; i < n; i++) {
+ for (ll j = 0; j < i; j++) {
+ if (hull[i] == hull[j]) return 3;
+ }
+ }
+ //hull is subset
+ for (pt p : hull) {
+ bool isP = false;
+ for (pt c : ps) isP |= c == p;
+ if (!isP) return 4;
+ }
+ //hull contains all points
+ for (pt p : hull) {
+ ll mi = 1;
+ for (ll i = 0; i < n; i++) {
+ mi = min(mi, cross(hull[i], hull[i + 1], p));
+ }
+ if (mi < 0) return 5; //outside
+ if (mi > 0) continue;
+ bool isCorner = 4;
+ for (pt c : hull) isCorner |= c == p;
+ if (!isCorner) return 6;
+ }
+ // hull is convex
+ if (n <= 2) return 0;
+ for (ll i = 0; i < n; i++) {
+ if (cross(hull[i], hull[i + 1], hull[(i + 2) % n]) <= 0) return 7;
+ }
+ return 0;
+ }
+}
+
+void stress_test(ll range) {
+ ll queries = 0;
+ for (int tries = 0; tries < 100'000; tries++) {
+ int n = Random::integer<int>(1, 100);
+ auto ps = Random::points<ll>(n, -range, range);
+ auto got = convexHull(ps);
+ if (isConvexHull(ps, got) > 0) cerr << "error" << FAIL;
+ queries += n;
+ }
+ cerr << "tested random queries: " << queries << endl;
+}
+
+constexpr int N = 2'000'000;
+void performance_test() {
+ timer t;
+ auto ps = Random::points<ll>(N, -1'000'000'000, 1'000'000'000);
+ t.start();
+ auto a = convexHull(ps);
+ t.stop();
+ hash_t hash = sz(a);
+ if (t.time > 500) cerr << "too slow: " << t.time << FAIL;
+ cerr << "tested performance: " << t.time << "ms (hash: " << hash << ")" << endl;
+}
+
+int main() {
+ stress_test(100);
+ stress_test(1'000'000'000);
+ performance_test();
+}
diff --git a/test/geometry/delaunay.cpp b/test/geometry/delaunay.cpp
new file mode 100644
index 0000000..7f8ec30
--- /dev/null
+++ b/test/geometry/delaunay.cpp
@@ -0,0 +1,144 @@
+#include "../util.h"
+using pt = complex<lll>;
+// Kreuzprodukt, 0, falls kollinear.
+auto cross(pt a, pt b) {return imag(conj(a) * b);}
+auto cross(pt p, pt a, pt b) {return cross(a - p, b - p);}
+#pragma GCC diagnostic ignored "-Wunused-variable"
+#include <geometry/delaunay.cpp>
+
+vector<pt> convexHull(vector<pt> pts){
+ sort(all(pts), [](const pt& a, const pt& b){
+ return real(a) == real(b) ? imag(a) < imag(b)
+ : real(a) < real(b);
+ });
+ pts.erase(unique(all(pts)), pts.end());
+ int k = 0;
+ vector<pt> h(2 * sz(pts));
+ auto half = [&](auto begin, auto end, int t) {
+ for (auto it = begin; it != end; it++) {
+ while (k > t && cross(h[k-2], h[k-1], *it) < 0) k--;//allow collinear points!
+ h[k++] = *it;
+ }};
+ half(all(pts), 1);// Untere Hülle.
+ half(next(pts.rbegin()), pts.rend(), k);// Obere Hülle.
+ h.resize(k);
+ return h;
+}
+
+lll area(const vector<pt>& poly) { //poly[0] == poly.back()
+ lll res = 0;
+ for (int i = 0; i + 1 < sz(poly); i++)
+ res += cross(poly[i], poly[i + 1]);
+ return res;
+}
+
+// Liegt p auf der Strecke a-b?
+bool pointInLineSegment(pt a, pt b, pt p) {
+ if (cross(a, b, p) != 0) return false;
+ auto dist = norm(a - b);
+ return norm(a - p) < dist && norm(b - p) < dist;
+}
+
+// Test auf Streckenschnitt zwischen a-b und c-d.
+// (nur intern)
+bool lineSegmentIntersection(pt a, pt b, pt c, pt d) {
+ if (cross(a, b, c) == 0 && cross(a, b, d) == 0) {
+ return pointInLineSegment(a,b,c) ||
+ pointInLineSegment(a,b,d) ||
+ pointInLineSegment(c,d,a) ||
+ pointInLineSegment(c,d,b);
+ }
+ return cross(a, b, c) * cross(a, b, d) < 0 &&
+ cross(c, d, a) * cross(c, d, b) < 0;
+}
+
+// 1 => c links von a->b
+// 0 => a, b und c kolliniear
+// -1 => c rechts von a->b
+int ccw(pt a, pt b, pt c) {
+ auto orien = cross(b - a, c - a);
+ return (orien > 0) - (orien < 0);
+}
+
+bool inOutCirc(pt a, pt b, pt c, pt p) {
+ lll p2 = norm(p);
+ lll A = norm(a)-p2;
+ lll B = norm(b)-p2;
+ lll C = norm(c)-p2;
+ return ccw(a, b, c) * (cross(p, a, b)*C + cross(p, b, c)*A + cross(p, c, a)*B) > 0;
+}
+
+
+void stress_test(ll range) {
+ ll queries = 0;
+ for (int tries = 0; tries < 100'000; tries++) {
+ int n = Random::integer<int>(3, 30);
+ auto ps = Random::points<lll>(n, -range, range);
+ bool skip = true;
+ for (int i = 2; i < n; i++) skip &= cross(ps[i-2], ps[i-1], ps[i]) == 0;
+ if (skip) continue;
+ for (int i = 0; i < n; i++) {
+ for (int j = 0; j < i; j++) {
+ skip |= ps[i] == ps[j];
+ }
+ }
+ if (skip) continue;
+
+ auto hull = convexHull(ps);
+ lll expectedArea = area(hull);
+ hull.pop_back();
+
+ auto got = delaunay(ps);
+ if (sz(got) % 3 != 0) cerr << "error: not triangles" << FAIL;
+ if (sz(got) / 3 + sz(hull) - 3 + 1 != 2 * sz(ps) - 4) cerr << "error: wrong number" << FAIL;
+
+ //all triangles should be oriented ccw
+ lll gotArea = 0;
+ for (int i = 0; i < sz(got); i += 3) gotArea += cross(got[i], got[i+1], got[i+2]);
+ if (gotArea != expectedArea) cerr << "error: wrong area" << FAIL;
+
+ for (int i = 0; i < sz(got); i++) {
+ int ii = i + 1;
+ if (i / 3 != ii / 3) ii -= 3;
+ for (int j = 0; j < i; j++) {
+ int jj = j + 1;
+ if (j / 3 != jj / 3) jj -= 3;
+
+ if (got[i] == got[j] && got[ii] == got[jj]) cerr << "error: dublicate" << FAIL;
+ if (lineSegmentIntersection(got[i], got[ii], got[j], got[jj])) cerr << "error: intersection" << FAIL;
+ }
+ bool seen = false;
+ for (pt p : ps) seen |= p == got[i];
+ if (!seen) cerr << "error: invalid point" << FAIL;
+ }
+ for (int i = 0; i < sz(got); i += 3) {
+ for (pt p : ps) {
+ if (p == got[i]) continue;
+ if (p == got[i+1]) continue;
+ if (p == got[i+2]) continue;
+ if (inOutCirc(got[i], got[i+1], got[i+2], p)) cerr << "error: not delaunay" << FAIL;
+ }
+ }
+ queries += n;
+ }
+ cerr << "tested random queries: " << queries << endl;
+}
+
+constexpr int N = 100'000;
+void performance_test() {
+ timer t;
+ auto ps = Random::points<lll>(N, -1'000'000'000, 1'000'000'000);
+ t.start();
+ auto got = delaunay(ps);
+ t.stop();
+ hash_t hash = sz(got);
+ if (t.time > 500) cerr << "too slow: " << t.time << FAIL;
+ cerr << "tested performance: " << t.time << "ms (hash: " << hash << ")" << endl;
+}
+
+int main() {
+ stress_test(10);
+ stress_test(10'000);
+ stress_test(1'000'000'000);
+ performance_test();
+}
diff --git a/test/geometry/formulas.cpp b/test/geometry/formulas.cpp
new file mode 100644
index 0000000..d63d431
--- /dev/null
+++ b/test/geometry/formulas.cpp
@@ -0,0 +1,127 @@
+#include "../util.h"
+constexpr ll EPS = 0;
+#define double ll
+#define polar polar<ll>
+#include <geometry/formulas.cpp>
+#undef polar
+#undef double
+
+void test_dot(ll range) {
+ ll queries = 0;
+ for (int tries = 0; tries < 100'000; tries++) {
+ auto p = Random::point<ll>(-range, range);
+ auto q = Random::point<ll>(-range, range);
+
+ ll expected = real(p) * real(q) + imag(p) * imag(q);
+ ll got = dot(p, q);
+
+ if (got != expected) cerr << "got: " << got << ", expected: " << expected << FAIL;
+ queries++;
+ }
+ cerr << "tested dot: " << queries << endl;
+}
+
+void test_norm(ll range) {
+ ll queries = 0;
+ for (int tries = 0; tries < 100'000; tries++) {
+ auto p = Random::point<ll>(-range, range);
+
+ ll expected = real(p) * real(p) + imag(p) * imag(p);
+ ll got = norm(p);
+
+ if (got != expected) cerr << "got: " << got << ", expected: " << expected << FAIL;
+ queries++;
+ }
+ cerr << "tested norm: " << queries << endl;
+}
+
+void test_cross(ll range) {
+ ll queries = 0;
+ for (int tries = 0; tries < 100'000; tries++) {
+ auto p = Random::point<ll>(-range, range);
+ auto q = Random::point<ll>(-range, range);
+
+ ll expected = real(p) * imag(q) - imag(p) * real(q);
+ ll got = cross(p, q);
+
+ if (got != expected) cerr << "got: " << got << ", expected: " << expected << FAIL;
+ queries++;
+ }
+ cerr << "tested cross1: " << queries << endl;
+
+ queries = 0;
+ for (int tries = 0; tries < 100'000; tries++) {
+ auto a = Random::point<ll>(-range, range);
+ auto b = Random::point<ll>(-range, range);
+ auto c = Random::point<ll>(-range, range);
+
+ ll expected = cross(b - a, c - a);
+ ll got = cross(a, b, c);
+
+ if (got != expected) cerr << "got: " << got << ", expected: " << expected << FAIL;
+ queries++;
+ }
+ cerr << "tested cross2: " << queries << endl;
+}
+
+void test_ccw(ll range) {
+ ll queries = 0;
+ for (int tries = 0; tries < 100'000; tries++) {
+ auto a = Random::point<ll>(-range, range);
+ auto b = Random::point<ll>(-range, range);
+ auto c = Random::point<ll>(-range, range);
+
+ ll expected = cross(a, b, c);
+ if (expected < 0) expected = -1;
+ if (expected > 0) expected = 1;
+ ll got = ccw(a, b, c);
+
+ if (got != expected) cerr << "got: " << got << ", expected: " << expected << FAIL;
+ queries++;
+ }
+ cerr << "tested ccw: " << queries << endl;
+}
+
+void test_isCoplanar(ll range) {(void) range;}// cant check this...
+
+void test_uniqueAngle(ll range) {
+ auto lessPt = [](pt a, pt b){
+ if (real(a) != real(b)) return real(a) < real(b);
+ return imag(a) < imag(b);
+ };
+ map<pt, pt, decltype(lessPt)> seen(lessPt);
+
+ ll queries = 0;
+ for (int tries = 0; tries < 100'000; tries++) {
+ pt expected = Random::point<ll>(-sqrt(range), sqrt(range));
+ ll g = abs(gcd(real(expected), imag(expected)));
+ if (g == 0) continue;
+ expected /= g;
+
+ pt rot = Random::point<ll>(-sqrt(range), sqrt(range));
+ if (norm(rot) == 0) continue;
+
+ pt got = uniqueAngle(expected * rot, pt(Random::integer<ll>(1, sqrt(range)), 0) * rot);
+ auto it = seen.emplace(got, expected).first;
+
+ if (it->second != expected) cerr << "error: inconsistent" << FAIL;
+ queries++;
+ }
+ cerr << "tested uniqueAngle: " << queries << " (" << sz(seen) << ")" << endl;
+}
+
+int main() {
+ test_dot(100);
+ test_dot(1'000'000'000);
+ test_norm(100);
+ test_norm(1'000'000'000);
+ test_cross(100);
+ test_cross(1'000'000'000);
+ test_ccw(100);
+ test_ccw(1'000'000'000);
+ test_isCoplanar(100);
+ test_isCoplanar(1'000'000'000);
+ test_uniqueAngle(100);
+ test_uniqueAngle(10'000);
+ test_uniqueAngle(1'000'000'000);
+}
diff --git a/test/geometry/linesAndSegments.cpp b/test/geometry/linesAndSegments.cpp
new file mode 100644
index 0000000..2943a67
--- /dev/null
+++ b/test/geometry/linesAndSegments.cpp
@@ -0,0 +1,240 @@
+#include "../util.h"
+constexpr double EPS = 1e-9;
+#define ll double
+double gcd(double x, double /**/) {return x;} //hacky
+#include <geometry/formulas.cpp>
+#undef ll
+#include <geometry/linesAndSegments.cpp>
+
+#include "../geometry.h"
+
+void stress_pointOnLine(ll range) {
+ ll queries = 0;
+ for (int tries = 0; tries < 100'000; tries++) {
+ auto [a, b] = Random::line(range);
+ pt p = Random::integerPoint(range);
+
+ bool expected = ccw(a, b, p) == 0;
+ bool got = pointOnLine(a, b, p);
+
+ if (got != expected) cerr << "error" << FAIL;
+ queries++;
+ }
+ cerr << "tested pointOnLine: " << queries << endl;
+}
+
+void stress_lineIntersection(ll range) {
+ ll queries = 0;
+ for (int tries = 0; tries < 100'000; tries++) {
+ auto [a, b] = Random::line(range);
+ auto [c, d] = Random::line(range);
+ if (ccw(a, b, c) == 0 && ccw(a, b, d) == 0) continue;
+
+ bool expected = ccw(0, a-b, c-d) == 0;
+ bool got = lineIntersection(a, b, c, d);
+
+ if (got != expected) cerr << "error" << FAIL;
+ queries++;
+ }
+ cerr << "tested lineIntersection: " << queries << endl;
+}
+
+void stress_lineIntersection2(ll range) {
+ ll queries = 0;
+ for (int tries = 0; tries < 100'000; tries++) {
+ auto [a, b] = Random::line(range);
+ auto [c, d] = Random::line(range);
+ if (ccw(0, a-b, c-d) == 0) continue;
+
+ auto got = lineIntersection2(a, b, c, d);
+
+ if (distToLine(a, b, got) > 1e-6) cerr << "error: 1" << FAIL;
+ if (distToLine(a, b, got) > 1e-6) cerr << "error: 2" << FAIL;
+ queries++;
+ }
+ cerr << "tested lineIntersection2: " << queries << endl;
+}
+
+void stress_distToLine(ll range) {
+ ll queries = 0;
+ for (int tries = 0; tries < 100'000; tries++) {
+ auto [a, b] = Random::line(range);
+ pt p = Random::integerPoint(range);
+
+ auto got = distToLine(a, b, p);
+ auto expected = abs(p - projectToLine(a, b, p));
+
+ if (float_error(got, expected) > 1e-6) cerr << "error" << FAIL;
+
+ queries++;
+ }
+ cerr << "tested distToLine: " << queries << endl;
+}
+
+void stress_projectToLine(ll range) {
+ ll queries = 0;
+ for (int tries = 0; tries < 100'000; tries++) {
+ auto [a, b] = Random::line(range);
+ pt p = Random::integerPoint(range);
+
+ auto got = projectToLine(a, b, p);
+
+ if (distToLine(a, b, got) > 1e-6) cerr << "error: 1" << FAIL;
+ if (dot((b-a)/abs(b-a), (got-p)/abs(got-p)) > 1e-6) cerr << "error: 2" << FAIL;
+
+ queries++;
+ }
+ cerr << "tested projectToLine: " << queries << endl;
+}
+
+void stress_sortLine(ll range) {
+ ll queries = 0;
+ for (int tries = 0; tries < 100'000; tries++) {
+ pt dir = 0;
+ while (norm(dir) == 0) dir = Random::integerPoint(range);
+ int n = Random::integer<int>(1, 30);
+ vector<pt> ps = Random::integerPoints(n, range);
+
+ sortLine(dir, ps);
+
+ for (int i = 1; i < n; i++) {
+ if (dot(dir, ps[i-1]) > dot(dir, ps[i])) cerr << "error" << FAIL;
+ }
+ queries+=n;
+ }
+ cerr << "tested sortLine: " << queries << endl;
+}
+
+void stress_pointOnSegment(ll range) {
+ ll queries = 0;
+ for (int tries = 0; tries < 100'000; tries++) {
+ auto [a, b] = Random::line(range);
+ pt p = Random::integerPoint(range);
+
+ bool expected = pointOnLine(a, b, p) && abs(a-p) <= abs(a-b) && abs(b-p) <= abs(a-b);
+ bool got = pointOnSegment(a, b, p);
+
+ if (got != expected) cerr << "error" << FAIL;
+ queries++;
+ }
+ cerr << "tested pointOnSegment: " << queries << endl;
+}
+
+void stress_distToSegment(ll range) {
+ ll queries = 0;
+ for (int tries = 0; tries < 100'000; tries++) {
+ auto [a, b] = Random::line(range);
+ pt p = Random::integerPoint(range);
+
+ double expected = min(abs(p-a), abs(p-b));
+ if (dot(b-a,p-a) >= 0 && dot(a-b,p-b) >= 0) expected = min(expected, distToLine(a, b, p));
+ double got = distToSegment(a, b, p);
+
+ if (float_error(got, expected) > 1e-6) cerr << "error" << FAIL;
+ queries++;
+ }
+ cerr << "tested distToSegment: " << queries << endl;
+}
+
+void stress_segmentIntersection(ll range) {
+ ll queries = 0;
+ for (int tries = 0; tries < 100'000; tries++) {
+ auto [a, b] = Random::line(range);
+ auto [c, d] = Random::line(range);
+
+ bool expected;
+ if (ccw(a, b, c) == 0 && ccw(a, b, d) == 0) {
+ expected = pointOnSegment(a,b,c) ||
+ pointOnSegment(a,b,d) ||
+ pointOnSegment(c,d,a) ||
+ pointOnSegment(c,d,b);
+ } else {
+ expected = ccw(a, b, c) * ccw(a, b, d) <= 0 &&
+ ccw(c, d, a) * ccw(c, d, b) <= 0;
+ }
+ bool got = segmentIntersection(a, b, c, d);
+
+ if (got != expected) cerr << "error" << FAIL;
+ queries++;
+ }
+ cerr << "tested segmentIntersection: " << queries << endl;
+}
+
+void stress_segmentIntersection2(ll range) {
+ ll queries = 0;
+ for (int tries = 0; tries < 100'000; tries++) {
+ auto [a, b] = Random::line(range);
+ auto [c, d] = Random::line(range);
+
+ auto got = segmentIntersection2(a, b, c, d);
+ auto tmp = segmentIntersection(a, b, c, d);
+
+ if (!got.empty() != tmp) cerr << "error: 1" << FAIL;
+ for (pt p : got) {
+ if (distToSegment(a, b, p) > 1e-6) cerr << "error: 2" << FAIL;
+ if (distToSegment(a, b, p) > 1e-6) cerr << "error: 3" << FAIL;
+ }
+ if (tmp) {
+ double gotDist = abs(got.front() - got.back());
+ double expectedDist = 0;
+ array<pt, 4> tmp2 = {a, b, c, d};
+ for (int i = 0; i < 4; i++) {
+ for (int j = 0; j < i; j++) {
+ if (!pointOnSegment(a, b, tmp2[i])) continue;
+ if (!pointOnSegment(c, d, tmp2[i])) continue;
+ if (!pointOnSegment(a, b, tmp2[j])) continue;
+ if (!pointOnSegment(c, d, tmp2[j])) continue;
+ expectedDist = max(expectedDist, abs(tmp2[i] - tmp2[j]));
+ }
+ }
+ if (float_error(gotDist, expectedDist) > 1e-6) cerr << "error: 4" << FAIL;
+ }
+ queries++;
+ }
+ cerr << "tested segmentIntersection2: " << queries << endl;
+}
+
+void stress_distBetweenSegments(ll range) {
+ ll queries = 0;
+ for (int tries = 0; tries < 100'000; tries++) {
+ auto [a, b] = Random::line(range);
+ auto [c, d] = Random::line(range);
+
+ double expected = 0;
+ if (!segmentIntersection(a, b, c, d)) {
+ expected = min({distToSegment(a, b, c), distToSegment(a, b, d),
+ distToSegment(c, d, a), distToSegment(c, d, b)});
+ }
+ double got = distBetweenSegments(a, b, c, d);
+
+ if (float_error(got, expected) > 1e-6) cerr << "error" << FAIL;
+ queries++;
+ }
+ cerr << "tested distBetweenSegments: " << queries << endl;
+}
+
+int main() {
+ stress_pointOnLine(100);
+ stress_pointOnLine(10'000);
+ stress_pointOnLine(1'000'000'000);
+ stress_lineIntersection(100);
+ stress_lineIntersection(1'000'000'000);
+ stress_lineIntersection2(100);
+ stress_lineIntersection2(1'000'000);
+ stress_distToLine(100);
+ stress_distToLine(1'000'000'000);
+ stress_projectToLine(100);
+ stress_projectToLine(1'000'000);
+ stress_sortLine(100);
+ stress_sortLine(1'000'000'000);
+ stress_pointOnSegment(100);
+ stress_pointOnSegment(1'000'000'000);
+ stress_distToSegment(100);
+ stress_distToSegment(1'000'000'000);
+ stress_segmentIntersection(100);
+ stress_segmentIntersection(1'000'000'000);
+ stress_segmentIntersection2(100);
+ stress_segmentIntersection2(1'000'000'000);
+ stress_distBetweenSegments(100);
+ stress_distBetweenSegments(1'000'000'000);
+}
diff --git a/test/geometry/polygon.cpp b/test/geometry/polygon.cpp
new file mode 100644
index 0000000..1dd46ca
--- /dev/null
+++ b/test/geometry/polygon.cpp
@@ -0,0 +1,296 @@
+#include "../util.h"
+constexpr ll EPS = 0;
+constexpr double INF = LD::INF;
+#define double ll
+#define polar polar<ll>
+#include <geometry/formulas.cpp>
+#undef polar
+#undef double
+double abs(pt p) {
+ return hypot(real(p), imag(p));
+}
+// Liegt p auf der Strecke a-b?
+bool pointOnLineSegment(pt a, pt b, pt p) {
+ if (cross(a, b, p) != 0) return false;
+ auto dist = norm(a - b);
+ return norm(a - p) <= dist && norm(b - p) <= dist;
+}
+// Entfernung von Punkt p zur Strecke a-b.
+double distToSegment(pt a, pt b, pt p) {
+ if (a == b) return abs(p - a);
+ if (dot(p - a, b - a) <= 0) return abs(p - a);
+ if (dot(p - b, b - a) >= 0) return abs(p - b);
+ return abs(cross(p - a, b - a)) / abs(b - a);
+}
+#pragma GCC diagnostic ignored "-Wunused-variable"
+#include <geometry/polygon.cpp>
+#include "../geometry.h"
+
+void test_area(ll range) {
+ int queries = 0;
+ for (int tries = 0; tries < 100'000; tries++) {
+ int n = Random::integer(3, 30);
+ auto ps = Random::polygon(n, range);
+ ps.push_back(ps[0]);
+
+ ll expected = 0;
+ for (int i = 0; i < n; i++) {
+ expected += cross(0, ps[i], ps[i+1]);
+ }
+ double got = area(ps) * 2;
+ if (got != expected) cerr << "got: " << got << ", expected: " << expected << FAIL;
+ queries += n;
+ }
+ cerr << "tested area: " << queries << endl;
+}
+
+bool ptLess(pt a, pt b) {
+ if (real(a) != real(b)) return real(a) < real(b);
+ return imag(a) < imag(b);
+}
+
+void test_windingNumber(ll range) {
+ int queries = 0;
+ for (int tries = 0; tries < 100'000; tries++) {
+ int n = Random::integer(3, 8);
+ auto ps = Random::polygon(n, range);
+ ps.push_back(ps[0]);
+
+ for (int i = 0; i < 100; i++) {
+ auto p = Random::point<ll>(-range, range);
+
+ ll expected = 0;
+ bool onBorder = false;
+ for (int j = 0; j < n; j++) {
+ int cur = details::lineSegmentIntersection(p, p + pt(1, 2'000'000'007), ps[j], ps[j+1]);
+ if (ptLess(ps[j], ps[j+1])) expected -= cur;
+ else expected += cur;
+ onBorder |= pointOnLineSegment(ps[j], ps[j+1], p);
+ }
+ if (onBorder) continue;
+ if (area(ps) < 0) expected = -expected;
+
+ bool got = windingNumber(p, ps);
+
+ if (got != expected) cerr << "got: " << got << ", expected: " << expected << FAIL;
+ }
+ queries += n;
+ }
+ cerr << "tested windingNumber: " << queries << endl;
+}
+
+void test_inside(ll range) {
+ int queries = 0;
+ for (int tries = 0; tries < 100'000; tries++) {
+ int n = Random::integer(3, 30);
+ auto ps = Random::polygon(n, range);
+ ps.push_back(ps[0]);
+
+ for (int i = 0; i < 100; i++) {
+ auto p = Random::point<ll>(-range, range);
+
+ ll count = 0;
+ bool onBorder = false;
+ for (int j = 0; j < n; j++) {
+ count += details::lineSegmentIntersection(p, p + pt(1, 2'000'000'007), ps[j], ps[j+1]);
+ onBorder |= pointOnLineSegment(ps[j], ps[j+1], p);
+ }
+ bool expected = (count % 2) && !onBorder;
+ bool got = inside(p, ps);
+
+ if (got != expected) cerr << "error" << FAIL;
+ }
+ queries += n;
+ }
+ cerr << "tested inside: " << queries << endl;
+}
+
+void test_insideConvex(ll range) {
+ int queries = 0;
+ for (int tries = 0; tries < 100'000; tries++) {
+ int n = Random::integer(3, 30);
+ auto ps = Random::convex(n, range);
+
+ for (int i = 0; i < 100; i++) {
+ auto p = Random::point<ll>(-range, range);
+
+ bool expected = true;
+ for (int j = 0; j < n; j++) expected &= cross(p, ps[j], ps[(j+1) % n]) > 0;
+
+ bool got = insideConvex(p, ps);
+
+ if (got != expected) {
+ for (pt pp : ps) cerr << pp << " ";
+ cerr << endl;
+ cerr << p << endl;
+ }
+
+ if (got != expected) cerr << "got: " << got << ", expected: " << expected << FAIL;
+ }
+ queries += n;
+ }
+ cerr << "tested insideConvex: " << queries << endl;
+}
+
+// convex hull without duplicates, h[0] != h.back()
+// apply comments if border counts as inside
+bool insideOrOnConvex(pt p, const vector<pt>& hull) {
+ int l = 0, r = sz(hull) - 1;
+ if (cross(hull[0], hull[r], p) > 0) return false;
+ while (l + 1 < r) {
+ int m = (l + r) / 2;
+ if (cross(hull[0], hull[m], p) >= 0) l = m;
+ else r = m;
+ }
+ return cross(hull[l], hull[r], p) >= 0;
+}
+
+void test_minkowski(ll range) {
+ int queries = 0;
+ for (int tries = 0; tries < 100'000; tries++) {
+ int n = Random::integer(3, 30);
+ auto A = Random::convex(n, range);
+ int m = Random::integer(3, 30);
+ auto B = Random::convex(n, range);
+
+ auto got = minkowski(A, B);
+ bool convex = true;
+ for (int i = 0; i < sz(got); i++) convex &= cross(got[i], got[(i+1) % sz(got)], got[(i+2) % sz(got)]) >= 0;
+ if (!convex) cerr << "error: not convex" << FAIL;
+
+ for (pt a : A) {
+ for (pt b : B) {
+ if (!insideOrOnConvex(a + b, got)) cerr << "error: not sum" << FAIL;
+ }
+ }
+ queries += n + m;
+ }
+ cerr << "tested minkowski: " << queries << endl;
+}
+
+double naive_dist(const vector<pt>& ps, const vector<pt>& qs) {
+ //check if intersect
+ double res = LD::INF;
+ bool intersect = true;
+ for (int i = 0; i < sz(qs); i++) {
+ bool sep = true;
+ for (pt p : ps) {
+ res = min(res, distToSegment(qs[i], qs[(i+1) % sz(qs)], p));
+ sep &= cross(qs[i], qs[(i+1) % sz(qs)], p) <= 0;
+ }
+ if (sep) intersect = false;
+ }
+ for (int i = 0; i < sz(ps); i++) {
+ bool sep = true;
+ for (pt q : qs) {
+ res = min(res, distToSegment(ps[i], ps[(i+1) % sz(ps)], q));
+ sep &= cross(ps[i], ps[(i+1) % sz(ps)], q) <= 0;
+ }
+ if (sep) intersect = false;
+ }
+ if (intersect) return 0;
+ return res;
+}
+
+void test_dist(ll range) {
+ int queries = 0;
+ int pos = 0;
+ for (int tries = 0; tries < 100'000; tries++) {
+ int n = Random::integer(3, 10);
+ auto A = Random::convex(n, range / 3);
+ int m = Random::integer(3, 10);
+ auto B = Random::convex(n, range / 3);
+
+ pt offset = Random::point<ll>(range / 3, 2 * range / 3);
+ for (pt& p : B) p += offset;
+
+ auto got = dist(A, B);
+ auto expected = naive_dist(A, B);
+
+ if (float_error(got, expected) > 1e-6) cerr << "got: " << got << ", expected: " << expected << FAIL;
+ if (got > 0) pos++;
+
+ queries += n + m;
+ }
+ cerr << "tested dist: " << queries << " (" << pos << ")" << endl;
+}
+
+void test_extremal(ll range) {
+ int queries = 0;
+ for (int tries = 0; tries < 100'000; tries++) {
+ int n = Random::integer(3, 30);
+ auto ps = Random::convex(n, range);
+ ps.push_back(ps[0]);
+
+ for (int i = 0; i < 100; i++) {
+ auto dir = Random::point<ll>(-range, range);
+ int tmp = extremal(ps, dir);
+ if (tmp < 0 || tmp >= n) cerr << "error: out of range" << FAIL;
+
+ auto got = ps[tmp];
+ bool extremal = true;
+ for (pt p : ps) extremal &= dot(dir, p) <= dot(dir, got);
+
+ if (!extremal) cerr << "error: not extremal" << FAIL;
+ queries += n;
+ }
+ }
+ cerr << "tested extremal: " << queries << endl;
+}
+
+void test_intersect(ll range) {
+ int queries = 0;
+ for (int tries = 0; tries < 100'000; tries++) {
+ int n = Random::integer(3, 10);
+ auto ps = Random::convex(n, range);
+ ps.push_back(ps[0]);
+
+ for (int i = 0; i < 100; i++) {
+ pt a = Random::point<ll>(-range, range);
+ pt b = a;
+ while (b == a) b = Random::point<ll>(-range, range);
+
+ auto got = intersectLine(ps, a, b);
+
+ vector<int> expected;
+ for (int j = 0; j < n; j++) {
+ if (cross(ps[j], a, b) > 0 && cross(ps[j+1], a, b) < 0) expected.push_back(j);
+ if (cross(ps[j], a, b) < 0 && cross(ps[j+1], a, b) > 0) expected.push_back(j);
+ if (cross(ps[j], a, b) == 0) {
+ if (cross(ps[j+1], a, b) != 0 ||
+ cross(ps[(j+n-1) % n], a, b) != 0) {
+ expected.push_back(j);
+ }
+ }
+ }
+ if (sz(expected) > 1 && expected[0] == expected[1]) expected.pop_back();
+
+ sort(all(got));
+ sort(all(expected));
+
+ if (got != expected) cerr << "error" << FAIL;
+
+ queries += n;
+ }
+ }
+ cerr << "tested intersect: " << queries << endl;
+}
+
+int main() {
+ test_area(100);
+ test_area(1'000'000'000);
+ test_windingNumber(100);
+ test_windingNumber(1'000'000'000);
+ test_inside(100);
+ test_inside(1'000'000'000);
+ test_insideConvex(100);
+ test_insideConvex(1'000'000'000);
+ test_minkowski(100);
+ test_minkowski(500'000'000);
+ test_dist(100);
+ test_dist(1'000'000'000);
+ test_extremal(100);
+ test_extremal(1'000'000'000);
+ test_intersect(100);
+ test_intersect(1'000'000'000);
+}
diff --git a/test/geometry/segmentIntersection.cpp b/test/geometry/segmentIntersection.cpp
new file mode 100644
index 0000000..9862be5
--- /dev/null
+++ b/test/geometry/segmentIntersection.cpp
@@ -0,0 +1,88 @@
+#include "../util.h"
+constexpr ll EPS = 0;
+#define double ll
+#define polar polar<ll>
+#include <geometry/formulas.cpp>
+#undef polar
+#undef double
+
+// Liegt p auf der Strecke a-b?
+bool pointOnLineSegment(pt a, pt b, pt p) {
+ if (cross(a, b, p) != 0) return false;
+ double dist = norm(a - b);
+ return norm(a - p) <= dist && norm(b - p) <= dist;
+}
+
+// Test auf Streckenschnitt zwischen a-b und c-d.
+bool lineSegmentIntersection(pt a, pt b, pt c, pt d) {
+ if (ccw(a, b, c) == 0 && ccw(a, b, d) == 0)
+ return pointOnLineSegment(a,b,c) ||
+ pointOnLineSegment(a,b,d) ||
+ pointOnLineSegment(c,d,a) ||
+ pointOnLineSegment(c,d,b);
+ return ccw(a, b, c) * ccw(a, b, d) <= 0 &&
+ ccw(c, d, a) * ccw(c, d, b) <= 0;
+}
+
+#include <geometry/segmentIntersection.cpp>
+
+vector<seg> randomSegs(int n, ll range) {
+ auto ps = Random::points<ll>(n, -range, range);
+ vector<seg> segs(n);
+ for (int i = 0; i < n; i++) {
+ pt b;
+ do {
+ b = Random::point<ll>(-pow(range, 0.8), pow(range, 0.8));
+ } while(norm(b) == 0);
+ segs[i] = {ps[i], ps[i] + b, i};
+ }
+ return segs;
+}
+
+bool naive(vector<seg>& segs) {
+ for (ll i = 0; i < sz(segs); i++) {
+ for (ll j = 0; j < i; j++) {
+ if (lineSegmentIntersection(segs[i].a, segs[i].b, segs[j].a, segs[j].b)) return true;
+ }
+ }
+ return false;
+}
+
+void stress_test(ll range) {
+ ll queries = 0;
+ ll intersection = 0;
+ ll notIntersection = 0;
+ for (int tries = 0; tries < 100'000; tries++) {
+ int n = Random::integer<int>(2, 100);
+ auto segs = randomSegs(n, range);
+ auto [a, b] = intersect(segs);
+ bool got = a >= 0;
+ if (got != (b >= 0)) cerr << "error: invalid ans" << FAIL;
+ auto expected = naive(segs);
+ if (got != expected) cerr << "error: intersection not found" << FAIL;
+ if (got && !lineSegmentIntersection(segs[a].a, segs[a].b, segs[b].a, segs[b].b)) cerr << "error: no intersection" << FAIL;
+ queries += n;
+ intersection += got;
+ notIntersection += !got;
+ }
+ cerr << "tested random queries: " << queries << "(" << intersection << ":" << notIntersection << ")" << endl;
+}
+
+constexpr int N = 1'000'000;
+void performance_test() {
+ timer t;
+ auto segs = randomSegs(N, 1'000'000'000);
+
+ t.start();
+ hash_t hash = intersect(segs).first;
+ t.stop();
+
+ if (t.time > 500) cerr << "too slow: " << t.time << FAIL;
+ cerr << "tested performance: " << t.time << "ms (hash: " << hash << ")" << endl;
+}
+
+int main() {
+ stress_test(100);
+ stress_test(1'000'000'000);
+ performance_test();
+}
diff --git a/test/geometry/sortAround.cpp b/test/geometry/sortAround.cpp
new file mode 100644
index 0000000..a27edc8
--- /dev/null
+++ b/test/geometry/sortAround.cpp
@@ -0,0 +1,83 @@
+#include "../util.h"
+constexpr ll EPS = 0;
+#define double ll
+#define polar polar<ll>
+#include <geometry/formulas.cpp>
+#undef polar
+#undef double
+#include <geometry/sortAround.cpp>
+
+//expected order:
+//1 8 7
+//2 . 6
+//3 4 5
+void test_tiny() {
+ vector<pt> expected = {
+ {-1, 1},
+ {-1, 0},
+ {-1,-1},
+ { 0,-1},
+ { 1,-1},
+ { 1, 0},
+ { 1, 1},
+ { 0, 1},
+ };
+ auto got = expected;
+ for (int i = 0; i < 100'000; i++) {
+ shuffle(all(got), Random::rng);
+ sortAround(0, got);
+ if (got != expected) cerr << "error" << FAIL;
+ }
+ cerr << "tested tiny" << endl;
+}
+
+void stress_test(ll range) {
+ ll queries = 0;
+ for (int tries = 0; tries < 100'000; tries++) {
+ int n = Random::integer<int>(2, 100);
+ auto ps = Random::points<ll>(n, -range, range);
+
+ auto contains = [&](pt p){
+ for (pt pp : ps) if (pp == p) return true;
+ return false;
+ };
+
+ pt c;
+ do {
+ c = Random::point<ll>(-range, range);
+ } while (contains(c));
+
+ sortAround(c, ps);
+
+ auto isLeft = [&](pt p){return real(p - c) < 0 || (real(p - c) == 0 && imag(p - c) < 0);};
+ auto isCCW = [&](pt a, pt b){return cross(c, a, b) > 0;};
+ if (!is_partitioned(all(ps), isLeft)) cerr << "error 1" << FAIL;
+ auto mid = partition_point(all(ps), isLeft);
+ if (!is_sorted(ps.begin(), mid, isCCW)) cerr << "error 2" << FAIL;
+ if (!is_sorted(mid, ps.end(), isCCW)) cerr << "error 3" << FAIL;
+ queries += n;
+ }
+ cerr << "tested random queries: " << queries << endl;
+}
+
+constexpr int N = 2'000'000;
+void performance_test() {
+ timer t;
+ auto ps = Random::points<ll>(N, -1'000'000'000, 1'000'000'000);
+
+ t.start();
+ sortAround(0, ps);
+ t.stop();
+
+ hash_t hash = 0;
+ for (pt p : ps) hash += real(p) * imag(p);
+ if (t.time > 500) cerr << "too slow: " << t.time << FAIL;
+ cerr << "tested performance: " << t.time << "ms (hash: " << hash << ")" << endl;
+}
+
+int main() {
+ test_tiny();
+ stress_test(100);
+ stress_test(1'000'000'000);
+ performance_test();
+}
diff --git a/test/geometry/triangle.cpp b/test/geometry/triangle.cpp
new file mode 100644
index 0000000..dc620ee
--- /dev/null
+++ b/test/geometry/triangle.cpp
@@ -0,0 +1,146 @@
+#include "../util.h"
+constexpr double EPS = 1e-6;
+#define ll double
+double gcd(double x, double /**/) {return x;} //hacky
+#include <geometry/formulas.cpp>
+#undef ll
+ll sgn(double x) {
+ return (x > EPS) - (x < -EPS);
+}
+#include <geometry/triangle.cpp>
+#include "../geometry.h"
+
+// Entfernung von Punkt p zur Geraden durch a-b. 2d und 3d
+double distToLine(pt a, pt b, pt p) {
+ return abs(cross(p - a, b - a)) / abs(b - a);
+}
+
+void test_centroid(ll range) {
+ int queries = 0;
+ for (int tries = 0; tries < 1'000'000; tries++) {
+ auto [a, b, c] = Random::triangle(range);
+
+ pt center = centroid(a, b, c);
+
+ if (distToLine(2.0*a, c+b, 2.0*center) > 1e-6) cerr << "error: 1" << FAIL;
+ if (distToLine(2.0*b, c+a, 2.0*center) > 1e-6) cerr << "error: 2" << FAIL;
+ if (distToLine(2.0*c, a+b, 2.0*center) > 1e-6) cerr << "error: 3" << FAIL;
+ queries++;
+ }
+ cerr << "tested centroid: " << queries << endl;
+}
+
+void test_area(ll range) {
+ int queries = 0;
+ for (int tries = 0; tries < 1'000'000; tries++) {
+ auto [a, b, c] = Random::triangle(range);
+
+ auto gotA = 2*area(a, b, c);
+ auto gotB = 2*area(abs(a-b), abs(b-c), abs(c-a));
+ auto expected = llround(gotA);
+
+ if (float_error(gotA, expected) > 1e-6) cerr << "error: 1" << FAIL;
+ if (float_error(gotB, expected) > 1e-3) cerr << "error: 2" << FAIL;
+ queries++;
+ }
+ cerr << "tested area: " << queries << endl;
+}
+
+void test_inCenter(ll range) {
+ int queries = 0;
+ for (int tries = 0; tries < 1'000'000; tries++) {
+ auto [a, b, c] = Random::triangle(range);
+
+ pt center = inCenter(a, b, c);
+
+ double da = distToLine(a, b, center);
+ double db = distToLine(b, c, center);
+ double dc = distToLine(c, a, center);
+
+ double avg = (da + db + dc) / 3.0;
+
+ if (float_error(da, avg) > 1e-6) cerr << "error: 1" << FAIL;
+ if (float_error(db, avg) > 1e-6) cerr << "error: 2" << FAIL;
+ if (float_error(dc, avg) > 1e-6) cerr << "error: 3" << FAIL;
+ queries++;
+ }
+ cerr << "tested inCenter: " << queries << endl;
+}
+
+void test_circumCenter(ll range) {
+ int queries = 0;
+ for (int tries = 0; tries < 1'000'000; tries++) {
+ auto [a, b, c] = Random::triangle(range);
+
+ pt center = circumCenter(a, b, c);
+
+ double da = abs(center - a);
+ double db = abs(center - b);
+ double dc = abs(center - c);
+
+ double avg = (da + db + dc) / 3.0;
+
+ if (float_error(da, avg) > 1e-6) cerr << "error: 1" << FAIL;
+ if (float_error(db, avg) > 1e-6) cerr << "error: 2" << FAIL;
+ if (float_error(dc, avg) > 1e-6) cerr << "error: 3" << FAIL;
+ queries++;
+ }
+ cerr << "tested circumCenter: " << queries << endl;
+}
+
+void test_insideOutCenter(ll range) {
+ int queries = 0;
+ for (int tries = 0; tries < 1'000'000; tries++) {
+ auto [a, b, c] = Random::triangle(range);
+ pt p = Random::integerPoint(range);
+
+ pt center = circumCenter(a, b, c);
+
+ double da = abs(center - a);
+ double db = abs(center - b);
+ double dc = abs(center - c);
+ double dp = abs(center - p);
+
+ double avg = (da + db + dc) / 3.0;
+
+ int expected = dp < avg ? 1 : -1;
+ if (float_error(dp, avg) < 1e-9) expected = 0;
+
+ if (insideOutCenter(a, b, c, p) != expected) cerr << "error" << FAIL;
+
+ queries++;
+ }
+ cerr << "tested insideOutCenter: " << queries << endl;
+}
+
+void test_similar(ll range) {
+ int queries = 0;
+ for (int tries = 0; tries < 1'000'000; tries++) {
+ auto [a, b, c] = Random::triangle(sqrt(range));
+ pt rot = Random::integerPoint(sqrt(range));
+ pt add = Random::integerPoint(range);
+
+ pt d = rot * a + add;
+ pt e = rot * b + add;
+ pt f = rot * c + add;
+
+ if (!similar(a, b, c, d, e, f)) cerr << "error" << FAIL;
+ queries++;
+ }
+ cerr << "tested similar: " << queries << endl;
+}
+
+int main() {
+ test_centroid(100);
+ test_centroid(1'000'000'000);
+ test_area(100);
+ test_area(1'000'000'000);
+ test_inCenter(100);
+ test_inCenter(1'000'000'000);
+ test_circumCenter(100);
+ test_circumCenter(1'000'000'000);
+ test_insideOutCenter(100);
+ test_insideOutCenter(1'000'000'000);
+ test_similar(100);
+ test_similar(1'000'000'000);
+}
generated by cgit v1.2.3 (git 2.39.1) at 2025-12-23 21:39:26 +0000