Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
17 changes: 11 additions & 6 deletions src/shapepy/analytic/bezier.py
Original file line number Diff line number Diff line change
Expand Up @@ -7,11 +7,11 @@
from functools import lru_cache
from typing import Iterable, Tuple, Union

from ..rbool import SubSetR1, WholeR1, from_any
from ..rbool import SubSetR1
from ..scalar.quadrature import inner
from ..scalar.reals import Math, Rational, Real
from ..tools import Is, To
from .polynomial import Polynomial
from .polynomial import Polynomial, scale_coefs, shift_coefs


@lru_cache(maxsize=None)
Expand Down Expand Up @@ -76,8 +76,13 @@ class Bezier(Polynomial):
"""

def __init__(
self, coefs: Iterable[Real], domain: Union[None, SubSetR1] = None
self,
coefs: Iterable[Real],
reparam: Tuple[Real, Real] = (0, 1),
*,
domain: Union[None, SubSetR1] = None,
):
domain = WholeR1() if domain is None else from_any(domain)
poly_coefs = tuple(bezier2polynomial(coefs))
super().__init__(poly_coefs, domain)
coefs = tuple(bezier2polynomial(coefs))
knota, knotb = reparam
coefs = shift_coefs(scale_coefs(coefs, knotb - knota), knota)
super().__init__(coefs, domain=domain)
69 changes: 43 additions & 26 deletions src/shapepy/analytic/polynomial.py
Original file line number Diff line number Diff line change
Expand Up @@ -5,14 +5,36 @@
from __future__ import annotations

from numbers import Real
from typing import Iterable, List, Union
from typing import Iterable, Iterator, List, Union

from ..rbool import IntervalR1, SubSetR1, WholeR1, from_any, move, scale
from ..rbool.tools import is_continuous
from ..scalar.reals import Math
from ..tools import Is, To
from .base import IAnalytic


def scale_coefs(coefs: Iterable[Real], amount: Real) -> Iterator[Real]:
"""Computes the polynomial p(t/A) from the coefficients of p(t)"""
if amount != 1:
inv = 1 / amount
coefs = (coef * inv**i for i, coef in enumerate(coefs))
yield from coefs


def shift_coefs(coefs: Iterable[Real], amount: Real) -> Iterator[Real]:
"""Computes the polynomial p(t-a) from the coefficients of p(t)"""
if amount != 0:
coefs = list(coefs)
for i, coef in enumerate(tuple(coefs)):
for j in range(i):
value = Math.binom(i, j) * (amount ** (i - j))
if (i + j) % 2:
value *= -1
coefs[j] += coef * value
yield from coefs


class Polynomial(IAnalytic):
"""
Defines a polynomial with coefficients
Expand All @@ -36,15 +58,20 @@ class Polynomial(IAnalytic):
5
"""

def __init__(self, coefs: Iterable[Real], domain: SubSetR1 = WholeR1()):
def __init__(
self, coefs: Iterable[Real], *, domain: Union[None, SubSetR1] = None
):
domain = WholeR1() if domain is None else from_any(domain)
if not is_continuous(domain):
raise ValueError(f"Domain {domain} is not continuous")
if not Is.iterable(coefs):
raise TypeError("Expected an iterable of coefficients")
coefs = tuple(coefs)
if len(coefs) == 0:
raise ValueError("Cannot receive an empty tuple")
degree = max((i for i, v in enumerate(coefs) if v * v > 0), default=0)
self.__coefs = coefs[: degree + 1]
self.__domain = from_any(domain)
self.__domain = domain

@property
def domain(self) -> SubSetR1:
Expand Down Expand Up @@ -81,26 +108,28 @@ def __add__(self, other: Union[Real, Polynomial]) -> Polynomial:
if not Is.instance(other, IAnalytic):
coefs = list(self)
coefs[0] += other
return Polynomial(coefs, self.domain)
return Polynomial(coefs, domain=self.domain)
if not Is.instance(other, Polynomial):
return NotImplemented
raise NotImplementedError
coefs = [0] * (1 + max(self.degree, other.degree))
for i, coef in enumerate(self):
coefs[i] += coef
for i, coef in enumerate(other):
coefs[i] += coef
return Polynomial(coefs, self.domain)
return Polynomial(coefs, domain=self.domain)

def __mul__(self, other: Union[Real, Polynomial]) -> Polynomial:
if not Is.instance(other, IAnalytic):
return Polynomial((other * coef for coef in self), self.domain)
return Polynomial(
(other * coef for coef in self), domain=self.domain
)
if not Is.instance(other, Polynomial):
return NotImplemented
raise NotImplementedError
coefs = [0 * self[0]] * (self.degree + other.degree + 1)
for i, coefi in enumerate(self):
for j, coefj in enumerate(other):
coefs[i + j] += coefi * coefj
return Polynomial(coefs, self.domain & other.domain)
return Polynomial(coefs, domain=self.domain & other.domain)

def eval(self, node: Real, derivate: int = 0) -> Real:
if node not in self.domain:
Expand Down Expand Up @@ -130,12 +159,12 @@ def eval(self, node: Real, derivate: int = 0) -> Real:

def derivate(self, times=1):
if self.degree < times:
return Polynomial([0 * self[0]], self.domain)
return Polynomial([0 * self[0]], domain=self.domain)
coefs = (
Math.factorial(n + times) // Math.factorial(n) * coef
for n, coef in enumerate(self[times:])
)
return Polynomial(coefs, self.domain)
return Polynomial(coefs, domain=self.domain)

def integrate(self, domain):
domain = from_any(domain)
Expand All @@ -161,21 +190,9 @@ def compose(self, function: IAnalytic) -> Polynomial:
if function.degree != 1:
raise ValueError("Only polynomials of degree = 1 are allowed")
shift_amount, scale_amount = tuple(function)
coefs = list(self)
domain = self.domain
if scale_amount != 1:
inv = 1 / scale_amount
coefs = list(coef * inv**i for i, coef in enumerate(self))
domain = scale(domain, scale_amount)
if shift_amount != 0:
for i, coef in enumerate(tuple(coefs)):
for j in range(i):
value = Math.binom(i, j) * (shift_amount ** (i - j))
if (i + j) % 2:
value *= -1
coefs[j] += coef * value
domain = move(domain, shift_amount)
return Polynomial(coefs, domain)
coefs = shift_coefs(scale_coefs(self, scale_amount), shift_amount)
domain = move(scale(self.domain, scale_amount), shift_amount)
return Polynomial(coefs, domain=domain)

def __repr__(self):
return str(self.domain) + ": " + self.__str__()
Expand Down
5 changes: 5 additions & 0 deletions src/shapepy/analytic/tools.py
Original file line number Diff line number Diff line change
Expand Up @@ -62,6 +62,11 @@ def find_minimum(
raise NotExpectedError(f"Invalid analytic: {type(analytic)}")


def is_constant(analytic: IAnalytic) -> bool:
"""Tells if the given analytic function is constant"""
return Is.instance(analytic, Polynomial) and analytic.degree == 0


class PolynomialFunctions:
"""Static class that stores static functions used for the generics
functions above. This class specifics for Polynomial"""
Expand Down
21 changes: 11 additions & 10 deletions src/shapepy/bool2d/boolean.py
Original file line number Diff line number Diff line change
Expand Up @@ -6,7 +6,6 @@
from __future__ import annotations

from copy import copy
from fractions import Fraction
from typing import Dict, Iterable, Tuple, Union

from shapepy.geometry.jordancurve import JordanCurve
Expand Down Expand Up @@ -370,12 +369,14 @@ def split_on_intersection(
jordansj = all_group_jordans[j]
for jordana in jordansi:
for jordanb in jordansj:
intersection |= jordana.piecewise & jordanb.piecewise
intersection |= (
jordana.parametrize() & jordanb.parametrize()
)
intersection.evaluate()
for jordans in all_group_jordans:
for jordan in jordans:
split_knots = intersection.all_knots[id(jordan.piecewise)]
jordan.piecewise.split(split_knots)
split_knots = intersection.all_knots[id(jordan.parametrize())]
jordan.parametrize().split(split_knots)

@staticmethod
def pursue_path(
Expand All @@ -396,14 +397,14 @@ def pursue_path(
We suppose there's no triple intersection
"""
matrix = []
all_segments = [tuple(jordan.piecewise) for jordan in jordans]
all_segments = [tuple(jordan.parametrize()) for jordan in jordans]
while True:
index_segment %= len(all_segments[index_jordan])
segment = all_segments[index_jordan][index_segment]
if (index_jordan, index_segment) in matrix:
break
matrix.append((index_jordan, index_segment))
last_point = segment(1)
last_point = segment(segment.knots[-1])
possibles = []
for i, jordan in enumerate(jordans):
if i == index_jordan:
Expand All @@ -415,7 +416,7 @@ def pursue_path(
continue
index_jordan = possibles[0]
for j, segj in enumerate(all_segments[index_jordan]):
if segj(0) == last_point:
if segj(segj.knots[0]) == last_point:
index_segment = j
break
return CyclicContainer(matrix)
Expand All @@ -434,7 +435,7 @@ def indexs_to_jordan(
"""
beziers = []
for index_jordan, index_segment in matrix_indexs:
new_bezier = jordans[index_jordan].piecewise[index_segment]
new_bezier = jordans[index_jordan].parametrize()[index_segment]
new_bezier = copy(new_bezier)
beziers.append(USegment(new_bezier))
new_jordan = JordanCurve(beziers)
Expand All @@ -448,7 +449,7 @@ def follow_path(
Returns a list of jordan curves which is the result
of the intersection between 'jordansa' and 'jordansb'
"""
assert all(map(Is.jordan, jordans))
assert all(Is.instance(j, JordanCurve) for j in jordans)
bez_indexs = []
for ind_jord, ind_seg in start_indexs:
indices_matrix = FollowPath.pursue_path(ind_jord, ind_seg, jordans)
Expand Down Expand Up @@ -480,7 +481,7 @@ def midpoints_one_shape(
"""
for i, jordan in enumerate(shapea.jordans):
for j, segment in enumerate(jordan.parametrize()):
mid_point = segment(Fraction(1, 2))
mid_point = segment((segment.knots[0] + segment.knots[-1]) / 2)
density = shapeb.density(mid_point)
mid_point_in = (float(density) > 0 and closed) or density == 1
if not inside ^ mid_point_in:
Expand Down
12 changes: 6 additions & 6 deletions src/shapepy/bool2d/density.py
Original file line number Diff line number Diff line change
Expand Up @@ -32,13 +32,13 @@ def half_density_jordan(
deltax = segment.xfunc - point.xcoord
deltay = segment.yfunc - point.ycoord
radius_square = deltax * deltax + deltay * deltay
minimal = find_minimum(radius_square, [0, 1])
minimal = find_minimum(radius_square, segment.domain)
if minimal < 1e-6:
place = where_minimum(radius_square, [0, 1])
place = where_minimum(radius_square, segment.domain)
if not Is.instance(place, SingleR1):
raise NotExpectedError(f"Not single value: {place}")
parameter = To.finite(place.internal)
angle = segment(parameter, 1).angle
angle = segment.eval(parameter, 1).angle
return line(angle)
raise NotExpectedError("Not found minimum < 1e-6")

Expand All @@ -61,10 +61,10 @@ def lebesgue_density_jordan(

segments = tuple(jordan.parametrize())
for i, segmenti in enumerate(segments):
if point == segmenti(0):
if point == segmenti(segmenti.knots[0]):
segmentj = segments[(i - 1) % len(segments)]
anglei = segmenti(0, 1).angle
anglej = segmentj(1, 1).angle
anglei = segmenti.eval(segmenti.knots[0], 1).angle
anglej = segmentj.eval(segmentj.knots[-1], 1).angle
return sector(anglei, ~anglej)

turns = IntegrateJordan.turns(jordan, point)
Expand Down
22 changes: 18 additions & 4 deletions src/shapepy/geometry/base.py
Original file line number Diff line number Diff line change
Expand Up @@ -5,9 +5,11 @@
from __future__ import annotations

from abc import ABC, abstractmethod
from typing import Iterable, Tuple
from typing import Iterable, Tuple, Union

from ..rbool import IntervalR1, WholeR1
from ..scalar.reals import Real
from ..tools import vectorize
from .box import Box
from .point import Point2D

Expand Down Expand Up @@ -72,19 +74,31 @@ class IParametrizedCurve(IGeometricCurve):
Class interface for parametrized curves
"""

@property
@abstractmethod
def domain(self) -> Union[IntervalR1, WholeR1]:
"""
The domain where the curve is defined.
"""
raise NotImplementedError

@property
@abstractmethod
def knots(self) -> Tuple[Real, ...]:
"""
The length of the curve
If the curve is not bounded, returns infinity
The subdivisions on the domain
"""
raise NotImplementedError

@abstractmethod
def __call__(self, node: Real, derivate: int = 0) -> Point2D:
def eval(self, node: Real, derivate: int = 0) -> Point2D:
"""Evaluates the curve at given node"""
raise NotImplementedError

@vectorize(1, 0)
def __call__(self, node: Real) -> Point2D:
return self.eval(node, 0)

def __and__(self, other: IParametrizedCurve):
return Future.intersect(self, other)

Expand Down
Loading