diff --git a/src/freeroute/geometry/float_point.py b/src/freeroute/geometry/float_point.py new file mode 100644 index 0000000..31c0028 --- /dev/null +++ b/src/freeroute/geometry/float_point.py @@ -0,0 +1,159 @@ +"""``FloatPoint`` — approximate (double-precision) point. + +Ports ``geometry/planar/FloatPoint.java``. + +**Arithmetic model: approximate (Python ``float``).** FreeRouting deliberately +keeps ``FloatPoint`` *outside* the exact ``Point`` hierarchy because float math +is inexact; it is used for distances, rounding, and speed-over-accuracy paths. +Its :meth:`side_of` has no reliable ``COLLINEAR`` result (numerical noise), so +exact side tests must go through :class:`~freeroute.geometry.point.Point`. +""" + +from __future__ import annotations + +import math + +from .side import Side + +__all__ = ["FloatPoint"] + + +class FloatPoint: + """A point in the plane with ``float`` coordinates.""" + + __slots__ = ("x", "y") + + def __init__(self, x: float, y: float) -> None: + self.x = float(x) + self.y = float(y) + + def __eq__(self, other: object) -> bool: + if not isinstance(other, FloatPoint): + return NotImplemented + return self.x == other.x and self.y == other.y + + def __hash__(self) -> int: + return hash((self.x, self.y)) + + def __repr__(self) -> str: + return f"FloatPoint({self.x!r}, {self.y!r})" + + # --- magnitudes --------------------------------------------------------- + + def size_square(self) -> float: + return self.x * self.x + self.y * self.y + + def size(self) -> float: + return math.sqrt(self.size_square()) + + def distance_square(self, other: FloatPoint) -> float: + dx = other.x - self.x + dy = other.y - self.y + return dx * dx + dy * dy + + def distance(self, other: FloatPoint) -> float: + return math.sqrt(self.distance_square(other)) + + def weighted_distance( + self, other: FloatPoint, horizontal_weight: float, vertical_weight: float + ) -> float: + dx = (self.x - other.x) * horizontal_weight + dy = (self.y - other.y) * vertical_weight + return math.sqrt(dx * dx + dy * dy) + + # --- arithmetic --------------------------------------------------------- + + def add(self, other: FloatPoint) -> FloatPoint: + return FloatPoint(self.x + other.x, self.y + other.y) + + def substract(self, other: FloatPoint) -> FloatPoint: + # Spelling matches the upstream method name ``substract``. + return FloatPoint(self.x - other.x, self.y - other.y) + + def middle_point(self, other: FloatPoint) -> FloatPoint: + return FloatPoint(0.5 * (self.x + other.x), 0.5 * (self.y + other.y)) + + def change_size(self, new_size: float) -> FloatPoint: + """A point on the ray from zero through this one, at distance ``new_size``.""" + if self.x == 0 and self.y == 0: + return self + length = math.sqrt(self.x * self.x + self.y * self.y) + return FloatPoint(self.x * new_size / length, self.y * new_size / length) + + def change_length(self, to_point: FloatPoint, new_length: float) -> FloatPoint: + dx = to_point.x - self.x + dy = to_point.y - self.y + if dx == 0 and dy == 0: + return to_point + length = math.sqrt(dx * dx + dy * dy) + return FloatPoint(self.x + dx * new_length / length, self.y + dy * new_length / length) + + def scalar_product(self, p1: FloatPoint, p2: FloatPoint) -> float: + """Scalar product of ``(p1 - self)`` and ``(p2 - self)``.""" + return (p1.x - self.x) * (p2.x - self.x) + (p1.y - self.y) * (p2.y - self.y) + + # --- orientation (approximate) ----------------------------------------- + + def side_of(self, p1: FloatPoint, p2: FloatPoint) -> Side: + """Approximate side of the directed line ``p1 -> p2``. + + ``COLLINEAR`` is unreliable here (float noise); use the exact + ``Point.side_of`` when collinearity matters. + """ + d21x = p2.x - p1.x + d21y = p2.y - p1.y + d01x = self.x - p1.x + d01y = self.y - p1.y + return Side.of(d21x * d01y - d21y * d01x) + + def is_contained_in_box(self, p1: FloatPoint, p2: FloatPoint, tolerance: float) -> bool: + min_x, max_x = (p1.x, p2.x) if p1.x < p2.x else (p2.x, p1.x) + if self.x < min_x - tolerance or self.x > max_x + tolerance: + return False + min_y, max_y = (p1.y, p2.y) if p1.y < p2.y else (p2.y, p1.y) + return min_y - tolerance <= self.y <= max_y + tolerance + + # --- rotation ----------------------------------------------------------- + + def turn_90_degree(self, factor: int) -> FloatPoint: + n = factor % 4 + if n == 0: + return FloatPoint(self.x, self.y) + if n == 1: + return FloatPoint(-self.y, self.x) + if n == 2: + return FloatPoint(-self.x, -self.y) + return FloatPoint(self.y, -self.x) + + def rotate(self, angle: float, pole: FloatPoint) -> FloatPoint: + if angle == 0: + return self + dx = self.x - pole.x + dy = self.y - pole.y + s = math.sin(angle) + c = math.cos(angle) + return FloatPoint(pole.x + dx * c - dy * s, pole.y + dx * s + dy * c) + + # --- rounding to integer points ---------------------------------------- + + def round(self): + """Round to the nearest :class:`~freeroute.geometry.point.IntPoint`.""" + from .point import IntPoint + + return IntPoint(_round_half_up(self.x), _round_half_up(self.y)) + + def bounding_box(self): + """Smallest :class:`~freeroute.geometry.box.IntBox` containing this point.""" + from .box import IntBox + + return IntBox( + math.floor(self.x), + math.floor(self.y), + math.ceil(self.x), + math.ceil(self.y), + ) + + +def _round_half_up(value: float) -> int: + """Round half up, matching Java ``Math.round`` (not Python banker's rounding).""" + return math.floor(value + 0.5) diff --git a/src/freeroute/geometry/limits.py b/src/freeroute/geometry/limits.py new file mode 100644 index 0000000..faf7ad3 --- /dev/null +++ b/src/freeroute/geometry/limits.py @@ -0,0 +1,19 @@ +"""Numeric limits from FreeRouting's ``geometry/planar/Limits.java``. + +In the Java source these bounds exist because coordinate products must fit in a +``long`` (or the mantissa of a ``double``) for the exact-orientation and +intersection tests. In this Python port the exact tests use unbounded ``int`` +arithmetic, so ``CRIT_INT`` is retained only to mirror the upstream promotion +rule from ``IntPoint`` to ``RationalPoint`` — it is not needed for correctness. +""" + +from __future__ import annotations + +#: 2^25 — upper bound such that the product of two values with absolute value at +#: most this fits in a double mantissa with room for one addition (Java only). +CRIT_INT: int = 33554432 + +#: 2^53 — largest double with all smaller integers exactly representable. +CRIT_DOUBLE: float = 9007199254740992.0 + +__all__ = ["CRIT_INT", "CRIT_DOUBLE"] diff --git a/src/freeroute/geometry/point.py b/src/freeroute/geometry/point.py new file mode 100644 index 0000000..4c9982e --- /dev/null +++ b/src/freeroute/geometry/point.py @@ -0,0 +1,236 @@ +"""Points: ``Point`` base, ``IntPoint``, ``RationalPoint``, and factories. + +Ports ``geometry/planar/{Point,IntPoint,RationalPoint}.java``. + +**Arithmetic model:** + +* :class:`IntPoint` — exact, Python ``int`` ``(x, y)``. +* :class:`RationalPoint` — exact, projective ``(x, y, z)`` Python ``int`` + representing the affine point ``(x/z, y/z)``. ``z == 0`` is a point at + infinity (the intersection of parallel lines). ``z`` is kept non-negative. + +Two points of different concrete type are compared/combined by promoting the +integer one to the rational representation, exactly as the Java double-dispatch +does — but with unbounded ``int`` instead of ``BigInteger``. +""" + +from __future__ import annotations + +import math + +from .float_point import FloatPoint +from .side import Side +from .vector import IntVector, RationalVector, Vector, _add_rational + +__all__ = ["Point", "IntPoint", "RationalPoint", "point", "rational_point", "ZERO_POINT"] + + +def _xyz(p: Point) -> tuple[int, int, int]: + if isinstance(p, IntPoint): + return (p.x, p.y, 1) + return (p.x, p.y, p.z) # RationalPoint + + +class Point: + """Abstract base for exact plane points.""" + + def to_float(self) -> FloatPoint: # pragma: no cover - overridden + raise NotImplementedError + + def is_infinite(self) -> bool: # pragma: no cover - overridden + raise NotImplementedError + + def translate_by(self, vector: Vector) -> Point: + if vector.is_zero(): + return self + return vector.add_to_point(self) + + def difference_by(self, other: Point) -> Vector: + """Return the vector ``self - other`` (exact).""" + if isinstance(self, IntPoint) and isinstance(other, IntPoint): + return IntVector(self.x - other.x, self.y - other.y) + sx, sy, sz = _xyz(self) + ox, oy, oz = _xyz(other) + x, y, z = _add_rational((sx, sy, sz), (-ox, -oy, oz)) + return RationalVector(x, y, z) + + def side_of(self, p1: Point, p2: Point) -> Side: + """Side of ``self`` relative to the directed line ``p1 -> p2``.""" + v1 = self.difference_by(p1) + v2 = p2.difference_by(p1) + return v1.side_of(v2) + + def side_of_line(self, line) -> Side: + return self.side_of(line.a, line.b) + + def compare_x(self, other: Point) -> int: + sx, _, sz = _xyz(self) + ox, _, oz = _xyz(other) + return _sign(sx * oz - ox * sz) + + def compare_y(self, other: Point) -> int: + _, sy, sz = _xyz(self) + _, oy, oz = _xyz(other) + return _sign(sy * oz - oy * sz) + + def compare_x_y(self, other: Point) -> int: + result = self.compare_x(other) + if result == 0: + result = self.compare_y(other) + return result + + def turn_90_degree(self, factor: int, pole: Point) -> Point: + v = self.difference_by(pole).turn_90_degree(factor) + return pole.translate_by(v) + + def mirror_vertical(self, pole: Point) -> Point: + v = self.difference_by(pole).mirror_at_y_axis() + return pole.translate_by(v) + + def mirror_horizontal(self, pole: Point) -> Point: + v = self.difference_by(pole).mirror_at_x_axis() + return pole.translate_by(v) + + +class IntPoint(Point): + """A point with exact integer coordinates.""" + + __slots__ = ("x", "y") + + def __init__(self, x: int, y: int) -> None: + self.x = int(x) + self.y = int(y) + + def __eq__(self, other: object) -> bool: + if not isinstance(other, IntPoint): + return NotImplemented + return self.x == other.x and self.y == other.y + + def __hash__(self) -> int: + return hash((self.x, self.y)) + + def __repr__(self) -> str: + return f"IntPoint({self.x}, {self.y})" + + def to_float(self) -> FloatPoint: + return FloatPoint(self.x, self.y) + + def is_infinite(self) -> bool: + return False + + def get_id_no(self) -> int: + return 31 * self.x + self.y + + def determinant(self, other: IntPoint) -> int: + """Exact determinant of the position vectors of ``self`` and ``other``.""" + return self.x * other.y - self.y * other.x + + def distance_square(self, other: IntPoint) -> float: + dx = other.x - self.x + dy = other.y - self.y + return float(dx * dx + dy * dy) + + def distance(self, other: IntPoint) -> float: + return math.sqrt(self.distance_square(other)) + + def surrounding_box(self): + from .box import IntBox + + return IntBox(self.x, self.y, self.x, self.y) + + def is_contained_in(self, box) -> bool: + return box.ll.x <= self.x <= box.ur.x and box.ll.y <= self.y <= box.ur.y + + +class RationalPoint(Point): + """A point with exact rational coordinates ``(x/z, y/z)``, ``z`` >= 0.""" + + __slots__ = ("x", "y", "z") + + def __init__(self, x: int, y: int, z: int) -> None: + if z < 0: + raise ValueError("RationalPoint: z is expected to be >= 0") + self.x = int(x) + self.y = int(y) + self.z = int(z) + + @staticmethod + def from_int_point(p: IntPoint) -> RationalPoint: + return RationalPoint(p.x, p.y, 1) + + def __eq__(self, other: object) -> bool: + if not isinstance(other, RationalPoint): + return NotImplemented + return self.x * other.z - other.x * self.z == 0 and self.y * other.z - other.y * self.z == 0 + + def __hash__(self) -> int: + if self.z != 0: + g = math.gcd(math.gcd(abs(self.x), abs(self.y)), self.z) or 1 + return hash((self.x // g, self.y // g, self.z // g)) + return hash((self.x, self.y, 0)) + + def __repr__(self) -> str: + return f"RationalPoint({self.x}, {self.y}, {self.z})" + + def to_float(self) -> FloatPoint: + if self.z == 0: + big = 3.4e38 + return FloatPoint(big, big) + return FloatPoint(self.x / self.z, self.y / self.z) + + def is_infinite(self) -> bool: + return self.z == 0 + + def is_contained_in(self, box) -> bool: + # Exact: compare numerators against z * boundary (z >= 0). + if self.x < box.ll.x * self.z: + return False + if self.y < box.ll.y * self.z: + return False + if self.x > box.ur.x * self.z: + return False + return self.y <= box.ur.y * self.z + + def surrounding_box(self): + from .box import IntBox + + fp = self.to_float() + return IntBox(math.floor(fp.x), math.floor(fp.y), math.ceil(fp.x), math.ceil(fp.y)) + + +ZERO_POINT = IntPoint(0, 0) + + +# --- factories (Point.get_instance) ----------------------------------------- + + +def point(x: int, y: int) -> IntPoint: + """Create an :class:`IntPoint`. + + FreeRouting promotes to ``RationalPoint`` past ``CRIT_INT`` to avoid ``long`` + overflow; unbounded Python ``int`` makes that unnecessary, so an integer + coordinate always yields an :class:`IntPoint`. + """ + return IntPoint(x, y) + + +def rational_point(x: int, y: int, z: int) -> Point: + """Create a point from projective coordinates, reducing to integer when exact. + + Mirrors ``Point.get_instance(BigInteger, BigInteger, BigInteger)``: a + negative denominator is normalized, and if ``z`` divides both numerators the + result collapses to an :class:`IntPoint`. + """ + if z < 0: + x, y, z = -x, -y, -z + if z != 0 and x % z == 0 and y % z == 0: + return IntPoint(x // z, y // z) + return RationalPoint(x, y, z) + + +def _sign(value: int) -> int: + if value > 0: + return 1 + if value < 0: + return -1 + return 0 diff --git a/src/freeroute/geometry/side.py b/src/freeroute/geometry/side.py new file mode 100644 index 0000000..cd0b050 --- /dev/null +++ b/src/freeroute/geometry/side.py @@ -0,0 +1,56 @@ +"""``Side`` and ``Signum`` — the sign-valued results of exact geometric tests. + +Ports ``geometry/planar/Side.java`` and ``datastructures/Signum.java``. Both are +three-valued signs; representing them as ``IntEnum`` with values ``+1 / -1 / 0`` +makes :meth:`Side.of` / :meth:`Signum.of` and :meth:`negate` trivial and keeps +the arithmetic exact (the input is always an exact integer determinant or +scalar product in this port, never a rounded double). +""" + +from __future__ import annotations + +from enum import IntEnum + +__all__ = ["Side", "Signum"] + + +class Side(IntEnum): + """Which side of a directed line/vector a point lies on. + + ``ON_THE_LEFT`` for a positive determinant, ``ON_THE_RIGHT`` for negative, + ``COLLINEAR`` for zero — matching ``Side.of`` in the Java source. + """ + + ON_THE_LEFT = 1 + ON_THE_RIGHT = -1 + COLLINEAR = 0 + + @staticmethod + def of(value: int | float) -> Side: + if value > 0: + return Side.ON_THE_LEFT + if value < 0: + return Side.ON_THE_RIGHT + return Side.COLLINEAR + + def negate(self) -> Side: + return Side(-self.value) + + +class Signum(IntEnum): + """Sign of a scalar quantity: ``POSITIVE`` / ``NEGATIVE`` / ``ZERO``.""" + + POSITIVE = 1 + NEGATIVE = -1 + ZERO = 0 + + @staticmethod + def of(value: int | float) -> Signum: + if value > 0: + return Signum.POSITIVE + if value < 0: + return Signum.NEGATIVE + return Signum.ZERO + + def negate(self) -> Signum: + return Signum(-self.value) diff --git a/src/freeroute/geometry/vector.py b/src/freeroute/geometry/vector.py new file mode 100644 index 0000000..507b585 --- /dev/null +++ b/src/freeroute/geometry/vector.py @@ -0,0 +1,274 @@ +"""Vectors: ``Vector`` base, ``IntVector``, ``RationalVector``. + +Ports ``geometry/planar/{Vector,IntVector,RationalVector}.java``. + +**Arithmetic model:** + +* :class:`IntVector` — exact, Python ``int`` ``(x, y)``. Determinants and + orientation are exact; Python's unbounded ``int`` removes the ``long`` / + ``double`` overflow concern that forces FreeRouting to cap coordinates at + ``CRIT_INT`` and use a ``double`` determinant in ``side_of`` (this port keeps + ``side_of`` exact). +* :class:`RationalVector` — exact, projective ``(x, y, z)`` Python ``int`` + representing the affine vector ``(x/z, y/z)``; ``z`` is kept non-negative. + Mirrors the upstream ``BigInteger`` triple. ``BigIntDirection`` is not needed + because ``int`` is unbounded, so normalized directions are always + :class:`~freeroute.geometry.direction.IntDirection`. +""" + +from __future__ import annotations + +import math + +from .float_point import FloatPoint +from .side import Side, Signum + +__all__ = ["Vector", "IntVector", "RationalVector"] + + +def _add_rational(a: tuple[int, int, int], b: tuple[int, int, int]) -> tuple[int, int, int]: + """Add two projective coordinate triples (``BigIntAux.add_rational_coordinates``).""" + if a[2] == b[2]: + return (a[0] + b[0], a[1] + b[1], a[2]) + return ( + a[0] * b[2] + b[0] * a[2], + a[1] * b[2] + b[1] * a[2], + a[2] * b[2], + ) + + +class Vector: + """Abstract base for plane vectors (translations of points).""" + + def _xyz(self) -> tuple[int, int, int]: # pragma: no cover - overridden + raise NotImplementedError + + def is_zero(self) -> bool: # pragma: no cover - overridden + raise NotImplementedError + + def is_multiple_of_45_degree(self) -> bool: + return self.is_orthogonal() or self.is_diagonal() + + def is_orthogonal(self) -> bool: # pragma: no cover - overridden + raise NotImplementedError + + def is_diagonal(self) -> bool: # pragma: no cover - overridden + raise NotImplementedError + + def side_of(self, other: Vector) -> Side: + """Side of ``self`` relative to the directed line from zero to ``other``. + + Exact: uses the numerators of the projective coordinates (the positive + denominators ``z`` do not change the sign of the determinant). + """ + sx, sy, _ = self._xyz() + ox, oy, _ = other._xyz() + return Side.of(ox * sy - oy * sx) + + def projection(self, other: Vector) -> Signum: + """Sign of the scalar product of ``self`` and ``other`` (exact).""" + sx, sy, _ = self._xyz() + ox, oy, _ = other._xyz() + return Signum.of(sx * ox + sy * oy) + + def scalar_product(self, other: Vector) -> float: + v1 = self.to_float() + v2 = other.to_float() + return v1.x * v2.x + v1.y * v2.y + + def add(self, other: Vector) -> Vector: + if isinstance(self, IntVector) and isinstance(other, IntVector): + return IntVector(self.x + other.x, self.y + other.y) + x, y, z = _add_rational(self._xyz(), other._xyz()) + return RationalVector(x, y, z) + + def negate(self) -> Vector: # pragma: no cover - overridden + raise NotImplementedError + + def to_float(self) -> FloatPoint: # pragma: no cover - overridden + raise NotImplementedError + + def to_normalized_direction(self): + from .direction import IntDirection + + x, y, _ = self._xyz() + g = math.gcd(abs(x), abs(y)) + if g > 1: + x //= g + y //= g + return IntDirection(x, y) + + def length_approx(self) -> float: + return self.to_float().size() + + def cos_angle(self, other: Vector) -> float: + result = self.scalar_product(other) + return result / (self.to_float().size() * other.to_float().size()) + + def angle_approx(self, other: Vector | None = None) -> float: + if other is None: + return IntVector(1, 0).angle_approx(self) + result = math.acos(max(-1.0, min(1.0, self.cos_angle(other)))) + if self.side_of(other) == Side.ON_THE_LEFT: + result = -result + return result + + def add_to_point(self, point): # pragma: no cover - overridden + raise NotImplementedError + + +# Standard zero vector. +def _zero() -> IntVector: + return IntVector(0, 0) + + +class IntVector(Vector): + """A vector with exact integer coordinates.""" + + __slots__ = ("x", "y") + + def __init__(self, x: int, y: int) -> None: + self.x = int(x) + self.y = int(y) + + def _xyz(self) -> tuple[int, int, int]: + return (self.x, self.y, 1) + + def __eq__(self, other: object) -> bool: + if not isinstance(other, IntVector): + return NotImplemented + return self.x == other.x and self.y == other.y + + def __hash__(self) -> int: + return hash((self.x, self.y)) + + def __repr__(self) -> str: + return f"IntVector({self.x}, {self.y})" + + def is_zero(self) -> bool: + return self.x == 0 and self.y == 0 + + def negate(self) -> IntVector: + return IntVector(-self.x, -self.y) + + def is_orthogonal(self) -> bool: + return self.x == 0 or self.y == 0 + + def is_diagonal(self) -> bool: + return abs(self.x) == abs(self.y) + + def determinant(self, other: IntVector) -> int: + """Exact determinant of ``self`` and ``other`` (``x*oy - y*ox``).""" + return self.x * other.y - self.y * other.x + + def turn_90_degree(self, factor: int) -> IntVector: + n = factor % 4 + if n == 0: + return IntVector(self.x, self.y) + if n == 1: + return IntVector(-self.y, self.x) + if n == 2: + return IntVector(-self.x, -self.y) + return IntVector(self.y, -self.x) + + def mirror_at_y_axis(self) -> IntVector: + return IntVector(-self.x, self.y) + + def mirror_at_x_axis(self) -> IntVector: + return IntVector(self.x, -self.y) + + def to_float(self) -> FloatPoint: + return FloatPoint(self.x, self.y) + + def add_to_point(self, point): + from .point import IntPoint + + if isinstance(point, IntPoint): + return IntPoint(point.x + self.x, point.y + self.y) + return point.translate_by(self) + + def change_length_approx(self, length: float) -> Vector: + from .point import ZERO_POINT + + new_point = self.to_float().change_size(length) + return new_point.round().difference_by(ZERO_POINT) + + +class RationalVector(Vector): + """A vector with exact rational coordinates ``(x/z, y/z)``, ``z`` >= 0.""" + + __slots__ = ("x", "y", "z") + + def __init__(self, x: int, y: int, z: int) -> None: + if z >= 0: + self.x, self.y, self.z = int(x), int(y), int(z) + else: + self.x, self.y, self.z = -int(x), -int(y), -int(z) + + @staticmethod + def from_int_vector(v: IntVector) -> RationalVector: + return RationalVector(v.x, v.y, 1) + + def _xyz(self) -> tuple[int, int, int]: + return (self.x, self.y, self.z) + + def __eq__(self, other: object) -> bool: + if not isinstance(other, RationalVector): + return NotImplemented + # (x1,y1,z1) == (x2,y2,z2) iff cross products vanish. + return self.x * other.z - other.x * self.z == 0 and self.y * other.z - other.y * self.z == 0 + + def __hash__(self) -> int: + if self.z != 0: + g = math.gcd(math.gcd(abs(self.x), abs(self.y)), self.z) + g = g or 1 + return hash((self.x // g, self.y // g, self.z // g)) + return hash((self.x, self.y, 0)) + + def __repr__(self) -> str: + return f"RationalVector({self.x}, {self.y}, {self.z})" + + def is_zero(self) -> bool: + return self.x == 0 and self.y == 0 + + def negate(self) -> RationalVector: + return RationalVector(-self.x, -self.y, self.z) + + def is_orthogonal(self) -> bool: + return self.x == 0 or self.y == 0 + + def is_diagonal(self) -> bool: + return abs(self.x) == abs(self.y) + + def turn_90_degree(self, factor: int) -> RationalVector: + n = factor % 4 + if n == 0: + return RationalVector(self.x, self.y, self.z) + if n == 1: + return RationalVector(-self.y, self.x, self.z) + if n == 2: + return RationalVector(-self.x, -self.y, self.z) + return RationalVector(self.y, -self.x, self.z) + + def mirror_at_y_axis(self) -> RationalVector: + return RationalVector(-self.x, self.y, self.z) + + def mirror_at_x_axis(self) -> RationalVector: + return RationalVector(self.x, -self.y, self.z) + + def to_float(self) -> FloatPoint: + if self.z == 0: + big = 1.0e38 + return FloatPoint(big, big) + return FloatPoint(self.x / self.z, self.y / self.z) + + def add_to_point(self, point): + from .point import IntPoint + + if isinstance(point, IntPoint): + new_x = self.z * point.x + self.x + new_y = self.z * point.y + self.y + from .point import RationalPoint + + return RationalPoint(new_x, new_y, self.z) + return point.translate_by(self) diff --git a/tests/geometry/test_point_vector.py b/tests/geometry/test_point_vector.py new file mode 100644 index 0000000..2760da1 --- /dev/null +++ b/tests/geometry/test_point_vector.py @@ -0,0 +1,161 @@ +"""Exact-arithmetic tests for points and vectors. + +Oracles are hand-derived from the geometry defined in FreeRouting's +``geometry/planar/{IntPoint,IntVector,RationalPoint,RationalVector}.java``. +Emphasis is on exactness: collinearity, large coordinates that would overflow a +Java ``long``, and rational equality. +""" + +from __future__ import annotations + +from freeroute.geometry import ( + IntPoint, + IntVector, + RationalPoint, + RationalVector, + Side, + Signum, + point, + rational_point, +) + +# --- IntVector determinant / side_of ---------------------------------------- + + +def test_determinant_exact(): + # determinant of (1,0),(0,1) = 1 + assert IntVector(1, 0).determinant(IntVector(0, 1)) == 1 + assert IntVector(0, 1).determinant(IntVector(1, 0)) == -1 + assert IntVector(2, 3).determinant(IntVector(4, 6)) == 0 # parallel + + +def test_determinant_exact_beyond_java_long(): + # 10**12 * 10**12 = 10**24 overflows a Java long (~9.2e18); Python int is exact. + a = IntVector(10**12, 0) + b = IntVector(0, 10**12) + assert a.determinant(b) == 10**24 + + +def test_side_of_collinear_left_right(): + # line from origin toward (1,0): +y is left, -y is right, on axis collinear. + right_dir = IntVector(1, 0) + assert IntVector(0, 1).side_of(right_dir) == Side.ON_THE_LEFT + assert IntVector(0, -1).side_of(right_dir) == Side.ON_THE_RIGHT + assert IntVector(5, 0).side_of(right_dir) == Side.COLLINEAR + + +def test_side_of_exact_near_collinear(): + # (10**9, 10**9) vs (10**9, 10**9 + 1): a hair off collinear, exact int detects it. + base = IntVector(10**9, 10**9) + assert IntVector(10**9, 10**9 + 1).side_of(base) != Side.COLLINEAR + assert IntVector(2 * 10**9, 2 * 10**9).side_of(base) == Side.COLLINEAR + + +def test_vector_add_and_negate(): + assert IntVector(1, 2).add(IntVector(3, 4)) == IntVector(4, 6) + assert IntVector(1, 2).negate() == IntVector(-1, -2) + + +def test_turn_90_degree(): + v = IntVector(1, 0) + assert v.turn_90_degree(1) == IntVector(0, 1) + assert v.turn_90_degree(2) == IntVector(-1, 0) + assert v.turn_90_degree(3) == IntVector(0, -1) + assert v.turn_90_degree(4) == IntVector(1, 0) + assert v.turn_90_degree(-1) == IntVector(0, -1) + + +def test_projection_signum(): + assert IntVector(1, 0).projection(IntVector(1, 0)) == Signum.POSITIVE + assert IntVector(1, 0).projection(IntVector(-1, 0)) == Signum.NEGATIVE + assert IntVector(1, 0).projection(IntVector(0, 1)) == Signum.ZERO + + +def test_orthogonal_diagonal(): + assert IntVector(3, 0).is_orthogonal() + assert IntVector(0, -4).is_orthogonal() + assert IntVector(2, 2).is_diagonal() + assert IntVector(-5, 5).is_diagonal() + assert IntVector(1, 2).is_multiple_of_45_degree() is False + + +# --- IntPoint --------------------------------------------------------------- + + +def test_point_difference_and_translate(): + p = IntPoint(3, 4) + q = IntPoint(1, 1) + d = p.difference_by(q) + assert isinstance(d, IntVector) + assert d == IntVector(2, 3) + assert q.translate_by(d) == p + + +def test_point_side_of_directed_line(): + p0 = IntPoint(0, 0) + p1 = IntPoint(2, 2) + assert IntPoint(1, 1).side_of(p0, p1) == Side.COLLINEAR + assert IntPoint(0, 1).side_of(p0, p1) == Side.ON_THE_LEFT + assert IntPoint(1, 0).side_of(p0, p1) == Side.ON_THE_RIGHT + + +def test_point_compare(): + assert IntPoint(1, 5).compare_x(IntPoint(3, 2)) == -1 + assert IntPoint(3, 5).compare_x(IntPoint(3, 2)) == 0 + assert IntPoint(3, 5).compare_y(IntPoint(3, 2)) == 1 + assert IntPoint(3, 2).compare_x_y(IntPoint(3, 5)) == -1 # x tie, y smaller + + +def test_point_turn_90_around_pole(): + pole = IntPoint(1, 1) + # (2,1) rotated 90 about (1,1) -> (1,2) + assert IntPoint(2, 1).turn_90_degree(1, pole) == IntPoint(1, 2) + + +# --- RationalPoint / factory ------------------------------------------------ + + +def test_rational_point_equality_scale_invariant(): + assert RationalPoint(3, 2, 2) == RationalPoint(6, 4, 4) + assert RationalPoint(3, 2, 2) != RationalPoint(3, 2, 3) + + +def test_rational_point_negative_denominator_normalized_via_factory(): + p = rational_point(3, 2, -2) # -> (−3/−2)=... factory flips sign -> (-3,-2,2) + assert isinstance(p, RationalPoint) + assert p == RationalPoint(-3, -2, 2) + + +def test_rational_point_reduces_to_int_when_divisible(): + p = rational_point(6, 4, 2) + assert isinstance(p, IntPoint) + assert p == IntPoint(3, 2) + + +def test_rational_point_infinite(): + inf = rational_point(1, 1, 0) + assert isinstance(inf, RationalPoint) + assert inf.is_infinite() + + +def test_rational_point_contained_in_box_exact(): + from freeroute.geometry import IntBox + + box = IntBox(0, 0, 2, 2) + assert RationalPoint(3, 2, 2).is_contained_in(box) # (1.5, 1.0) inside + assert not RationalPoint(5, 2, 2).is_contained_in(box) # (2.5, 1.0) outside + + +def test_int_and_rational_difference_promotes(): + p = IntPoint(2, 2) + q = RationalPoint(3, 2, 2) # (1.5, 1.0) + d = p.difference_by(q) + assert isinstance(d, RationalVector) + # (2,2) - (1.5,1.0) = (0.5, 1.0) == (1,2,2) + assert d == RationalVector(1, 2, 2) + + +def test_factory_point_is_intpoint_even_for_huge_coords(): + p = point(10**15, -(10**15)) + assert isinstance(p, IntPoint) + assert (p.x, p.y) == (10**15, -(10**15))