Add exact planar points and vectors
Ports the point/vector foundation of geometry/planar: Side and Signum (three-valued signs), Limits, FloatPoint (approximate), and the exact IntPoint/IntVector plus projective RationalPoint/RationalVector. Arithmetic model per type is documented in each module. The key simplification over the Java source: Python's unbounded int makes the exact orientation determinants and rational (x,y,z) coordinates trivial, so side_of is kept exact (upstream uses a double for speed) and no BigInteger or CRIT_INT overflow promotion is needed. Rational points use the projective triple with z=0 denoting the point at infinity. Tests cover collinearity, determinants beyond Java long range, rational equality/reduction, and integer/rational promotion.
This commit is contained in:
parent
12b0a231f0
commit
33cb196bc8
159
src/freeroute/geometry/float_point.py
Normal file
159
src/freeroute/geometry/float_point.py
Normal file
@ -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)
|
||||
19
src/freeroute/geometry/limits.py
Normal file
19
src/freeroute/geometry/limits.py
Normal file
@ -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"]
|
||||
236
src/freeroute/geometry/point.py
Normal file
236
src/freeroute/geometry/point.py
Normal file
@ -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
|
||||
56
src/freeroute/geometry/side.py
Normal file
56
src/freeroute/geometry/side.py
Normal file
@ -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)
|
||||
274
src/freeroute/geometry/vector.py
Normal file
274
src/freeroute/geometry/vector.py
Normal file
@ -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)
|
||||
161
tests/geometry/test_point_vector.py
Normal file
161
tests/geometry/test_point_vector.py
Normal file
@ -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))
|
||||
Loading…
x
Reference in New Issue
Block a user