geo-code/m01/angles.cpp

// GEO01 stretch S4: angles from coordinates, and the orientation test.
//
// Build and run (Developer PowerShell for Visual Studio 2026):
//     cl /std:c++latest /O2 /EHsc /nologo angles.cpp && .\angles.exe
// or with g++ or clang++:
//     g++ -std=c++23 -O2 angles.cpp -o angles && ./angles
//
// What it shows: measured in floating point, a triangle's three angles add to 180 degrees only
// approximately, within about 1e-13 degrees. That is evidence for Theorem 1.3, never a proof of it.
// The sign of the cross product you compute along the way is the orientation test, the predicate
// the computational-geometry bridge (B4) is built on.
#include <algorithm>
#include <array>
#include <cmath>
#include <cstdio>
#include <numbers>
#include <random>

struct Vec {
    double x, y;
};
Vec operator-(Vec a, Vec b) { return {a.x - b.x, a.y - b.y}; }
double dot(Vec u, Vec v) { return u.x * v.x + u.y * v.y; }
// The 2D cross product: the signed area of the parallelogram spanned by u and v.
double cross(Vec u, Vec v) { return u.x * v.y - u.y * v.x; }
// Orientation of the path a -> b -> c: positive for a left (counterclockwise) turn,
// negative for a right turn, zero when the three points are collinear.
double orient(Vec a, Vec b, Vec c) { return cross(b - a, c - a); }

// The non-reflex angle at vertex b, between the rays b->a and b->c, in degrees (0 to 180).
double angle_at(Vec a, Vec b, Vec c) {
    // TODO: return the angle at b between rays b->a and b->c, in degrees.
    return 0.0;
}

int main() {
    // 1. Three angles whose values you know.
    const Vec o{0, 0}, e{1, 0}, n{0, 1}, w{-1, 0}, ne{1, 1};
    std::printf("right angle      %8.3f  (want 90)\n", angle_at(e, o, n));
    std::printf("45 degrees       %8.3f  (want 45)\n", angle_at(e, o, ne));
    std::printf("straight angle   %8.3f  (want 180)\n", angle_at(e, o, w));
    std::printf("zero angle       %8.3f  (want 0)\n", angle_at(e, o, Vec{2, 0}));
    std::printf("tiny angle       %.3e  (want about 5.73e-07)\n", angle_at(Vec{1, 1e-8}, o, e));

    // 2. A million random triangles: how far from 180 do the angle sums stray?
    std::mt19937_64 rng(2026);
    std::uniform_real_distribution<double> coord(-1000.0, 1000.0);
    double worst = 0.0;
    for (int i = 0; i < 1'000'000; ++i) {
        const std::array<Vec, 3> t{Vec{coord(rng), coord(rng)}, Vec{coord(rng), coord(rng)}, Vec{coord(rng), coord(rng)}};
        const double sum = angle_at(t[2], t[0], t[1]) + angle_at(t[0], t[1], t[2]) + angle_at(t[1], t[2], t[0]);
        worst = std::max(worst, std::abs(sum - 180.0));
    }
    std::printf("worst |angle sum - 180| over 1e6 random triangles: %.2e degrees\n", worst);

    // 3. The orientation test.
    std::printf("orient(o, e, n) = %+.1f  (left turn)\n", orient(o, e, n));
    std::printf("orient(o, n, e) = %+.1f  (right turn)\n", orient(o, n, e));
    std::printf("orient(o, e, (2,0)) = %+.1f  (collinear)\n", orient(o, e, Vec{2, 0}));
}