From 533799de9ae6a56df5b17836fe7b395a3d7a8466 Mon Sep 17 00:00:00 2001 From: Eivind Fonn Date: Fri, 2 Jan 2026 15:27:23 +0100 Subject: [PATCH 01/22] Typing: volume_factory.py --- src/splipy/volume_factory.py | 11 +++++------ 1 file changed, 5 insertions(+), 6 deletions(-) diff --git a/src/splipy/volume_factory.py b/src/splipy/volume_factory.py index 3bb06a4..75624cd 100644 --- a/src/splipy/volume_factory.py +++ b/src/splipy/volume_factory.py @@ -459,9 +459,6 @@ def loft(*srfs: Surface | Sequence[Surface]) -> Volume: # ensure all centers have the same dimension (typically some are 2D points and others are 3D points) max_dim = max(s.dimension for s in surfaces) x = [np.pad(xi, (0, max_dim - len(xi)), mode="constant") for xi in x] - dist = [0] - for x1, x0 in zip(x[1:], x[:-1]): - dist.append(dist[-1] + np.linalg.norm(x1 - x0)) # clone input, so we don't change those references # make sure everything has the same dimension since we need to compute length @@ -469,6 +466,8 @@ def loft(*srfs: Surface | Sequence[Surface]) -> Volume: if len(surfaces) == 2: return edge_surfaces(surfaces) if len(surfaces) == 3: + dist = curve_length_parametrization(x) + # can't do cubic spline interpolation, so we'll do quadratic basis3 = BSplineBasis(3) else: @@ -502,13 +501,13 @@ def loft(*srfs: Surface | Sequence[Surface]) -> Volume: Nw_inv = np.linalg.inv(Nw) # compute interpolation points in physical space - x = np.zeros((m1, m2, n, dim)) + cp = np.zeros((m1, m2, n, dim)) for i in range(n): tmp = np.tensordot(Nv, surfaces[i].controlpoints, axes=(1, 1)) - x[:, :, i, :] = np.tensordot(Nu, tmp, axes=(1, 1)) + cp[:, :, i, :] = np.tensordot(Nu, tmp, axes=(1, 1)) # solve interpolation problem - cp = np.tensordot(Nw_inv, x, axes=(1, 2)) + cp = np.tensordot(Nw_inv, cp, axes=(1, 2)) cp = np.tensordot(Nv_inv, cp, axes=(1, 2)) cp = np.tensordot(Nu_inv, cp, axes=(1, 2)) From 6b0ca9529376bb6cade86be053cfeae75554f51a Mon Sep 17 00:00:00 2001 From: Eivind Fonn Date: Fri, 2 Jan 2026 15:31:53 +0100 Subject: [PATCH 02/22] Typing: surface_factory.py --- src/splipy/surface_factory.py | 14 ++++++++------ 1 file changed, 8 insertions(+), 6 deletions(-) diff --git a/src/splipy/surface_factory.py b/src/splipy/surface_factory.py index 97f1e35..332960b 100644 --- a/src/splipy/surface_factory.py +++ b/src/splipy/surface_factory.py @@ -679,6 +679,9 @@ def thicken(curve: Curve, amount: Scalar | Callable[..., float]) -> Surface: return edge_curves(right, left) # dimension=3, we will create a surrounding tube + # callable amount is not supported + if callable(amount): + raise TypeError("Callable amount not supported for three-dimensional curves") return sweep(curve, curve_factory.circle(r=amount)) @@ -764,9 +767,6 @@ def loft(*in_curves: Curve | Sequence[Curve]) -> Surface: # ensure all centers have the same dimension (typically some are 2D points and others are 3D points) max_dim = max(c.dimension for c in curves) x = [np.pad(xi, (0, max_dim - len(xi)), mode="constant") for xi in x] - dist = [0] - for x1, x0 in zip(x[1:], x[:-1]): - dist.append(dist[-1] + np.linalg.norm(x1 - x0)) # clone input, so we don't change those references # make sure everything has the same dimension since we need to compute length @@ -774,6 +774,8 @@ def loft(*in_curves: Curve | Sequence[Curve]) -> Surface: if len(curves) == 2: return edge_curves(curves) if len(curves) == 3: + dist = curve_length_parametrization(x) + # can't do cubic spline interpolation, so we'll do quadratic basis2 = BSplineBasis(3) else: @@ -801,12 +803,12 @@ def loft(*in_curves: Curve | Sequence[Curve]) -> Surface: Nv_inv = np.linalg.inv(Nv) # compute interpolation points in physical space - x = np.zeros((m, n, curves[0][0].size)) + cp: FloatArray = np.zeros((m, n, curves[0][0].size)) for i in range(n): - x[:, i, :] = Nu @ curves[i].controlpoints + cp[:, i, :] = Nu @ curves[i].controlpoints # solve interpolation problem - cp = np.tensordot(Nv_inv, x, axes=(1, 1)) + cp = np.tensordot(Nv_inv, cp, axes=(1, 1)) cp = np.tensordot(Nu_inv, cp, axes=(1, 1)) # re-order controlpoints so they match up with Surface constructor From 7c1584697d15fde99ce4c51aaab319709b6ea1f9 Mon Sep 17 00:00:00 2001 From: Eivind Fonn Date: Fri, 2 Jan 2026 15:32:28 +0100 Subject: [PATCH 03/22] Typing: svg.py --- src/splipy/io/svg.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/splipy/io/svg.py b/src/splipy/io/svg.py index 95620c2..cb02f39 100644 --- a/src/splipy/io/svg.py +++ b/src/splipy/io/svg.py @@ -74,7 +74,7 @@ class SVG(MasterIO): margin: float all_objects: list[SplineObject] - all_kwargs: list[tuple[Any]] + all_kwargs: list[dict[str, Any]] center: tuple[float, float] offset: tuple[float, float] From 400bc28b331a72276083f914125b2c5f1d511a7c Mon Sep 17 00:00:00 2001 From: Eivind Fonn Date: Fri, 2 Jan 2026 23:26:31 +0100 Subject: [PATCH 04/22] More convenient type aliases --- src/splipy/typing.py | 30 ++++++++++++++++++++++++++++-- 1 file changed, 28 insertions(+), 2 deletions(-) diff --git a/src/splipy/typing.py b/src/splipy/typing.py index 88b6b5c..28684ac 100644 --- a/src/splipy/typing.py +++ b/src/splipy/typing.py @@ -1,16 +1,42 @@ from __future__ import annotations -from typing import Literal, SupportsFloat, TypedDict, TypeVar +from collections.abc import Sequence +from typing import Literal, TypedDict, TypeVar import numpy as np import numpy.typing as npt B = TypeVar("B", bound=npt.NBitBase) +# Anything that can be converted to an array. Use this type for parameters that +# are (multidimensional) lists of points, such as controlpoints. Generally, only +# use this for public-facing functions. Use np.asarray to convert to an actual +# array, and then use FloatArray internally in Splipy. Avoid using np.array, +# since that will make a copy of objects that are already arrays. If you need a +# copy to mutate, be explicit and use np.asarray().copy(). type ArrayLike = npt.ArrayLike -type Scalar = SupportsFloat + +# Alias to signal that we expect control points - makes no difference to type +# checking but it's helpful for humans. +type ControlPoints = ArrayLike + +# Alias to signal that we expect points (which are not precisely the same thing +# as control points) - makes no difference to type checking but it's helpful for +# humans. +type Points = ArrayLike + +# Anything that acts and behaves as a float. Note that np.floating and +# np.integer encompasses floats of many different sizes. +type Scalar = float | np.floating | int | np.integer + type FloatArray = npt.NDArray[np.floating] type IntArray = npt.NDArray[np.integer] + +# Again, these are identical, but we can signal intent to a human reader based +# on which we use. +type Knots = Sequence[Scalar] | FloatArray | IntArray +type Point = Sequence[Scalar] | FloatArray | IntArray + type Direction = Literal["u", "v", "w", "U", "V", "W"] | int type SectionElement = Literal[-1, 0] | None type Section = tuple[SectionElement, ...] From 47d18f1ba4abec79b42267bddb9b7b9e82114176 Mon Sep 17 00:00:00 2001 From: Eivind Fonn Date: Fri, 2 Jan 2026 23:28:32 +0100 Subject: [PATCH 05/22] Typing: zero errors in the factory files --- src/splipy/basis.py | 14 +- src/splipy/curve.py | 4 +- src/splipy/curve_factory.py | 294 +++++++++++++++++++--------------- src/splipy/io/threedm.py | 3 +- src/splipy/splineobject.py | 34 ++-- src/splipy/surface_factory.py | 63 ++++---- src/splipy/utils/__init__.py | 29 +++- src/splipy/utils/curve.py | 8 +- src/splipy/volume_factory.py | 34 ++-- 9 files changed, 286 insertions(+), 197 deletions(-) diff --git a/src/splipy/basis.py b/src/splipy/basis.py index 9c3e723..69d1442 100644 --- a/src/splipy/basis.py +++ b/src/splipy/basis.py @@ -15,7 +15,9 @@ from .utils import ensure_listlike_old if TYPE_CHECKING: - from .typing import ArrayLike, FloatArray, Scalar + from splipy.typing import Knots + + from .typing import FloatArray, Scalar __all__ = ["BSplineBasis"] @@ -41,7 +43,7 @@ class BSplineBasis: def __init__( self, order: int = 2, - knots: ArrayLike | None = None, + knots: Knots | None = None, periodic: int = -1, ) -> None: """Construct a B-Spline basis with a given order and knot vector. @@ -151,7 +153,7 @@ def greville(self, index: int | None = None) -> float | FloatArray: @overload def evaluate( self, - t: ArrayLike | Scalar, + t: Knots | Scalar, d: int = 0, from_right: bool = ..., ) -> npt.NDArray[np.double]: ... @@ -159,7 +161,7 @@ def evaluate( @overload def evaluate( self, - t: ArrayLike | Scalar, + t: Knots | Scalar, d: int = 0, from_right: bool = ..., sparse: Literal[False] = ..., @@ -168,7 +170,7 @@ def evaluate( @overload def evaluate( self, - t: ArrayLike | Scalar, + t: Knots | Scalar, d: int = 0, from_right: bool = ..., sparse: Literal[True] = ..., @@ -176,7 +178,7 @@ def evaluate( def evaluate( self, - t: ArrayLike | Scalar, + t: Knots | Scalar, d: int = 0, from_right: bool = True, sparse: bool = False, diff --git a/src/splipy/curve.py b/src/splipy/curve.py index 48502b7..35dc79e 100644 --- a/src/splipy/curve.py +++ b/src/splipy/curve.py @@ -12,7 +12,7 @@ from .utils import ensure_listlike, is_singleton if TYPE_CHECKING: - from collections.abc import Sequence + from collections.abc import Callable, Sequence from .typing import ArrayLike, Direction, FloatArray, Scalar @@ -494,7 +494,7 @@ def closest_point(self, pt: ArrayLike, t0: Scalar = None) -> tuple[FloatArray, f break return self(t), t - def error(self, target: Curve) -> tuple[FloatArray, float]: + def error(self, target: Curve | Callable[[FloatArray], FloatArray]) -> tuple[FloatArray, float]: """Computes the L2 (squared and per knot span) and max error between this curve and a target curve diff --git a/src/splipy/curve_factory.py b/src/splipy/curve_factory.py index e476a5c..e34a994 100644 --- a/src/splipy/curve_factory.py +++ b/src/splipy/curve_factory.py @@ -2,22 +2,34 @@ from __future__ import annotations -import copy import inspect -from math import ceil, cos, pi, sin, sqrt -from typing import TYPE_CHECKING +from enum import IntEnum +from math import ceil, pi, sqrt +from typing import TYPE_CHECKING, Any, Literal, overload import numpy as np import scipy.sparse as sp import scipy.sparse.linalg as splinalg from numpy.linalg import norm +from splipy.utils.curve import curve_length_parametrization + from . import state from .basis import BSplineBasis from .curve import Curve -from .utils import flip_and_move_plane_geometry, rotate_local_x_axis +from .utils import ( + flip_and_move_plane_geometry, + knot_vector, + normalize_points, + rotate_local_x_axis, + with_repeated_knots, +) if TYPE_CHECKING: + from collections.abc import Callable + + from splipy.typing import FloatArray, Knots, Point, Points + from .typing import Scalar __all__ = [ @@ -38,7 +50,11 @@ ] -class Boundary: +type P2C0 = Literal["p2c0", "c0p2", "P2c0", "c0P2", "p2C0", "C0p2", "P2C0", "C0P2"] +type P4C1 = Literal["p4c1", "c1p4", "P4c1", "c1P4", "p4C1", "C1p4", "P4C1", "C1P4"] + + +class Boundary(IntEnum): """Enumeration representing different boundary conditions used in :func:`interpolate`.""" @@ -61,7 +77,7 @@ class Boundary: """Use `TANGENT` for the start and `NATURAL` for the end.""" -def line(a, b, relative=False) -> Curve: +def line(a: Point, b: Point, relative: bool = False) -> Curve: """Create a line between two points. :param array-like a: Start point @@ -70,12 +86,22 @@ def line(a, b, relative=False) -> Curve: :return: Linear curve from *a* to *b* :rtype: Curve """ + a = np.asarray(a) + b = np.asarray(b) if relative: - b = tuple(ai + bi for ai, bi in zip(a, b)) + return Curve(controlpoints=[a, a + b]) return Curve(controlpoints=[a, b]) -def polygon(*points, **keywords): +@overload +def polygon(*points: Point, t: Knots | None = None, relative: bool = False) -> Curve: ... + + +@overload +def polygon(points: Points, /, *, t: Knots | None = None, relative: bool = False) -> Curve: ... + + +def polygon(*points_in: Point | Points, t: Knots | None = None, relative: bool = False) -> Curve: """Create a linear interpolation between input points. :param [array-like] points: The points to interpolate @@ -85,33 +111,20 @@ def polygon(*points, **keywords): :return: Linear curve through the input points :rtype: Curve """ - if len(points) == 1: - points = points[0] - - knot = keywords.get("t", []) - if len(knot) == 0: # establish knot vector based on eucledian length between points - knot = [0, 0] - prevPt = points[0] - for pt in points[1:]: - dist = 0 - for x0, x1 in zip(prevPt, pt): # loop over (x,y) and maybe z-coordinate - dist += (x1 - x0) ** 2 - knot.append(knot[-1] + sqrt(dist)) - prevPt = pt - knot.append(knot[-1]) + points = normalize_points(*points_in) + + if t is None: # establish knot vector based on Euclidean length between points + knot = curve_length_parametrization(points, reps=2) else: # use knot vector given as input argument - knot = [knot[0]] + list(knot) + [knot[-1]] + knot = with_repeated_knots(np.asarray(t), reps=2) - relative = keywords.get("relative", False) if relative: - points = list(points) - for i in range(1, len(points)): - points[i] = [x0 + x1 for (x0, x1) in zip(points[i - 1], points[i])] + points = np.cumsum(points, axis=0) return Curve(BSplineBasis(2, knot), points) -def n_gon(n=5, r=1, center=(0, 0, 0), normal=(0, 0, 1)): +def n_gon(n: int = 5, r: Scalar = 1, center: Point = (0, 0, 0), normal: Point = (0, 0, 1)) -> Curve: """Create a regular polygon of *n* equal sides centered at the origin. :param int n: Number of sides and vertices @@ -128,20 +141,21 @@ def n_gon(n=5, r=1, center=(0, 0, 0), normal=(0, 0, 1)): if n < 3: raise ValueError("regular polygons need at least 3 sides") - cp = [] dt = 2 * pi / n - knot = [-1] - for i in range(n): - cp.append([r * cos(i * dt), r * sin(i * dt)]) - knot.append(i) - knot += [n, n + 1] + cp = r * np.column_stack([np.cos(np.arange(n) * dt), np.sin(np.arange(n) * dt)]) + knot = np.arange(-1, n + 2) basis = BSplineBasis(2, knot, 0) - result = Curve(basis, cp) return flip_and_move_plane_geometry(result, center, normal) -def circle(r: Scalar = 1, center=(0, 0, 0), normal=(0, 0, 1), type="p2C0", xaxis=(1, 0, 0)) -> Curve: +def circle( + r: Scalar = 1, + center: Point = (0, 0, 0), + normal: Point = (0, 0, 1), + type: P2C0 | P4C1 = "p2c0", + xaxis: Point = (1, 0, 0), +) -> Curve: """Create a circle. :param float r: Radius @@ -155,40 +169,48 @@ def circle(r: Scalar = 1, center=(0, 0, 0), normal=(0, 0, 1), type="p2C0", xaxis """ if r <= 0: raise ValueError("radius needs to be positive") + tp = type.lower() - if type == "p2C0" or type == "C0p2": + if tp == "p2c0" or tp == "c0p2": w = 1.0 / sqrt(2) - controlpoints = [ - [1, 0, 1], - [w, w, w], - [0, 1, 1], - [-w, w, w], - [-1, 0, 1], - [-w, -w, w], - [0, -1, 1], - [w, -w, w], - ] + controlpoints = np.array( + [ + [1, 0, 1], + [w, w, w], + [0, 1, 1], + [-w, w, w], + [-1, 0, 1], + [-w, -w, w], + [0, -1, 1], + [w, -w, w], + ], + dtype=float, + ) knot = np.array([-1, 0, 0, 1, 1, 2, 2, 3, 3, 4, 4, 5]) / 4.0 * 2 * pi result = Curve(BSplineBasis(3, knot, 0), controlpoints, True) - elif type.lower() == "p4c1" or type.lower() == "c1p4": + + elif tp == "p4c1" or tp == "c1p4": w = 2 * sqrt(2) / 3 a = 1.0 / 2 / sqrt(2) b = 1.0 / 6 * (4 * sqrt(2) - 1) - controlpoints = [ - [1, -a, 1], - [1, a, 1], - [b, b, w], - [a, 1, 1], - [-a, 1, 1], - [-b, b, w], - [-1, a, 1], - [-1, -a, 1], - [-b, -b, w], - [-a, -1, 1], - [a, -1, 1], - [b, -b, w], - ] + controlpoints = np.array( + [ + [1, -a, 1], + [1, a, 1], + [b, b, w], + [a, 1, 1], + [-a, 1, 1], + [-b, b, w], + [-1, a, 1], + [-1, -a, 1], + [-b, -b, w], + [-a, -1, 1], + [a, -1, 1], + [b, -b, w], + ], + dtype=float, + ) knot = np.array([-1, -1, 0, 0, 0, 1, 1, 1, 2, 2, 2, 3, 3, 3, 4, 4, 4, 5, 5]) / 4.0 * 2 * pi result = Curve(BSplineBasis(5, knot, 1), controlpoints, True) else: @@ -199,7 +221,14 @@ def circle(r: Scalar = 1, center=(0, 0, 0), normal=(0, 0, 1), type="p2C0", xaxis return flip_and_move_plane_geometry(result, center, normal) -def ellipse(r1=1, r2=1, center=(0, 0, 0), normal=(0, 0, 1), type="p2C0", xaxis=(1, 0, 0)) -> Curve: +def ellipse( + r1: Scalar = 1, + r2: Scalar = 1, + center: Point = (0, 0, 0), + normal: Point = (0, 0, 1), + type: P2C0 | P4C1 = "p2C0", + xaxis: Point = (1, 0, 0), +) -> Curve: """Create an ellipse :param float r1: Radius along xaxis @@ -218,7 +247,7 @@ def ellipse(r1=1, r2=1, center=(0, 0, 0), normal=(0, 0, 1), type="p2C0", xaxis=( return flip_and_move_plane_geometry(result, center, normal) -def circle_segment_from_three_points(x0, x1, x2): +def circle_segment_from_three_points(x0: Point, x1: Point, x2: Point) -> Curve: """circle_segment_from_three_points(x0, x1, x2) Create a circle segment going from the point x0 to x2 through x1 @@ -270,7 +299,13 @@ def circle_segment_from_three_points(x0, x1, x2): return result -def circle_segment(theta, r=1, center=(0, 0, 0), normal=(0, 0, 1), xaxis=(1, 0, 0)): +def circle_segment( + theta: Scalar, + r: Scalar = 1, + center: Point = (0, 0, 0), + normal: Point = (0, 0, 1), + xaxis: Point = (1, 0, 0), +) -> Curve: """Create a circle segment starting parallel to the rotated x-axis. :param float theta: Angle in radians @@ -284,7 +319,7 @@ def circle_segment(theta, r=1, center=(0, 0, 0), normal=(0, 0, 1), xaxis=(1, 0, :raises ValueError: If theta is not in the range *[-2pi, 2pi]* """ # error test input - if abs(theta) > 2 * pi: + if np.abs(theta) > 2 * pi: raise ValueError("theta needs to be in range [-2pi,2pi]") if r <= 0: raise ValueError("radius needs to be positive") @@ -292,28 +327,21 @@ def circle_segment(theta, r=1, center=(0, 0, 0), normal=(0, 0, 1), xaxis=(1, 0, return circle(r, center, normal) # build knot vector - knot_spans = int(ceil(abs(theta) / (2 * pi / 3))) - knot = [0] - for i in range(knot_spans + 1): - knot += [i] * 2 - knot += [knot_spans] # knot vector [0,0,0,1,1,2,2,..,n,n,n] - knot = np.array(knot) / float(knot[-1]) * theta # set parametric space to [0,theta] + knot_spans = int(ceil(np.abs(theta) / (2 * pi / 3))) + knot = knot_vector(0, theta, num_intervals=knot_spans, interior_reps=2, endpoint_reps=3) n = (knot_spans - 1) * 2 + 3 # number of control points needed - cp = [] - t = 0 # current angle dt = float(theta) / knot_spans / 2 # angle step - - # build control points - for i in range(n): - w = 1 - (i % 2) * (1 - cos(dt)) # weights = 1 and cos(dt) every other i - x = r * cos(t) - y = r * sin(t) - cp += [[x, y, w]] - t += dt + cp = np.column_stack( + [ + r * np.cos(np.arange(n) * dt), + r * np.sin(np.arange(n) * dt), + 1 - (np.arange(n) % 2) * (1 - np.cos(dt)), # weights = 1 and cos(dt) every other i + ] + ) if theta < 0: - cp.reverse() + cp = cp[::-1, ...] result = Curve(BSplineBasis(3, np.flip(knot, 0)), cp, True) else: result = Curve(BSplineBasis(3, knot), cp, True) @@ -321,7 +349,7 @@ def circle_segment(theta, r=1, center=(0, 0, 0), normal=(0, 0, 1), xaxis=(1, 0, return flip_and_move_plane_geometry(result, center, normal) -def interpolate(x, basis, t=None): +def interpolate(x: Points, basis: BSplineBasis, t: Knots | None = None) -> Curve: """Perform general spline interpolation on a provided basis. :param matrix-like x: Matrix *X[i,j]* of interpolation points *xi* with @@ -333,21 +361,20 @@ def interpolate(x, basis, t=None): :rtype: Curve """ # wrap input into an array - x = np.array(x) + points = np.asarray(x) # evaluate all basis functions in the interpolation points - if t is None: - t = basis.greville() + t = basis.greville() if t is None else np.asarray(t) N = basis.evaluate(t, sparse=True) # solve interpolation problem - cp = splinalg.spsolve(N, x) - cp = cp.reshape(x.shape) + cp = splinalg.spsolve(N, points) + cp = cp.reshape(points.shape) return Curve(basis, cp) -def least_square_fit(x, basis, t): +def least_square_fit(x: Points, basis: BSplineBasis, t: Knots) -> Curve: """Perform a least-square fit of a point cloud onto a spline basis :param matrix-like x: Matrix *X[i,j]* of interpolation points *xi* with @@ -363,12 +390,17 @@ def least_square_fit(x, basis, t): N = basis.evaluate(t) # solve interpolation problem - controlpoints, _, _, _ = np.linalg.lstsq(N, x, rcond=None) + controlpoints, _, _, _ = np.linalg.lstsq(N, np.asarray(x), rcond=None) return Curve(basis, controlpoints) -def cubic_curve(x, boundary=Boundary.FREE, t=None, tangents=None): +def cubic_curve( + x: Points, + boundary: Boundary = Boundary.FREE, + t: Knots | None = None, + tangents: Points | None = None, +) -> Curve: """Perform cubic spline interpolation on a provided basis. The valid boundary conditions are enumerated in :class:`Boundary`. The @@ -389,6 +421,7 @@ def cubic_curve(x, boundary=Boundary.FREE, t=None, tangents=None): :return: Interpolated curve :rtype: Curve """ + x = np.asarray(x) # if periodic input is not closed, make sure we do it now if boundary == Boundary.PERIODIC and not ( @@ -404,44 +437,40 @@ def cubic_curve(x, boundary=Boundary.FREE, t=None, tangents=None): # augment interpolation knot by euclidian distance to end t = list(t) + [t[-1] + norm(np.array(x[0, :]) - np.array(x[-2, :]))] - len(x) if t is None: - t = [0.0] - for x0, x1 in zip(x[:-1, :], x[1:, :]): - # eucledian distance between two consecutive points - dist = norm(np.array(x1) - np.array(x0)) - t.append(t[-1] + dist) + t = curve_length_parametrization(x) # modify knot vector for chosen boundary conditions - knot = [t[0]] * 3 + list(t) + [t[-1]] * 3 + t_array = np.asarray(t) + knot = with_repeated_knots(t_array, reps=4) if boundary == Boundary.FREE: - del knot[-5] - del knot[4] + knot = np.delete(knot, [4, -5]) elif boundary == Boundary.HERMITE: - knot = sorted(list(knot) + list(t[1:-1])) + knot = np.sort(np.concatenate((knot, t_array[1:-1]))) # create the interpolation basis and interpolation matrix on this if boundary == Boundary.PERIODIC: # C2-periodic knots - knot[0] = t[0] + t[-4] - t[-1] - knot[1] = t[0] + t[-3] - t[-1] - knot[2] = t[0] + t[-2] - t[-1] - knot[-3] = t[-1] + t[1] - t[0] - knot[-2] = t[-1] + t[2] - t[0] - knot[-1] = t[-1] + t[3] - t[0] + knot[0] = t_array[0] + t_array[-4] - t_array[-1] + knot[1] = t_array[0] + t_array[-3] - t_array[-1] + knot[2] = t_array[0] + t_array[-2] - t_array[-1] + knot[-3] = t_array[-1] + t_array[1] - t_array[0] + knot[-2] = t_array[-1] + t_array[2] - t_array[0] + knot[-1] = t_array[-1] + t_array[3] - t_array[0] basis = BSplineBasis(4, knot, 2) - # do not duplicate the interpolation at the sem (start=end is the same point) + # do not duplicate the interpolation at the seam (start=end is the same point) # identical points equal singular interpolation matrix - t = t[:-1] + t_array = t_array[:-1] x = x[:-1, :] else: basis = BSplineBasis(4, knot) - N = basis(t, sparse=True) # left-hand-side matrix + N = basis(t_array, sparse=True) # left-hand-side matrix # add derivative boundary conditions if applicable if boundary in [Boundary.TANGENT, Boundary.HERMITE, Boundary.TANGENTNATURAL]: + assert tangents is not None if boundary == Boundary.TANGENT: dn = basis([t[0], t[-1]], d=1) elif boundary == Boundary.TANGENTNATURAL: @@ -470,7 +499,7 @@ def cubic_curve(x, boundary=Boundary.FREE, t=None, tangents=None): return Curve(basis, cp) -def bezier(pts, quadratic=False, relative=False): +def bezier(pts: Points, quadratic: bool = False, relative: bool = False) -> Curve: """Generate a cubic or quadratic bezier curve from a set of control points :param [array-like] pts: list of control-points. In addition to a starting @@ -485,19 +514,24 @@ def bezier(pts, quadratic=False, relative=False): """ p = 3 if quadratic else 4 + points = np.asarray(pts) + # compute number of intervals - n = int((len(pts) - 1) / (p - 1)) + n = int((len(points) - 1) / (p - 1)) # generate uniform knot vector of repeated integers - knot = list(range(n + 1)) * (p - 1) + [0, n] - knot.sort() + knot = knot_vector(num_intervals=n, interior_reps=p - 1, endpoint_reps=p) + if relative: - pts = copy.deepcopy(pts) - for i in range(1, len(pts)): - pts[i] = [x0 + x1 for (x0, x1) in zip(pts[i - 1], pts[i])] - return Curve(BSplineBasis(p, knot), pts) + points = np.cumsum(points, axis=0) + return Curve(BSplineBasis(p, knot), points) -def manipulate(crv, f, normalized=False, vectorized=False): +def manipulate( + crv: Curve, + f: Callable, + normalized: bool = False, + vectorized: bool = False, +) -> Curve: """Create a new curve based on an expression-evaluation of an existing one :param Curve crv: original curve on which f is to be applied :param function f: expression of the physical point *x*, the velocity @@ -540,9 +574,12 @@ def move_3_to_right_fast(x, v): t = np.array(b.greville()) n = len(crv) + argv: list[Any] + destination: FloatArray + if vectorized: x = crv(t) - arg_names = inspect.getargspec(f).args + arg_names = inspect.getfullargspec(f).args argc = len(arg_names) argv = [0] * argc for j in range(argc): @@ -567,11 +604,12 @@ def move_3_to_right_fast(x, v): a[:] = [acc / norm(acc) for acc in a] argv[j] = a destination = f(*argv) + else: destination = np.zeros((len(crv), crv.dimension)) for t1, i in zip(t, range(len(t))): x = crv(t1) - arg_names = inspect.getargspec(f).args + arg_names = inspect.getfullargspec(f).args argc = len(arg_names) argv = [0] * argc for j in range(argc): @@ -602,7 +640,13 @@ def move_3_to_right_fast(x, v): return Curve(b, controlpoints) -def fit(x, t0, t1, rtol=1e-4, atol=0.0): +def fit( + x: Callable[..., FloatArray], + t0: Scalar, + t1: Scalar, + rtol: Scalar = 1e-4, + atol: Scalar = 0.0, +) -> Curve: """Computes an interpolation for a parametric curve up to a specified tolerance. The method will iteratively refine parts where needed resulting in a non-uniform knot vector with as optimized knot locations as possible. @@ -698,7 +742,7 @@ def move_along_tangent(t): return crv -def fit_points(x, t=[], rtol=1e-4, atol=0.0): +def fit_points(x: Points, t: Knots | None = None, rtol: Scalar = 1e-4, atol: Scalar = 0.0) -> Curve: """Computes an approximation for a list of points up to a specified tolerance. The method will iteratively refine parts where needed resulting in a non-uniform knot vector with as optimized knot locations as possible. The target curve is the @@ -714,5 +758,5 @@ def fit_points(x, t=[], rtol=1e-4, atol=0.0): :return: Curve Non-uniform cubic B-spline curve """ - linear = polygon(x, t=t) if len(t) > 0 else polygon(x) + linear = polygon(x, t=t) return fit(linear, linear.start(0), linear.end(0), rtol=rtol, atol=atol) diff --git a/src/splipy/io/threedm.py b/src/splipy/io/threedm.py index 03ddb1f..b6c9cee 100644 --- a/src/splipy/io/threedm.py +++ b/src/splipy/io/threedm.py @@ -34,6 +34,7 @@ from types import TracebackType from splipy.splineobject import SplineObject + from splipy.typing import Point # The rhino3dm type hints are incomplete, hence we have some shims. @@ -85,7 +86,7 @@ def read(self) -> list[SplineObject]: result.append(self.read_surface(nsrf)) if type(geom) is Line: - result.append(curve_factory.line(geom.From, geom.To)) + result.append(curve_factory.line(cast("Point", geom.From), cast("Point", geom.To))) continue if type(geom) is PolylineCurve: geom = geom.ToPolyline() diff --git a/src/splipy/splineobject.py b/src/splipy/splineobject.py index cad5002..11748ce 100644 --- a/src/splipy/splineobject.py +++ b/src/splipy/splineobject.py @@ -23,6 +23,8 @@ ) if TYPE_CHECKING: + from splipy.typing import ControlPoints, Point + from .typing import ArrayLike, Direction, FloatArray, Scalar, SectionElement, SectionKwargs __all__ = ["SplineObject"] @@ -67,7 +69,7 @@ class SplineObject: @staticmethod def construct_subclass( bases: Sequence[BSplineBasis], - controlpoints: ArrayLike, + controlpoints: ControlPoints, rational: bool, raw: bool = True, ) -> SplineObject: @@ -85,7 +87,7 @@ def construct_self(cls, bases: Sequence[BSplineBasis], controlpoints: FloatArray def __init__( self, bases: Sequence[BSplineBasis | None], - controlpoints: ArrayLike | None = None, + controlpoints: ControlPoints | None = None, rational: bool = False, raw: bool = False, ) -> None: @@ -123,7 +125,7 @@ def __init__( self.controlpoints = cps else: - self.controlpoints = np.array(controlpoints, dtype=np.float64) + self.controlpoints = np.asarray(controlpoints, dtype=np.float64) self.dimension = self.controlpoints.shape[-1] - rational self.rational = rational @@ -368,7 +370,7 @@ def get_antiderivative_spline( def get_antiderivative_spline( self, direction: Direction | None = None, constant: ArrayLike | None = None - ) -> SplineObject: + ) -> SplineObject | list[SplineObject]: """Compute the antiderivative (integral) of the spline object in a given parametric direction. The antiderivative is computed by inverting the derivative operator on @@ -1086,12 +1088,12 @@ def translate(self, x: ArrayLike) -> Self: return self @overload - def scale(self, arg: ArrayLike, /) -> Self: ... + def scale(self, arg: Point, /) -> Self: ... @overload def scale(self, arg: Scalar, *args: Scalar) -> Self: ... - def scale(self, arg: ArrayLike | Scalar, *args: Scalar) -> Self: + def scale(self, arg: Point | Scalar, *args: Scalar) -> Self: """Scale, or magnify the object by a given amount. In case of one input argument, the scaling is uniform. @@ -1621,46 +1623,46 @@ def shape(self) -> tuple[int, ...]: """The dimensions of the control point array.""" return self.controlpoints.shape[:-1] - def __iadd__(self, x: ArrayLike) -> Self: + def __iadd__(self, x: Point) -> Self: self.translate(x) return self - def __isub__(self, x: ArrayLike) -> Self: + def __isub__(self, x: Point) -> Self: self.translate(-np.asarray(x, dtype=np.float64)) # can't do -x if x is a list, so we rewrap it here return self - def __imul__(self, x: ArrayLike | Scalar) -> Self: + def __imul__(self, x: Point | Scalar) -> Self: self.scale(x) return self - def __itruediv__(self, x: ArrayLike | Scalar) -> Self: + def __itruediv__(self, x: Point | Scalar) -> Self: self.scale(1.0 / np.asarray(x, dtype=np.float64)) return self __ifloordiv__ = __itruediv__ # integer division (should not distinguish) - def __add__(self, x: ArrayLike) -> Self: + def __add__(self, x: Point) -> Self: new_obj = copy.deepcopy(self) new_obj += x return new_obj - def __radd__(self, x: ArrayLike) -> Self: + def __radd__(self, x: Point) -> Self: return self + x - def __sub__(self, x: ArrayLike) -> Self: + def __sub__(self, x: Point) -> Self: new_obj = copy.deepcopy(self) new_obj -= x return new_obj - def __mul__(self, x: ArrayLike | Scalar) -> Self: + def __mul__(self, x: Point | Scalar) -> Self: new_obj = copy.deepcopy(self) new_obj *= x return new_obj - def __rmul__(self, x: ArrayLike | Scalar) -> Self: + def __rmul__(self, x: Point | Scalar) -> Self: return self * x - def __truediv__(self, x: ArrayLike | Scalar) -> Self: + def __truediv__(self, x: Point | Scalar) -> Self: new_obj = copy.deepcopy(self) new_obj /= x return new_obj diff --git a/src/splipy/surface_factory.py b/src/splipy/surface_factory.py index 332960b..7a6bee9 100644 --- a/src/splipy/surface_factory.py +++ b/src/splipy/surface_factory.py @@ -21,7 +21,7 @@ if TYPE_CHECKING: from collections.abc import Callable, Sequence - from splipy.typing import ArrayLike, FloatArray, Scalar + from splipy.typing import FloatArray, Knots, Point, Scalar __all__ = [ "square", @@ -41,7 +41,7 @@ ] -def square(size: Scalar = 1, lower_left: ArrayLike = (0, 0)) -> Surface: +def square(size: Scalar = 1, lower_left: Point = (0, 0)) -> Surface: """Create a square with parametric origin at *(0,0)*. :param float size: Size(s), either a single scalar or a tuple of scalars per axis @@ -57,10 +57,10 @@ def square(size: Scalar = 1, lower_left: ArrayLike = (0, 0)) -> Surface: def disc( r: Scalar = 1, - center: ArrayLike = (0, 0, 0), - normal: ArrayLike = (0, 0, 1), + center: Point = (0, 0, 0), + normal: Point = (0, 0, 1), type: Literal["radial", "square"] = "radial", - xaxis: ArrayLike = (1, 0, 0), + xaxis: Point = (1, 0, 0), ) -> Surface: """Create a circular disc. The *type* parameter distinguishes between different parametrizations. @@ -73,7 +73,6 @@ def disc( :return: The disc :rtype: Surface """ - r = float(r) if type == "radial": c1 = curve_factory.circle(r, center=center, normal=normal, xaxis=xaxis) c2 = flip_and_move_plane_geometry(c1 * 0, center, normal) @@ -83,17 +82,19 @@ def disc( return result if type == "square": w = 1 / sqrt(2) - cp: list[list[float]] = [ - [-r * w, -r * w, 1], - [0, -r, w], - [r * w, -r * w, 1], - [-r, 0, w], - [0, 0, 1], - [r, 0, w], - [-r * w, r * w, 1], - [0, r, w], - [r * w, r * w, 1], - ] + cp = np.array( + [ + [-r * w, -r * w, 1], + [0, -r, w], + [r * w, -r * w, 1], + [-r, 0, w], + [0, 0, 1], + [r, 0, w], + [-r * w, r * w, 1], + [0, r, w], + [r * w, r * w, 1], + ] + ) basis1 = BSplineBasis(3) basis2 = BSplineBasis(3) result = Surface(basis1, basis2, cp, True) @@ -103,9 +104,9 @@ def disc( def sphere( r: Scalar = 1, - center: ArrayLike = (0, 0, 0), - zaxis: ArrayLike = (0, 0, 1), - xaxis: ArrayLike = (1, 0, 0), + center: Point = (0, 0, 0), + zaxis: Point = (0, 0, 1), + xaxis: Point = (1, 0, 0), ) -> Surface: """Create a spherical shell. @@ -125,7 +126,7 @@ def sphere( return flip_and_move_plane_geometry(result, center, zaxis) -def extrude(curve: Curve, amount: ArrayLike) -> Surface: +def extrude(curve: Curve, amount: Point) -> Surface: """Extrude a curve by sweeping it to a given height. :param Curve curve: Curve to extrude @@ -144,7 +145,7 @@ def extrude(curve: Curve, amount: ArrayLike) -> Surface: return Surface(curve.bases[0], BSplineBasis(2), cp, curve.rational) -def revolve(curve: Curve, theta: Scalar = 2 * pi, axis: ArrayLike = (0, 0, 1)) -> Surface: +def revolve(curve: Curve, theta: Scalar = 2 * pi, axis: Point = (0, 0, 1)) -> Surface: """Revolve a surface by sweeping a curve in a rotational fashion around the *z* axis. @@ -194,9 +195,9 @@ def revolve(curve: Curve, theta: Scalar = 2 * pi, axis: ArrayLike = (0, 0, 1)) - def cylinder( r: Scalar = 1, h: Scalar = 1, - center: ArrayLike = (0, 0, 0), - axis: ArrayLike = (0, 0, 1), - xaxis: ArrayLike = (1, 0, 0), + center: Point = (0, 0, 0), + axis: Point = (0, 0, 1), + xaxis: Point = (1, 0, 0), ) -> Surface: """Create a cylinder shell with no top or bottom @@ -214,9 +215,9 @@ def cylinder( def torus( minor_r: Scalar = 1, major_r: Scalar = 3, - center: ArrayLike = (0, 0, 0), - normal: ArrayLike = (0, 0, 1), - xaxis: ArrayLike = (1, 0, 0), + center: Point = (0, 0, 0), + normal: Point = (0, 0, 1), + xaxis: Point = (1, 0, 0), ) -> Surface: """Create a torus (doughnut) by revolving a circle of size *minor_r* around the *z* axis with radius *major_r*. @@ -591,7 +592,7 @@ def finitestrain_patch(bottom: Curve, right: Curve, top: Curve, left: Curve) -> return srf -def thicken(curve: Curve, amount: Scalar | Callable[..., float]) -> Surface: +def thicken(curve: Curve, amount: Scalar | Callable[..., Scalar]) -> Surface: """Generate a surface by adding thickness to a curve. - For 2D curves this will generate a 2D planar surface with the curve @@ -821,7 +822,7 @@ def loft(*in_curves: Curve | Sequence[Curve]) -> Surface: def interpolate( x: FloatArray, bases: Sequence[BSplineBasis], - u: Sequence[ArrayLike] | None = None, + u: Sequence[Knots] | None = None, ) -> Surface: """Interpolate a surface on a set of regular gridded interpolation points `x`. @@ -851,7 +852,7 @@ def interpolate( return Surface(bases[0], bases[1], cp.transpose(1, 0, 2).reshape((np.prod(surf_shape), dim))) -def least_square_fit(x: FloatArray, bases: Sequence[BSplineBasis], u: Sequence[ArrayLike]) -> Surface: +def least_square_fit(x: FloatArray, bases: Sequence[BSplineBasis], u: Sequence[Knots]) -> Surface: """Perform a least-square fit of a point cloud `x` onto a spline basis. The points can be either a matrix (in which case the first index is diff --git a/src/splipy/utils/__init__.py b/src/splipy/utils/__init__.py index c7f9250..1eeb4bf 100644 --- a/src/splipy/utils/__init__.py +++ b/src/splipy/utils/__init__.py @@ -9,7 +9,16 @@ if TYPE_CHECKING: from splipy.splineobject import SplineObject - from splipy.typing import Direction, FloatArray, Section, SectionElement, SectionKwargs + from splipy.typing import ( + Direction, + FloatArray, + IntArray, + Point, + Points, + Section, + SectionElement, + SectionKwargs, + ) def knot_vector( @@ -40,6 +49,15 @@ def iter() -> Iterator[float]: return np.fromiter(iter(), dtype=np.float64, count=count) +def with_repeated_knots( + knots: FloatArray | IntArray, + reps: int = 1, +) -> FloatArray: + if reps == 1: + return knots + return np.pad(knots, pad_width=reps - 1, mode="edge") + + def is_right_hand(patch, tol=1e-3): param = tuple((a + b) / 2 for a, b in zip(patch.start(), patch.end())) @@ -187,6 +205,15 @@ def ensure_listlike_old(x, dups=1): return [] +def normalize_points(*points: Point | Points) -> FloatArray: + """Utility function for functions that accept a single sequence of points or + multiple points as parameters. + """ + if len(points) == 1: + return np.asarray(points[0]) + return np.asarray(points) + + def rotate_local_x_axis(xaxis=(1, 0, 0), normal=(0, 0, 1)): # rotate xaxis vector back to reference domain (r=1, around origin) theta = atan2(normal[1], normal[0]) diff --git a/src/splipy/utils/curve.py b/src/splipy/utils/curve.py index 14de1ef..85c014b 100644 --- a/src/splipy/utils/curve.py +++ b/src/splipy/utils/curve.py @@ -1,11 +1,13 @@ from __future__ import annotations +from itertools import repeat + __doc__ = "Implementation of various curve utilities" import numpy as np -def curve_length_parametrization(pts, normalize=False): +def curve_length_parametrization(pts, normalize=False, reps=1): """Calculate knots corresponding to a curvelength parametrization of a set of points. @@ -14,10 +16,12 @@ def curve_length_parametrization(pts, normalize=False): :return: The parametrization :rtype: [float] """ - knots = [0.0] + knots = [0.0] * reps for i in range(1, len(pts)): knots.append(knots[-1] + np.linalg.norm(pts[i] - pts[i - 1])) + knots.extend(repeat(knots[-1], reps - 1)) + if normalize: length = knots[-1] knots = [k / length for k in knots] diff --git a/src/splipy/volume_factory.py b/src/splipy/volume_factory.py index 75624cd..d016d62 100644 --- a/src/splipy/volume_factory.py +++ b/src/splipy/volume_factory.py @@ -19,7 +19,7 @@ from collections.abc import Sequence from splipy.curve import Curve - from splipy.typing import ArrayLike, FloatArray, Scalar + from splipy.typing import FloatArray, Point, Scalar __all__ = [ "cube", @@ -34,7 +34,7 @@ ] -def cube(size: Scalar = 1, lower_left: ArrayLike = (0, 0, 0)) -> Volume: +def cube(size: Scalar = 1, lower_left: Point = (0, 0, 0)) -> Volume: """Create a cube with parmetric origin at *(0,0,0)*. :param float size: Size(s), either a single scalar or a tuple of scalars per axis @@ -49,7 +49,9 @@ def cube(size: Scalar = 1, lower_left: ArrayLike = (0, 0, 0)) -> Volume: def sphere( - r: Scalar = 1, center: ArrayLike = (0, 0, 0), type: Literal["radial", "square"] = "radial" + r: Scalar = 1, + center: Point = (0, 0, 0), + type: Literal["radial", "square"] = "radial", ) -> Volume: """Create a solid sphere @@ -124,7 +126,7 @@ def sphere( raise ValueError("invalid type argument") -def revolve(surf: Surface, theta: Scalar = 2 * pi, axis: ArrayLike = (0, 0, 1)) -> Volume: +def revolve(surf: Surface, theta: Scalar = 2 * pi, axis: Point = (0, 0, 1)) -> Volume: """Revolve a volume by sweeping a surface in a rotational fashion around an axis. @@ -170,9 +172,9 @@ def revolve(surf: Surface, theta: Scalar = 2 * pi, axis: ArrayLike = (0, 0, 1)) def torus( minor_r: Scalar = 1, major_r: Scalar = 3, - center: ArrayLike = (0, 0, 0), - normal: ArrayLike = (0, 0, 1), - xaxis: ArrayLike = (1, 0, 0), + center: Point = (0, 0, 0), + normal: Point = (0, 0, 1), + xaxis: Point = (1, 0, 0), type: Literal["radial", "square"] = "radial", ) -> Volume: """Create a torus (doughnut) by revolving a circle of size *minor_r* @@ -200,9 +202,9 @@ def torus( def cylinder( r: Scalar = 1, h: Scalar = 1, - center: ArrayLike = (0, 0, 0), - axis: ArrayLike = (0, 0, 1), - xaxis: ArrayLike = (1, 0, 0), + center: Point = (0, 0, 0), + axis: Point = (0, 0, 1), + xaxis: Point = (1, 0, 0), type: Literal["radial", "square"] = "radial", ) -> Volume: """Create a solid cylinder @@ -219,7 +221,7 @@ def cylinder( return extrude(surface_factory.disc(r, center, axis, xaxis=xaxis, type=type), h * np.array(axis)) -def extrude(surf: Surface, amount: ArrayLike) -> Volume: +def extrude(surf: Surface, amount: Point) -> Volume: """Extrude a surface by sweeping it to a given height. :param Surface surf: Surface to extrude @@ -518,7 +520,9 @@ def loft(*srfs: Surface | Sequence[Surface]) -> Volume: def interpolate( - x: FloatArray, bases: Sequence[BSplineBasis], u: Sequence[FloatArray] | None = None + x: FloatArray, + bases: Sequence[BSplineBasis], + u: Sequence[FloatArray] | None = None, ) -> Volume: """Interpolate a volume on a set of regular gridded interpolation points `x`. @@ -548,7 +552,11 @@ def interpolate( return Volume(bases[0], bases[1], bases[2], cp.transpose(2, 1, 0, 3).reshape((np.prod(vol_shape), dim))) -def least_square_fit(x: FloatArray, bases: Sequence[BSplineBasis], u: Sequence[FloatArray]) -> Volume: +def least_square_fit( + x: FloatArray, + bases: Sequence[BSplineBasis], + u: Sequence[FloatArray], +) -> Volume: """Perform a least-square fit of a point cloud `x` onto a spline basis. The points can be either a matrix (in which case the first index is From bdc7a33280b293ac115606d402beede54449c8ff Mon Sep 17 00:00:00 2001 From: Eivind Fonn Date: Sat, 3 Jan 2026 22:08:45 +0100 Subject: [PATCH 06/22] Typing: utils/curve.py --- src/splipy/surface_factory.py | 5 +++-- src/splipy/utils/curve.py | 24 ++++++++++++++++-------- src/splipy/volume_factory.py | 5 +++-- 3 files changed, 22 insertions(+), 12 deletions(-) diff --git a/src/splipy/surface_factory.py b/src/splipy/surface_factory.py index 7a6bee9..c96c638 100644 --- a/src/splipy/surface_factory.py +++ b/src/splipy/surface_factory.py @@ -781,10 +781,11 @@ def loft(*in_curves: Curve | Sequence[Curve]) -> Surface: basis2 = BSplineBasis(3) else: # create knot vector from the euclidian length between the curves - dist = curve_length_parametrization([c.center() for c in curves]) + dist = curve_length_parametrization(np.asarray([c.center() for c in curves])) + knot = curve_length_parametrization(np.asarray([c.center() for c in curves]), reps=4) # using "free" boundary condition by setting N'''(u) continuous at second to last and second knot - knot = [dist[0]] * 4 + dist[2:-2] + [dist[-1]] * 4 + knot = np.delete(knot, [4, -5]) basis2 = BSplineBasis(4, knot) n = len(curves) diff --git a/src/splipy/utils/curve.py b/src/splipy/utils/curve.py index 85c014b..3ae6612 100644 --- a/src/splipy/utils/curve.py +++ b/src/splipy/utils/curve.py @@ -2,34 +2,42 @@ from itertools import repeat +from splipy.curve import Curve +from splipy.typing import FloatArray, Points + __doc__ = "Implementation of various curve utilities" import numpy as np -def curve_length_parametrization(pts, normalize=False, reps=1): +def curve_length_parametrization(pts: Points, normalize: bool = False, reps: int = 1) -> FloatArray: """Calculate knots corresponding to a curvelength parametrization of a set of points. :param numpy.array pts: A set of points :param bool normalize: Whether to normalize the parametrization + :param reps int: How many repetitions of the first and last knot to return. :return: The parametrization :rtype: [float] """ - knots = [0.0] * reps - for i in range(1, len(pts)): - knots.append(knots[-1] + np.linalg.norm(pts[i] - pts[i - 1])) + points = np.asarray(pts) + npts = points.shape[0] + knots = np.zeros((npts + (reps - 1) * 2,), dtype=float) + + distances = np.linalg.norm(points[1:, ...] - points[:-1, ...], axis=1) + distances = np.cumsum(distances) + knots[reps:reps-1+npts] = distances - knots.extend(repeat(knots[-1], reps - 1)) + if reps > 1: + knots[-reps+1:] = knots[-reps] if normalize: - length = knots[-1] - knots = [k / length for k in knots] + knots /= knots[-1] return knots -def get_curve_points(curve): +def get_curve_points(curve: Curve) -> FloatArray: """Evaluate the curve in all its knots. :param curve: The curve diff --git a/src/splipy/volume_factory.py b/src/splipy/volume_factory.py index d016d62..ae2eb4b 100644 --- a/src/splipy/volume_factory.py +++ b/src/splipy/volume_factory.py @@ -474,10 +474,11 @@ def loft(*srfs: Surface | Sequence[Surface]) -> Volume: basis3 = BSplineBasis(3) else: # create knot vector from the euclidian length between the surfaces - dist = curve_length_parametrization([s.center() for s in surfaces]) + dist = curve_length_parametrization(np.asarray([s.center() for s in surfaces])) + knot = curve_length_parametrization(np.asarray([s.center() for s in surfaces]), reps=4) # using "free" boundary condition by setting N'''(u) continuous at second to last and second knot - knot = [dist[0]] * 4 + dist[2:-2] + [dist[-1]] * 4 + knot = np.delete(knot, [4, -5]) basis3 = BSplineBasis(4, knot) n = len(surfaces) From ad5df7ab72ccd220ab5aa4ae201a6e9c9f0a8a82 Mon Sep 17 00:00:00 2001 From: Eivind Fonn Date: Sat, 3 Jan 2026 22:20:07 +0100 Subject: [PATCH 07/22] Typing: curve.py --- src/splipy/basis.py | 7 ++++--- src/splipy/curve.py | 8 ++++++-- src/splipy/typing.py | 2 ++ 3 files changed, 12 insertions(+), 5 deletions(-) diff --git a/src/splipy/basis.py b/src/splipy/basis.py index 69d1442..19165cf 100644 --- a/src/splipy/basis.py +++ b/src/splipy/basis.py @@ -17,7 +17,7 @@ if TYPE_CHECKING: from splipy.typing import Knots - from .typing import FloatArray, Scalar + from .typing import FloatArray, Int, Scalar __all__ = ["BSplineBasis"] @@ -133,12 +133,13 @@ def greville_single(self, index: int) -> float: return float(np.sum(self.knots[index + 1 : index + self.order]) / (self.order - 1)) @overload - def greville(self, index: int) -> float: ... + def greville(self, index: Int) -> float: ... @overload def greville(self) -> FloatArray: ... - def greville(self, index: int | None = None) -> float | FloatArray: + + def greville(self, index: Int | None = None) -> float | FloatArray: """Fetch greville points, also known as knot averages: .. math:: \\sum_{j=i+1}^{i+p-1} \\frac{t_j}{p-1} diff --git a/src/splipy/curve.py b/src/splipy/curve.py index 35dc79e..f40f7f5 100644 --- a/src/splipy/curve.py +++ b/src/splipy/curve.py @@ -6,6 +6,8 @@ import numpy as np import scipy.sparse.linalg as splinalg +from splipy.typing import Point + from . import state from .basis import BSplineBasis from .splineobject import SplineObject @@ -433,7 +435,7 @@ def rebuild(self, p: int, n: int) -> Curve: # return new resampled curve return Curve(basis, controlpoints) - def _closest_point_linear_curve(self, pt: ArrayLike) -> tuple[FloatArray, float]: + def _closest_point_linear_curve(self, pt: Point) -> tuple[FloatArray, float]: """Computes the closest point on a linear curve to a given point. :param array-like pt: point to which the closest point on the curve is sought :return: the closest point on the curve and its parametric location @@ -455,7 +457,7 @@ def _closest_point_linear_curve(self, pt: ArrayLike) -> tuple[FloatArray, float] t = t1 return self(t), t - def closest_point(self, pt: ArrayLike, t0: Scalar = None) -> tuple[FloatArray, float]: + def closest_point(self, pt: Point, t0: Scalar | None = None) -> tuple[FloatArray, float]: """Computes the closest point on this curve to a given point. This is done by newton iteration and is using the state variables `controlpoint_absolute_tolerance` and `controlpoint_relative_tolerance` to determine convergence; but limited to 15 iterations. @@ -469,6 +471,8 @@ def closest_point(self, pt: ArrayLike, t0: Scalar = None) -> tuple[FloatArray, f if t0 is None: dist = [np.linalg.norm(cp - pt) for cp in self.controlpoints] + + pt = np.asarray(pt) i = np.argmin(dist) t0 = self.bases[0].greville(i) t = t0 diff --git a/src/splipy/typing.py b/src/splipy/typing.py index 28684ac..6b89108 100644 --- a/src/splipy/typing.py +++ b/src/splipy/typing.py @@ -29,6 +29,8 @@ # np.integer encompasses floats of many different sizes. type Scalar = float | np.floating | int | np.integer +type Int = int | np.integer + type FloatArray = npt.NDArray[np.floating] type IntArray = npt.NDArray[np.integer] From 333ed7a7195227ceeb34a4406fb8afa11ce4bc49 Mon Sep 17 00:00:00 2001 From: Eivind Fonn Date: Sat, 3 Jan 2026 23:44:24 +0100 Subject: [PATCH 08/22] Typing: splinemodel.py --- src/splipy/basis.py | 2 +- src/splipy/io/ofoam.py | 17 +- src/splipy/splinemodel.py | 317 +++++++++++++++++++++++------------ src/splipy/splineobject.py | 6 +- src/splipy/utils/__init__.py | 8 +- 5 files changed, 226 insertions(+), 124 deletions(-) diff --git a/src/splipy/basis.py b/src/splipy/basis.py index 19165cf..a5b30e9 100644 --- a/src/splipy/basis.py +++ b/src/splipy/basis.py @@ -120,7 +120,7 @@ def greville_all(self) -> FloatArray: operator = np.ones((self.order - 1,), dtype=np.float64) / (self.order - 1) return np.convolve(self.knots[1 : -1 - (self.periodic + 1)], operator, mode="valid") - def greville_single(self, index: int) -> float: + def greville_single(self, index: Int) -> float: """Fetch a greville point, also known as a knot averages: .. math:: \\sum_{j=i+1}^{i+p-1} \\frac{t_j}{p-1} diff --git a/src/splipy/io/ofoam.py b/src/splipy/io/ofoam.py index 9d64256..3a5e406 100644 --- a/src/splipy/io/ofoam.py +++ b/src/splipy/io/ofoam.py @@ -3,6 +3,8 @@ from itertools import groupby from operator import itemgetter from pathlib import Path +from types import TracebackType +from typing import Self import numpy as np @@ -12,10 +14,10 @@ class OpenFOAM: target: Path - def __init__(self, target): + def __init__(self, target: str) -> None: self.target = Path(target) - def __enter__(self): + def __enter__(self) -> Self: # Create the target directory if it does not exist if not self.target.exists(): self.target.mkdir(parents=True, exist_ok=True) @@ -25,10 +27,15 @@ def __enter__(self): return self - def __exit__(self, exc_type, exc_value, traceback): + def __exit__( + self, + exc_type: type[BaseException], + exc_value: BaseException, + traceback: TracebackType, + ) -> None: pass - def _header(self, cls, obj, note=None): + def _header(self, cls: str, obj: str, note: str | None = None) -> str: s = "FoamFile\n{\n" s += " version 2.0;\n" s += " format ascii;\n" @@ -39,7 +46,7 @@ def _header(self, cls, obj, note=None): s += "}\n" return s - def write(self, model): + def write(self, model: SplineModel) -> None: assert isinstance(model, SplineModel), "OpenFOAM.write only supports SplineModel objects" # Only linear volumes in 3D, please diff --git a/src/splipy/splinemodel.py b/src/splipy/splinemodel.py index ff12089..066ae95 100644 --- a/src/splipy/splinemodel.py +++ b/src/splipy/splinemodel.py @@ -2,12 +2,17 @@ from collections import Counter, OrderedDict, namedtuple from collections.abc import Callable, Iterator +from dataclasses import dataclass from itertools import chain, islice, permutations, product from operator import itemgetter from pathlib import Path -from typing import Any +from re import S +from typing import Any, Sequence, Unpack, cast, overload import numpy as np +import numpy.typing as npt + +from splipy.typing import ControlPoints, FloatArray, IntArray, Scalar, Section, Int, SectionElement, SectionKwargs from . import state from .splineobject import SplineObject @@ -27,17 +32,38 @@ from collections.abc import MutableMapping -def _section_to_index(section): +def _section_to_index(section: Section) -> tuple[slice | int]: """Replace all `None` in `section` with `slice(None)`, so that it works as a numpy array indexing tuple. """ - return tuple(slice(None) if s is None else s for s in section) + return tuple(slice(None) if s is None else s for s in section) # type: ignore[return-value] + + +class FaceScalar(np.record): + nodes: np.ndarray[tuple[int], np.dtype[np.int64]] + owner: np.int64 + neighbor: np.int64 + name: object -face_t = np.dtype([("nodes", int, (4,)), ("owner", int, ()), ("neighbor", int, ()), ("name", object, ())]) +class FaceArray(np.recarray[Any, np.dtype[FaceScalar]]): + nodes: np.ndarray[Any, np.dtype[np.int64]] + owner: np.ndarray[Any, np.dtype[np.int64]] + neighbor: np.ndarray[Any, np.dtype[np.int64]] + name: np.ndarray[Any, np.dtype[np.object_]] -class VertexDict(MutableMapping): +face_t = np.dtype( + [ + ("nodes", int, (4,)), + ("owner", int, ()), + ("neighbor", int, ()), + ("name", object, ()), + ], +) + + +class VertexDict[T](MutableMapping[FloatArray, T]): """A dictionary where the keys are numpy arrays, and where equality is computed in an approximate sense for floating point numbers. @@ -47,20 +73,19 @@ class VertexDict(MutableMapping): rtol: float atol: float - _keys: list[np.ndarray | None] - _values: list[Any] + _keys: list[FloatArray | None] + _values: list[T | None] lut: dict[tuple[int, ...], list[tuple[int, float]]] - def __init__(self, rtol=1e-5, atol=1e-8): - # List of (key, value) pairs + def __init__(self, rtol: float = 1e-5, atol: float = 1e-8) -> None: self.rtol = rtol self.atol = atol self._keys = [] self._values = [] self.lut = {} - def _bounds(self, key): + def _bounds(self, key: Scalar) -> tuple[Scalar, Scalar]: if key >= self.atol: return ( (key - self.atol) / (1 + self.rtol), @@ -78,14 +103,14 @@ def _bounds(self, key): (key + self.atol) / (1 - self.rtol), ) - def _candidate(self, key): + def _candidate(self, key: FloatArray) -> int: """Return the internal index for the first stored mapping that matches the given key. :param numpy.array key: The key to look for :raises KeyError: If the key is not found """ - candidates = None + candidates: set[int] | None = None for coord, k in np.ndenumerate(key): lut = self.lut.setdefault(coord, []) minval, maxval = self._bounds(k) @@ -95,12 +120,14 @@ def _candidate(self, key): candidates = {i for i, _ in lut[lo:hi]} else: candidates &= {i for i, _ in lut[lo:hi]} + if not candidates: + raise KeyError(key) for c in candidates: if self._keys[c] is not None: return c raise KeyError(key) - def _insert(self, key, value): + def _insert(self, key: FloatArray, value: T) -> None: newindex = len(self._values) for coord, v in np.ndenumerate(key): lut = self.lut.setdefault(coord, []) @@ -108,7 +135,7 @@ def _insert(self, key, value): self._keys.append(key) self._values.append(value) - def __setitem__(self, key, value): + def __setitem__(self, key: FloatArray, value: T) -> None: """Assign a key to a value.""" try: c = self._candidate(key) @@ -116,15 +143,17 @@ def __setitem__(self, key, value): except KeyError: self._insert(key, value) - def __getitem__(self, key): + def __getitem__(self, key: FloatArray) -> T: """Gets the value assigned to a key. :raises KeyError: If the key is not found """ c = self._candidate(key) - return self._values[c] + value = self._values[c] + assert value is not None + return value - def __delitem__(self, key): + def __delitem__(self, key: FloatArray) -> None: """Deletes an assignment.""" try: i = self._candidate(key) @@ -133,18 +162,16 @@ def __delitem__(self, key): self._keys[i] = None self._values[i] = None - def __iter__(self): + def __iter__(self) -> Iterator[FloatArray]: """Iterate over all keys. .. note:: This generates all the stored keys, not all matching keys. """ - yield from self._keys + for k in self._keys: + if k is not None: + yield k - def items(self): - """Return a list of key, value pairs.""" - yield from self._values - - def __len__(self): + def __len__(self) -> int: """Returns the number of stored assignments.""" return len(self._values) @@ -154,7 +181,6 @@ class OrientationError(RuntimeError): :class:`splipy.SplineModel.Orientation` indicating an inability to match two objects. """ - pass @@ -162,7 +188,6 @@ class TwinError(RuntimeError): """A `TwinError` is raised when two objects with identical interfaces are added, but different interiors. """ - pass @@ -178,7 +203,11 @@ class Orientation: direction `d` *in the reference system* should be reversed. """ - def __init__(self, perm, flip): + perm: Sequence[int] + perm_inv: Sequence[int] + flip: Sequence[bool] + + def __init__(self, perm: Sequence[int], flip: Sequence[bool]) -> None: """Initialize an Orientation object. .. warning:: This constructor is for internal use. Use @@ -190,7 +219,7 @@ def __init__(self, perm, flip): self.perm_inv = tuple(perm.index(d) for d in range(len(perm))) @classmethod - def compute(cls, cpa, cpb=None): + def compute(cls, cpa: SplineObject, cpb: SplineObject | None = None) -> Orientation: """Compute and return a new orientation object representing the mapping between `cpa` (the reference system) and `cpb` (the mapped system). @@ -257,10 +286,10 @@ def compute(cls, cpa, cpb=None): raise OrientationError("Non-matching objects") @property - def pardim(self): + def pardim(self) -> int: return len(self.perm) - def __mul__(self, other): + def __mul__(self, other: Orientation) -> Orientation: """Compose two mappings. If `ort_left` maps system `A` (reference) to system `B`, and @@ -275,13 +304,13 @@ def __mul__(self, other): return Orientation(perm, flip) - def map_array(self, array): + def map_array[T: np.generic](self, array: npt.NDArray[T]) -> npt.NDArray[T]: """Map an array in the mapped system to the reference system.""" array = array.transpose(*self.perm) flips = tuple(slice(None, None, -1) if f else slice(None) for f in self.flip) return array[flips] - def map_section(self, section): + def map_section(self, section: Section) -> Section: """Map a section in the mapped system to the reference system. The input is a section tuple as described in @@ -291,11 +320,8 @@ def map_section(self, section): """ permuted = tuple(section[d] for d in self.perm) - flipped = () - for ( - s, - f, - ) in zip(permuted, self.flip): + flipped: Section = () + for s, f in zip(permuted, self.flip): # Flipping only applies to indexed directions, not variable ones if f and s is not None: flipped += (0 if s == -1 else -1,) @@ -304,7 +330,7 @@ def map_section(self, section): return flipped - def view_section(self, section): + def view_section(self, section: Section) -> Orientation: """Reduce a mapping to a lower dimension. The input is a section tuple as described in @@ -332,7 +358,7 @@ def view_section(self, section): return self.__class__(new_perm, new_flip) @property - def ifem_format(self): + def ifem_format(self) -> int: """Compute the orientation in IFEM format. For one-dimensional objects, this is a single binary digit indicating @@ -390,7 +416,22 @@ class TopologicalNode: of any kind. """ - def __init__(self, obj, lower_nodes, index): + obj: SplineObject + index: int + lower_nodes: list[tuple[TopologicalNode, ...]] + higher_nodes: dict[int, list[TopologicalNode]] + owner: TopologicalNode | None + + name: str | None + cell_numbers: IntArray | None + cp_numbers: IntArray | None + + def __init__( + self, + obj: SplineObject, + lower_nodes: list[tuple[TopologicalNode, ...]], + index: int + ) -> None: """Initialize a `TopologicalNode` object associated with the given `SplineObject` and lower order nodes. @@ -419,26 +460,26 @@ def __init__(self, obj, lower_nodes, index): node._transfer_ownership(self) @property - def pardim(self): + def pardim(self) -> int: return self.obj.pardim @property - def nhigher(self): + def nhigher(self) -> int: return len(self.higher_nodes[self.pardim + 1]) @property - def super_owner(self): + def super_owner(self) -> TopologicalNode: """Return the highest owning node.""" owner = self while owner.owner is not None: owner = owner.owner return owner - def assign_higher(self, node): + def assign_higher(self, node: TopologicalNode) -> None: """Add a link to a node of higher dimension.""" self.higher_nodes.setdefault(node.pardim, []).append(node) - def view(self, other_obj=None): + def view(self, other_obj: SplineObject | None = None) -> NodeView: """Return a `NodeView` object of this node. The returned view has an orientation that matches that of the input @@ -451,7 +492,7 @@ def view(self, other_obj=None): orientation = Orientation.compute(self.obj, other_obj) if other_obj else Orientation.compute(self.obj) return NodeView(self, orientation) - def _transfer_ownership(self, new_owner): + def _transfer_ownership(self, new_owner: TopologicalNode) -> None: """Transfers ownership of this node to a new owner. This operation is transitive, so all child nodes owned by this node, or who are owner-less will also be transferred. @@ -465,7 +506,7 @@ def _transfer_ownership(self, new_owner): if child.owner is self or child.owner is None: child._transfer_ownership(new_owner) - def generate_cp_numbers(self, start=0): + def generate_cp_numbers(self, start: Int = 0) -> Int: """Generate a control point numbering starting at `start`. Return the next unused index.""" assert self.owner is None @@ -488,7 +529,7 @@ def generate_cp_numbers(self, start=0): self.assign_cp_numbers(numbers) return start + nowned - def assign_cp_numbers(self, numbers): + def assign_cp_numbers(self, numbers: IntArray) -> None: """Directly assign control point numbers.""" self.cp_numbers = numbers @@ -500,17 +541,22 @@ def assign_cp_numbers(self, numbers): # orientations not matching up. node.assign_cp_numbers(numbers[_section_to_index(section)]) - def read_cp_numbers(self): + def read_cp_numbers(self) -> None: """Read control point numbers for unowned control points from child nodes.""" + assert self.cp_numbers is not None + for node, section in zip(self.lower_nodes[-1], sections(self.pardim, self.pardim - 1)): if node.owner is not self: # The two sections may not agree on orientation, so we fix this here. - ori = Orientation.compute(self.obj.section(*section), node.obj) + sub = self.obj.section(*section) + assert isinstance(sub, SplineObject) + ori = Orientation.compute(sub, node.obj) + assert node.cp_numbers is not None self.cp_numbers[_section_to_index(section)] = ori.map_array(node.cp_numbers) assert (self.cp_numbers != -1).all() - def generate_cell_numbers(self, start=0): + def generate_cell_numbers(self, start: Int = 0) -> Int: """Generate a cell numbering starting at `start`. Return the next unused index.""" assert self.owner is None @@ -520,16 +566,19 @@ def generate_cell_numbers(self, start=0): self.cell_numbers = np.reshape(np.arange(start, start + nelems, dtype=int), shape) return start + nelems - def faces(self): + def faces(self) -> list[FaceArray]: """Return all faces owned by this node, as a list of numpy arrays with dtype `face_t`.""" assert self.pardim == 3 assert self.obj.order() == (2, 2, 2) + assert self.cp_numbers is not None + assert self.cell_numbers is not None + shape = [len(kvec) - 1 for kvec in self.obj.knots()] ncells = np.prod(shape) - retval = [] + retval: list[FaceArray] = [] - def mkindex(dim, z, a, b): - rval = [a, b] if dim != 1 else [b, a] + def mkindex(dim: int, z: slice | int, a: slice, b: slice) -> tuple[slice | int, ...]: + rval: list[slice | int] = [a, b] if dim != 1 else [b, a] rval.insert(dim, z) return tuple(rval) @@ -549,7 +598,7 @@ def mkindex(dim, z, a, b): faces["nodes"][:, 3] = self.cp_numbers[mkindex(d, np.s_[1:-1], np.s_[:-1], np.s_[1:])].flatten() faces["owner"] = self.cell_numbers[mkindex(d, np.s_[:-1], np.s_[:], np.s_[:])].flatten() faces["neighbor"] = self.cell_numbers[mkindex(d, np.s_[1:], np.s_[:], np.s_[:])].flatten() - retval.append(faces) + retval.append(faces.view(FaceArray)) # Go through the two boundaries for bdnode, bdindex in zip(islice(lower, 2), (0, -1)): @@ -588,16 +637,18 @@ def mkindex(dim, z, a, b): # Get the spline object on that interface as oriented from the neighbor's perspective nb_sec = section_from_index(3, 2, nb_index) nb_obj = neighbor.obj.section(*nb_sec) + assert isinstance(nb_obj, SplineObject) # Compute the relative orientation ori = Orientation.compute(bdnode.obj, nb_obj) # Get the neighbor cell numbers from the neighbor's # perspective, and map them to our system + assert neighbor.cell_numbers is not None cellidxs = neighbor.cell_numbers[_section_to_index(nb_sec)] faces["neighbor"] = ori.map_array(cellidxs).flatten() - retval.append(faces) + retval.append(faces.view(FaceArray)) for faces in retval: assert ((faces["owner"] < faces["neighbor"]) | (faces["neighbor"] == -1)).all() @@ -613,7 +664,10 @@ class NodeView: persistent. """ - def __init__(self, node, orientation=None): + node: TopologicalNode + orientation: Orientation + + def __init__(self, node: TopologicalNode, orientation: Orientation) -> None: """Initialize a `NodeView` object with the given node and orientation. .. warning:: This constructor is for internal use. @@ -622,18 +676,18 @@ def __init__(self, node, orientation=None): self.orientation = orientation @property - def pardim(self): + def pardim(self) -> int: return self.node.pardim @property - def name(self): + def name(self) -> str | None: return self.node.name @name.setter - def name(self, value): + def name(self, value: str | None) -> None: self.node.name = value - def section(self, *args, **kwargs): + def section(self, *args: SectionElement, **kwargs: Unpack[SectionKwargs]) -> NodeView: """Return a section. See :func:`splipy.SplineObject.section` for more details on the input arguments. @@ -657,32 +711,32 @@ def section(self, *args, **kwargs): return NodeView(node, ref_ori * my_ori) - def corner(self, i): + def corner(self, i: int) -> NodeView: """Return the i'th corner.""" return self.section(*section_from_index(self.pardim, 0, i)) @property - def corners(self): + def corners(self) -> tuple[NodeView, ...]: """A tuple of all corners.""" - return tuple(self.section(s) for s in sections(self.pardim, 0)) + return tuple(self.section(*s) for s in sections(self.pardim, 0)) - def edge(self, i): + def edge(self, i: int) -> NodeView: """Return the i'th edge.""" return self.section(*section_from_index(self.pardim, 1, i)) @property - def edges(self): + def edges(self) -> tuple[NodeView, ...]: """A tuple of all edges.""" - return tuple(self.section(s) for s in sections(self.pardim, 1)) + return tuple(self.section(*s) for s in sections(self.pardim, 1)) - def face(self, i): + def face(self, i: int) -> NodeView: """Return the i'th face.""" return self.section(*section_from_index(self.pardim, 2, i)) @property - def faces(self): + def faces(self) -> tuple[NodeView, ...]: """A tuple of all faces.""" - return tuple(self.section(s) for s in sections(self.pardim, 2)) + return tuple(self.section(*s) for s in sections(self.pardim, 2)) class ObjectCatalogue: @@ -690,7 +744,13 @@ class ObjectCatalogue: at most `pardim` parametric directions. """ - def __init__(self, pardim): + pardim: int + count: int + internal: OrderedDict[tuple[TopologicalNode, ...], list[TopologicalNode]] + lower: ObjectCatalogue | VertexDict[TopologicalNode] + callbacks: dict[str, list[Callable[[TopologicalNode], None]]] + + def __init__(self, pardim: int) -> None: """Initialize a catalogue for objects of parametric dimension `pardim`. """ @@ -710,11 +770,11 @@ def __init__(self, pardim): # Callbacks for events self.callbacks = {} - def add_callback(self, event: str, callback: Callable[[TopologicalNode], None]): + def add_callback(self, event: str, callback: Callable[[TopologicalNode], None]) -> None: """Add a callback function to be called on a given event.""" self.callbacks.setdefault(event, []).append(callback) - def lookup(self, obj, add=False, raise_on_twins=()): + def lookup(self, obj: SplineObject, add: bool = False, raise_on_twins: Sequence[int] = ()) -> NodeView: """Obtain the `NodeView` object corresponding to a given object. If the keyword argument `add` is true, this function may generate one @@ -738,21 +798,25 @@ def lookup(self, obj, add=False, raise_on_twins=()): """ # Pass lower-dimensional objects through to the lower levels if self.pardim > obj.pardim: + assert isinstance(self.lower, ObjectCatalogue) return self.lower.lookup(obj, add=add, raise_on_twins=raise_on_twins) # Special case for points: self.lower is a mapping from array to node if self.pardim == 0: + lower = cast("VertexDict[TopologicalNode]", self.lower) cps = obj.controlpoints if obj.rational: cps = cps[..., :-1] if add: node = TopologicalNode(obj, [], index=self.count) self.count += 1 - rval = self.lower.setdefault(cps, node).view() + rval = lower.setdefault(cps, node).view() for cb in self.callbacks.get("add", []): cb(node) return rval - return self.lower[cps].view() + return lower[cps].view() + + assert isinstance(self.lower, ObjectCatalogue) # Get all nodes of lower dimension (points, vertices, etc.) # This involves a recursive call to self.lower.__call__ @@ -811,7 +875,7 @@ def lookup(self, obj, add=False, raise_on_twins=()): raise KeyError("No such object found") return self._add(obj, lower_nodes) - def add(self, obj, raise_on_twins=()): + def add(self, obj: SplineObject, raise_on_twins: Sequence[int] = ()) -> NodeView: """Add new nodes to the graph to accommodate the given object, then return the corresponding `NodeView` object. @@ -836,7 +900,7 @@ def add(self, obj, raise_on_twins=()): """ return self.lookup(obj, add=True, raise_on_twins=raise_on_twins) - def _add(self, obj, lower_nodes): + def _add(self, obj: SplineObject, lower_nodes: list[tuple[TopologicalNode, ...]]) -> NodeView: node = TopologicalNode(obj, lower_nodes, index=self.count) self.count += 1 # Assign the new node to each possible permutation of lower-order @@ -852,16 +916,19 @@ def _add(self, obj, lower_nodes): __call__ = add __getitem__ = lookup - def top_nodes(self): + def top_nodes(self) -> list[TopologicalNode]: """Return all nodes of the highest parametric dimension.""" return self.nodes(self.pardim) - def nodes(self, pardim): + def nodes(self, pardim: int) -> list[TopologicalNode]: """Return all nodes of a given parametric dimension.""" if self.pardim == pardim: if self.pardim > 0: + assert isinstance(self.lower, ObjectCatalogue) return list(uniquify(chain.from_iterable(self.internal.values()))) - return list(uniquify(self.lower.values())) + lower = cast("VertexDict[TopologicalNode]", self.lower) + return list(uniquify(lower.values())) + assert isinstance(self.lower, ObjectCatalogue) return self.lower.nodes(pardim) @@ -870,7 +937,19 @@ def nodes(self, pardim): class SplineModel: - def __init__(self, pardim=3, dimension=3, objs=[], force_right_hand=False): + pardim: int + dimension: int + force_right_hand: bool + cataloge: ObjectCatalogue + names: dict[str, SplineObject] + + def __init__( + self, + pardim: int = 3, + dimension: int = 3, + objs: Sequence[SplineObject] | None = None, + force_right_hand: bool = False, + ) -> None: self.pardim = pardim self.dimension = dimension @@ -882,45 +961,49 @@ def __init__(self, pardim=3, dimension=3, objs=[], force_right_hand=False): self.catalogue = ObjectCatalogue(pardim) self.names = {} - self.add(objs) + if objs is not None: + self.add(objs) - def add_callback(self, event: str, callback: Callable[[TopologicalNode], None]): - catalogue = self.catalogue + def add_callback(self, event: str, callback: Callable[[TopologicalNode], None]) -> None: + catalogue: ObjectCatalogue | VertexDict[TopologicalNode] = self.catalogue while isinstance(catalogue, ObjectCatalogue): catalogue.add_callback(event, callback) catalogue = catalogue.lower - def add(self, obj, name=None, raise_on_twins=True): + def add( + self, + obj: SplineObject | Sequence[SplineObject], + name: str | None = None, + raise_on_twins: bool | Sequence[int] = True, + ) -> None: if raise_on_twins is True: raise_on_twins = tuple(range(self.pardim + 1)) elif raise_on_twins is False: raise_on_twins = () - if isinstance(obj, SplineObject): - obj = [obj] - self._validate(obj) - self._generate(obj, raise_on_twins=raise_on_twins) + self._validate([obj] if isinstance(obj, SplineObject) else obj) + self._generate([obj] if isinstance(obj, SplineObject) else obj, raise_on_twins=raise_on_twins) if name and isinstance(obj, SplineObject): self.names[name] = obj - def __getitem__(self, obj): + def __getitem__(self, obj: SplineObject) -> NodeView: return self.catalogue[obj] def objects(self) -> Iterator[SplineObject]: for node in self.catalogue.top_nodes(): yield node.obj - def boundary(self, name=None): + def boundary(self, name: str | None = None) -> Iterator[TopologicalNode]: for node in self.catalogue.nodes(self.pardim - 1): if node.nhigher == 1 and (name is None or name == node.name): yield node - def assign_boundary(self, name): + def assign_boundary(self, name: str) -> None: """Give a name to all unnamed boundary nodes.""" for node in self.boundary(): if node.name is None: node.name = name - def _validate(self, objs): + def _validate(self, objs: Sequence[SplineObject]) -> None: if any(p.dimension != self.dimension for p in objs): raise ValueError("Patches with different dimension added") if any(p.pardim > self.pardim for p in objs): @@ -931,10 +1014,10 @@ def _validate(self, objs): indices = ", ".join(map(str, left_inds)) raise ValueError(f"Possibly left-handed patches detected, indexes {indices}") - def _generate(self, objs, **kwargs): + def _generate(self, objs: Sequence[SplineObject], raise_on_twins: Sequence[int]) -> None: for i, p in enumerate(objs): try: - self.catalogue.add(p, **kwargs) + self.catalogue.add(p, raise_on_twins=raise_on_twins) except OrientationError as err: # TODO(Eivind): Mutating exceptions is fishy. if len(err.args) > 1: @@ -944,55 +1027,66 @@ def _generate(self, objs, **kwargs): ) raise err - def generate_cp_numbers(self): - index = 0 + def generate_cp_numbers(self) -> None: + index: Int = 0 for node in self.catalogue.top_nodes(): index = node.generate_cp_numbers(index) self.ncps = index for node in self.catalogue.top_nodes(): node.read_cp_numbers() - def generate_cell_numbers(self): - index = 0 + def generate_cell_numbers(self) -> None: + index: Int = 0 for node in self.catalogue.top_nodes(): index = node.generate_cell_numbers(index) self.ncells = index - def cps(self): + def cps(self) -> ControlPoints: cps = np.zeros((self.ncps, self.dimension)) for node in self.catalogue.top_nodes(): + assert node.cp_numbers is not None indices = node.cp_numbers.reshape(-1) values = node.obj.controlpoints.reshape(-1, self.dimension) cps[indices] = values return cps - def faces(self): + def faces(self) -> FaceArray: assert self.pardim == 3 faces = list(chain.from_iterable(node.faces() for node in self.catalogue.top_nodes())) - return np.hstack(faces) + return np.hstack(faces).view(FaceArray) - def summary(self): - c = self.catalogue + def summary(self) -> None: + c: ObjectCatalogue | VertexDict[TopologicalNode] = self.catalogue while isinstance(c, ObjectCatalogue): print(f"Dim {c.pardim}: {len(c.top_nodes())}") c = c.lower - def write_ifem(self, filename): + def write_ifem(self, filename: str) -> None: IFEMWriter(self).write(filename) -IFEMConnection = namedtuple("IFEMConnection", ["master", "slave", "midx", "sidx", "orient"]) +@dataclass +class IFEMConnection: + master: int + slave: int + midx: int + sidx: int + orient: int class IFEMWriter: - def __init__(self, model): + model: SplineModel + nodes: list[TopologicalNode] + node_ids: dict[TopologicalNode, int] + + def __init__(self, model: SplineModel) -> None: self.model = model # List the nodes so that the order is deterministic self.nodes = list(model.catalogue.top_nodes()) self.node_ids = {node: i for i, node in enumerate(self.nodes)} - def connections(self): + def connections(self) -> Iterator[IFEMConnection]: p = self.model.pardim # For every object in the model... @@ -1037,7 +1131,7 @@ def connections(self): orient=orientation.ifem_format, ) - def write(self, filename): + def write(self, filename: str) -> None: lines = [ "", "", @@ -1067,11 +1161,12 @@ def write(self, filename): ) for name in names: - entries = {} + entries: dict[int, set[int]] = {} for node in self.model.catalogue.nodes(self.model.pardim - 1): if node.name != name: continue parent = node.owner + assert parent is not None sub_idx = next( idx for idx, sub in enumerate(parent.lower_nodes[self.model.pardim - 1]) if sub is node ) diff --git a/src/splipy/splineobject.py b/src/splipy/splineobject.py index 11748ce..4187964 100644 --- a/src/splipy/splineobject.py +++ b/src/splipy/splineobject.py @@ -451,8 +451,8 @@ def get_antiderivative_spline( delta_knot = old_knots[i + p] - old_knots[i] # Build index slices: [..., i, ...] and [..., i+1, ...] - idx_current = [slice(None)] * self.pardim - idx_next = [slice(None)] * self.pardim + idx_current: list[slice | int] = [slice(None)] * self.pardim + idx_next: list[slice | int] = [slice(None)] * self.pardim idx_current[d] = i idx_next[d] = i + 1 @@ -1038,7 +1038,7 @@ def reparam( return self - def translate(self, x: ArrayLike) -> Self: + def translate(self, x: Point) -> Self: """Translate (i.e. move) the object by a given distance. :param array-like x: The vector to translate by. diff --git a/src/splipy/utils/__init__.py b/src/splipy/utils/__init__.py index 1eeb4bf..99bf178 100644 --- a/src/splipy/utils/__init__.py +++ b/src/splipy/utils/__init__.py @@ -100,7 +100,7 @@ def rotation_matrix(theta, axis): ) -def sections(src_dim, tgt_dim): +def sections(src_dim: int, tgt_dim: int) -> Iterator[Section]: """Generate all boundary sections from a source dimension to a target dimension. For example, `sections(3,1)` generates all edges on a volume. @@ -116,10 +116,10 @@ def sections(src_dim, tgt_dim): args = [None] * src_dim for f, i in zip(fixed, indices[::-1]): args[f] = i - yield args + yield tuple(args) -def section_from_index(src_dim, tgt_dim, i): +def section_from_index(src_dim, tgt_dim, i) -> Section: """Return the i'th section from a source dimension to a target dimension. See :func:`splipy.Utils.sections` for more information. @@ -127,7 +127,7 @@ def section_from_index(src_dim, tgt_dim, i): for j, s in enumerate(sections(src_dim, tgt_dim)): if i == j: return s - return None + raise ValueError(f"No such section: {i} for dimensions {src_dim} and {tgt_dim}") def section_to_index(section): From f42a15dc4517f867d71c93fdb0ba1cb6a3b73bc4 Mon Sep 17 00:00:00 2001 From: Eivind Fonn Date: Sat, 3 Jan 2026 23:52:27 +0100 Subject: [PATCH 09/22] Typing: ofoam.py --- src/splipy/io/ofoam.py | 14 +++++++------- src/splipy/splinemodel.py | 2 +- 2 files changed, 8 insertions(+), 8 deletions(-) diff --git a/src/splipy/io/ofoam.py b/src/splipy/io/ofoam.py index 3a5e406..5f83fb2 100644 --- a/src/splipy/io/ofoam.py +++ b/src/splipy/io/ofoam.py @@ -8,7 +8,7 @@ import numpy as np -from splipy.splinemodel import SplineModel +from splipy.splinemodel import FaceArray, SplineModel class OpenFOAM: @@ -58,7 +58,7 @@ def write(self, model: SplineModel) -> None: model.generate_cp_numbers() model.generate_cell_numbers() faces = model.faces() - ninternal = sum(faces["name"] is None) + ninternal = sum(faces.name == None) note = ( f"nPoints: {model.ncps} nCells: {model.ncells} nFaces: {len(faces)} nInternalFaces: {ninternal}" ) @@ -69,11 +69,11 @@ def write(self, model: SplineModel) -> None: # - All faces in the same boundary must be contiguous # - Low number owners before high number owners # - Low number neighbors before high number neighbors - faces = list(faces) - faces = sorted(faces, key=itemgetter("neighbor")) - faces = sorted(faces, key=itemgetter("owner")) - faces = sorted(faces, key=lambda x: (x["name"] is not None, x["name"])) - faces = np.array(faces) + faces_list = list(faces) + faces_list = sorted(faces_list, key=itemgetter("neighbor")) + faces_list = sorted(faces_list, key=itemgetter("owner")) + faces_list = sorted(faces_list, key=lambda x: (x["name"] is not None, x["name"])) + faces = np.array(faces_list).view(FaceArray) # Write the points file (vertex coordinates) with (self.target / "points").open("w") as f: diff --git a/src/splipy/splinemodel.py b/src/splipy/splinemodel.py index 066ae95..b91014c 100644 --- a/src/splipy/splinemodel.py +++ b/src/splipy/splinemodel.py @@ -1041,7 +1041,7 @@ def generate_cell_numbers(self) -> None: index = node.generate_cell_numbers(index) self.ncells = index - def cps(self) -> ControlPoints: + def cps(self) -> FloatArray: cps = np.zeros((self.ncps, self.dimension)) for node in self.catalogue.top_nodes(): assert node.cp_numbers is not None From 00041be42c2d1ec96c9d414561c07b4176091af0 Mon Sep 17 00:00:00 2001 From: Eivind Fonn Date: Sun, 4 Jan 2026 10:34:55 +0100 Subject: [PATCH 10/22] Typing: stl.py --- src/splipy/io/stl.py | 17 +++++++++-------- 1 file changed, 9 insertions(+), 8 deletions(-) diff --git a/src/splipy/io/stl.py b/src/splipy/io/stl.py index 500571b..d417ff3 100644 --- a/src/splipy/io/stl.py +++ b/src/splipy/io/stl.py @@ -148,17 +148,18 @@ def __enter__(self) -> Self: def write(self, obj: SplineModel | Surface | Volume, n: int | Sequence[int] | None = None) -> None: if isinstance(obj, SplineModel): if obj.pardim == 3: # volume model - for surface in obj.boundary(): - self.write_surface(surface.obj, n) + for node in obj.boundary(): + assert isinstance(node.obj, Surface) + self.write_surface(node.obj, n) elif obj.pardim == 2: # surface model - for surface in obj.objects(): - assert isinstance(surface, Surface) - self.write_surface(surface, n) + for comp in obj.objects(): + assert isinstance(comp, Surface) + self.write_surface(comp, n) elif isinstance(obj, Volume): - for surface in obj.faces(): - if surface is not None: # happens with periodic volumes - self.write_surface(surface, n) + for surf in obj.faces(): + if surf is not None: # happens with periodic volumes + self.write_surface(surf, n) elif isinstance(obj, Surface): self.write_surface(obj, n) From cd94127e14482ee808517ccf1fb0bbd7d5c0c674 Mon Sep 17 00:00:00 2001 From: Eivind Fonn Date: Sun, 4 Jan 2026 10:37:34 +0100 Subject: [PATCH 11/22] Remove bisect --- src/splipy/splinemodel.py | 4 +- src/splipy/utils/bisect.py | 101 ------------------------------------- 2 files changed, 2 insertions(+), 103 deletions(-) delete mode 100644 src/splipy/utils/bisect.py diff --git a/src/splipy/splinemodel.py b/src/splipy/splinemodel.py index b91014c..fa39404 100644 --- a/src/splipy/splinemodel.py +++ b/src/splipy/splinemodel.py @@ -1,5 +1,6 @@ from __future__ import annotations +import bisect from collections import Counter, OrderedDict, namedtuple from collections.abc import Callable, Iterator from dataclasses import dataclass @@ -17,7 +18,6 @@ from . import state from .splineobject import SplineObject from .utils import ( - bisect, check_section, is_right_hand, section_from_index, @@ -76,7 +76,7 @@ class VertexDict[T](MutableMapping[FloatArray, T]): _keys: list[FloatArray | None] _values: list[T | None] - lut: dict[tuple[int, ...], list[tuple[int, float]]] + lut: dict[tuple[int, ...], list[tuple[int, Scalar]]] def __init__(self, rtol: float = 1e-5, atol: float = 1e-8) -> None: self.rtol = rtol diff --git a/src/splipy/utils/bisect.py b/src/splipy/utils/bisect.py deleted file mode 100644 index ce44e98..0000000 --- a/src/splipy/utils/bisect.py +++ /dev/null @@ -1,101 +0,0 @@ -"""Bisection algorithms.""" - -from __future__ import annotations - - -def insort_right(a, x, lo=0, hi=None, *, key=None): - """Insert item x in list a, and keep it sorted assuming a is sorted. - - If x is already in a, insert it to the right of the rightmost x. - - Optional args lo (default 0) and hi (default len(a)) bound the - slice of a to be searched. - """ - lo = bisect_right(a, x, lo, hi) if key is None else bisect_right(a, key(x), lo, hi, key=key) - a.insert(lo, x) - - -def bisect_right(a, x, lo=0, hi=None, *, key=None): - """Return the index where to insert item x in list a, assuming a is sorted. - - The return value i is such that all e in a[:i] have e <= x, and all e in - a[i:] have e > x. So if x already appears in the list, a.insert(i, x) will - insert just after the rightmost x already there. - - Optional args lo (default 0) and hi (default len(a)) bound the - slice of a to be searched. - """ - - if lo < 0: - raise ValueError("lo must be non-negative") - if hi is None: - hi = len(a) - # Note, the comparison uses "<" to match the - # __lt__() logic in list.sort() and in heapq. - if key is None: - while lo < hi: - mid = (lo + hi) // 2 - if x < a[mid]: - hi = mid - else: - lo = mid + 1 - else: - while lo < hi: - mid = (lo + hi) // 2 - if x < key(a[mid]): - hi = mid - else: - lo = mid + 1 - return lo - - -def insort_left(a, x, lo=0, hi=None, *, key=None): - """Insert item x in list a, and keep it sorted assuming a is sorted. - - If x is already in a, insert it to the left of the leftmost x. - - Optional args lo (default 0) and hi (default len(a)) bound the - slice of a to be searched. - """ - - lo = bisect_left(a, x, lo, hi) if key is None else bisect_left(a, key(x), lo, hi, key=key) - a.insert(lo, x) - - -def bisect_left(a, x, lo=0, hi=None, *, key=None): - """Return the index where to insert item x in list a, assuming a is sorted. - - The return value i is such that all e in a[:i] have e < x, and all e in - a[i:] have e >= x. So if x already appears in the list, a.insert(i, x) will - insert just before the leftmost x already there. - - Optional args lo (default 0) and hi (default len(a)) bound the - slice of a to be searched. - """ - - if lo < 0: - raise ValueError("lo must be non-negative") - if hi is None: - hi = len(a) - # Note, the comparison uses "<" to match the - # __lt__() logic in list.sort() and in heapq. - if key is None: - while lo < hi: - mid = (lo + hi) // 2 - if a[mid] < x: - lo = mid + 1 - else: - hi = mid - else: - while lo < hi: - mid = (lo + hi) // 2 - if key(a[mid]) < x: - lo = mid + 1 - else: - hi = mid - return lo - - -# Create aliases -bisect = bisect_right -insort = insort_right From 951e608a59c4d9b7c707a637f9147a21ca05ddc9 Mon Sep 17 00:00:00 2001 From: Eivind Fonn Date: Sun, 4 Jan 2026 10:41:27 +0100 Subject: [PATCH 12/22] Typing: io/master.py --- src/splipy/io/grdecl.py | 2 +- src/splipy/io/master.py | 21 +++++++++++++++++---- src/splipy/io/stl.py | 9 +++++++-- src/splipy/io/threedm.py | 5 +++-- 4 files changed, 28 insertions(+), 9 deletions(-) diff --git a/src/splipy/io/grdecl.py b/src/splipy/io/grdecl.py index 0b53c4a..692e55a 100644 --- a/src/splipy/io/grdecl.py +++ b/src/splipy/io/grdecl.py @@ -168,7 +168,7 @@ def get_discontinuous_z(self) -> FloatArray: return self.Xz -class GRDECL(MasterIO): +class GRDECL: filename: str attribute: dict[str, FloatArray | IntArray] fstream: TextIO diff --git a/src/splipy/io/master.py b/src/splipy/io/master.py index 7a6d592..bd84219 100644 --- a/src/splipy/io/master.py +++ b/src/splipy/io/master.py @@ -1,18 +1,31 @@ from __future__ import annotations +from types import TracebackType +from typing import Self, Sequence + +from splipy.splinemodel import SplineModel +from splipy.splineobject import SplineObject class MasterIO: - def __init__(self, filename): + def __init__(self, filename: str) -> None: """Create an IO object attached to a file. :param str filename: The file to read from or write to """ raise NotImplementedError() - def __enter__(self): + def __enter__(self) -> Self: raise NotImplementedError() - def write(self, obj): + def __exit__( + self, + exc_type: type[BaseException], + exc_value: BaseException, + traceback: TracebackType, + ) -> None: + pass + + def write(self, obj: SplineObject | Sequence[SplineObject] | SplineModel) -> None: """Write one or more objects to the file. :param obj: The object(s) to write @@ -20,7 +33,7 @@ def write(self, obj): """ raise NotImplementedError() - def read(self): + def read(self) -> list[SplineObject]: """Reads all the objects from the file. :return: Objects diff --git a/src/splipy/io/stl.py b/src/splipy/io/stl.py index d417ff3..766eea4 100644 --- a/src/splipy/io/stl.py +++ b/src/splipy/io/stl.py @@ -9,6 +9,7 @@ import numpy as np from splipy.splinemodel import SplineModel +from splipy.splineobject import SplineObject from splipy.surface import Surface from splipy.utils import ensure_listlike from splipy.volume import Volume @@ -145,7 +146,7 @@ def __enter__(self) -> Self: self.writer = ASCII_STL_Writer(Path(self.filename).open("w")) return self - def write(self, obj: SplineModel | Surface | Volume, n: int | Sequence[int] | None = None) -> None: + def write(self, obj: SplineObject | Sequence[SplineObject] | SplineModel, n: int | Sequence[int] | None = None) -> None: if isinstance(obj, SplineModel): if obj.pardim == 3: # volume model for node in obj.boundary(): @@ -164,9 +165,13 @@ def write(self, obj: SplineModel | Surface | Volume, n: int | Sequence[int] | No elif isinstance(obj, Surface): self.write_surface(obj, n) - else: + elif isinstance(obj, SplineObject): raise ValueError("Unsopported object for STL format") + else: + for sub in obj: + self.write(sub) + def write_surface(self, surface: Surface, n: int | Sequence[int] | None = None) -> None: # choose evaluation points as one of three cases: # 1. specified with input diff --git a/src/splipy/io/threedm.py b/src/splipy/io/threedm.py index b6c9cee..0116ff0 100644 --- a/src/splipy/io/threedm.py +++ b/src/splipy/io/threedm.py @@ -1,6 +1,6 @@ from __future__ import annotations -from typing import TYPE_CHECKING, Any, Protocol, Self, cast +from typing import TYPE_CHECKING, Any, Protocol, Self, Sequence, cast import numpy as np from rhino3dm import ( @@ -26,6 +26,7 @@ from rhino3dm import Surface as threedmSurface # name conflict with splipy from splipy import BSplineBasis, Curve, Surface, curve_factory +from splipy.splinemodel import SplineModel from .master import MasterIO @@ -63,7 +64,7 @@ def __init__(self, filename: str) -> None: def __enter__(self) -> Self: return self - def write(self, _: SplineObject) -> None: + def write(self, _: SplineObject | Sequence[SplineObject] | SplineModel) -> None: raise OSError("Writing to 3DM not supported") def read(self) -> list[SplineObject]: From e7e6b60b0d5eb03411ce05a862c117c0b84ec1e5 Mon Sep 17 00:00:00 2001 From: Eivind Fonn Date: Sun, 4 Jan 2026 10:43:20 +0100 Subject: [PATCH 13/22] Typing: state.py --- src/splipy/state.py | 12 +++++++++++- 1 file changed, 11 insertions(+), 1 deletion(-) diff --git a/src/splipy/state.py b/src/splipy/state.py index 56bb6c0..006fb2d 100644 --- a/src/splipy/state.py +++ b/src/splipy/state.py @@ -4,6 +4,7 @@ import sys from contextlib import contextmanager +from typing import Iterator, TypedDict, Unpack states = [ "controlpoint_relative_tolerance", @@ -35,8 +36,17 @@ """Since splipy insists on finite parametric domains, we define 'unbounded' here""" +class StateKwargs(TypedDict, total=False): + controlpoint_absolute_tolerance: float + controlpoint_relative_tolerance: float + parametric_absolute_tolerance: float + parametric_relative_tolerance: float + knot_tolerance: float + unlimited: float + + @contextmanager -def state(**kwargs): +def state(**kwargs: Unpack[StateKwargs]) -> Iterator[None]: """A context manager for running code in a modified state. This takes an arbitrary number of keyword arguments, which correspond to From 439c57c27bab4b821c12fee517187224242896a3 Mon Sep 17 00:00:00 2001 From: Eivind Fonn Date: Sun, 4 Jan 2026 11:02:57 +0100 Subject: [PATCH 14/22] Typing: utils/__init__.py --- src/splipy/basis.py | 4 +- src/splipy/splineobject.py | 2 +- src/splipy/surface.py | 5 +- src/splipy/utils/__init__.py | 110 ++++++++++++++++++----------------- 4 files changed, 60 insertions(+), 61 deletions(-) diff --git a/src/splipy/basis.py b/src/splipy/basis.py index a5b30e9..e1578be 100644 --- a/src/splipy/basis.py +++ b/src/splipy/basis.py @@ -12,7 +12,7 @@ from scipy.sparse import csr_matrix from . import state -from .utils import ensure_listlike_old +from .utils import ensure_listlike if TYPE_CHECKING: from splipy.typing import Knots @@ -229,7 +229,7 @@ def evaluate_old(self, t, d=0, from_right=True, sparse=False): # type: ignore[n :rtype: numpy.array """ # for single-value input, wrap it into a list so it don't crash on the loop below - t = ensure_listlike_old(t) + t = ensure_listlike(t) self.snap(t) p = self.order # knot vector order diff --git a/src/splipy/splineobject.py b/src/splipy/splineobject.py index 4187964..eaa3125 100644 --- a/src/splipy/splineobject.py +++ b/src/splipy/splineobject.py @@ -1139,7 +1139,7 @@ def scale(self, arg: Point | Scalar, *args: Scalar) -> Self: return self - def rotate(self, theta: float, normal: ArrayLike = (0, 0, 1)) -> Self: + def rotate(self, theta: Scalar, normal: ArrayLike = (0, 0, 1)) -> Self: """Rotate the object around an axis. :param float theta: Angle to rotate about, measured in radians diff --git a/src/splipy/surface.py b/src/splipy/surface.py index 73c4b36..6483a88 100644 --- a/src/splipy/surface.py +++ b/src/splipy/surface.py @@ -144,11 +144,8 @@ def derivative( u = np.atleast_1d(np.asarray(u, dtype=np.float64)) v = np.atleast_1d(np.asarray(v, dtype=np.float64)) - # u = ensure_listlike_old(u) - # v = ensure_listlike_old(v) + result: FloatArray = np.zeros((len(u), len(v), self.dimension)) - # dNus = [self.bases[0].evaluate(u, d, above) for d in range(derivs[0]+1)] - # dNvs = [self.bases[1].evaluate(v, d, above) for d in range(derivs[1]+1)] dNus = [self.bases[0].evaluate(u, d, above[0]) for d in range(np.sum(derivs) + 1)] dNvs = [self.bases[1].evaluate(v, d, above[1]) for d in range(np.sum(derivs) + 1)] diff --git a/src/splipy/utils/__init__.py b/src/splipy/utils/__init__.py index 99bf178..7eb49d7 100644 --- a/src/splipy/utils/__init__.py +++ b/src/splipy/utils/__init__.py @@ -3,10 +3,12 @@ from collections.abc import Iterator, Sequence, Sized from itertools import combinations, product, repeat from math import atan2, sqrt -from typing import TYPE_CHECKING, SupportsFloat, TypeVar, Unpack +from typing import TYPE_CHECKING, Any, Iterable, Literal, SupportsFloat, TypeVar, Unpack, cast import numpy as np +from splipy.typing import Int, Knots, Scalar + if TYPE_CHECKING: from splipy.splineobject import SplineObject from splipy.typing import ( @@ -49,16 +51,13 @@ def iter() -> Iterator[float]: return np.fromiter(iter(), dtype=np.float64, count=count) -def with_repeated_knots( - knots: FloatArray | IntArray, - reps: int = 1, -) -> FloatArray: +def with_repeated_knots(knots: FloatArray, reps: int = 1) -> FloatArray: if reps == 1: return knots return np.pad(knots, pad_width=reps - 1, mode="edge") -def is_right_hand(patch, tol=1e-3): +def is_right_hand(patch: SplineObject, tol: float = 1e-3) -> bool: param = tuple((a + b) / 2 for a, b in zip(patch.start(), patch.end())) if patch.dimension == patch.pardim == 3: @@ -72,7 +71,7 @@ def is_right_hand(patch, tol=1e-3): dw = dw / np.linalg.norm(dw) # Compare cross product - return np.dot(dw, np.cross(du, dv)) >= tol + return np.dot(dw, np.cross(du, dv)) >= tol # type: ignore[no-any-return] if patch.dimension == patch.pardim == 2: du = patch.derivative(*param, d=(1, 0)) @@ -82,13 +81,14 @@ def is_right_hand(patch, tol=1e-3): du = du / np.linalg.norm(du) dv = dv / np.linalg.norm(dv) - return np.cross(du, dv) >= tol + return np.cross(du, dv) >= tol # type: ignore[return-value] raise ValueError("Right-handedness only defined for 2D or 3D patches in 2D or 3D space, respectively") -def rotation_matrix(theta, axis): - axis = axis / np.sqrt(np.dot(axis, axis)) +def rotation_matrix(theta: Scalar, axis: Point) -> FloatArray: + axis = np.asarray(axis, dtype=float) + axis /= np.linalg.norm(axis) a = np.cos(theta / 2) b, c, d = -axis * np.sin(theta / 2) return np.array( @@ -96,7 +96,8 @@ def rotation_matrix(theta, axis): [a * a + b * b - c * c - d * d, 2 * (b * c - a * d), 2 * (b * d + a * c)], [2 * (b * c + a * d), a * a + c * c - b * b - d * d, 2 * (c * d - a * b)], [2 * (b * d - a * c), 2 * (c * d + a * b), a * a + d * d - b * b - c * c], - ] + ], + dtype=float, ) @@ -110,16 +111,17 @@ def sections(src_dim: int, tgt_dim: int) -> Iterator[Section]: """ # Enumerate all combinations of fixed directions nfixed = src_dim - tgt_dim + pool: list[Literal[0, -1]] = [0, -1] for fixed in combinations(range(src_dim), r=nfixed): # Enumerate all {0,-1}^n over the fixed directions - for indices in product([0, -1], repeat=nfixed): - args = [None] * src_dim + for indices in product(pool, repeat=nfixed): + args: list[Literal[0, -1] | None] = [None] * src_dim for f, i in zip(fixed, indices[::-1]): args[f] = i yield tuple(args) -def section_from_index(src_dim, tgt_dim, i) -> Section: +def section_from_index(src_dim: int, tgt_dim: int, i: int) -> Section: """Return the i'th section from a source dimension to a target dimension. See :func:`splipy.Utils.sections` for more information. @@ -130,14 +132,14 @@ def section_from_index(src_dim, tgt_dim, i) -> Section: raise ValueError(f"No such section: {i} for dimensions {src_dim} and {tgt_dim}") -def section_to_index(section): +def section_to_index(section: Section) -> int: """Return the index corresponding to a section.""" src_dim = len(section) tgt_dim = sum(1 for s in section if s is None) for i, t in enumerate(sections(src_dim, tgt_dim)): if tuple(section) == tuple(t): return i - return None + assert False def check_section(*args: SectionElement, pardim: int, **kwargs: Unpack[SectionKwargs]) -> Section: @@ -152,6 +154,7 @@ def check_section(*args: SectionElement, pardim: int, **kwargs: Unpack[SectionKw while len(args_list) < pardim: args_list.append(None) for k in set(kwargs.keys()) & set("uvw"): + k = cast("Literal['u', 'v', 'w']", k) index = "uvw".index(k) args_list[index] = kwargs[k] return tuple(args_list) @@ -167,14 +170,7 @@ def check_direction(direction: Direction, pardim: int) -> int: raise ValueError("Invalid direction") -def ensure_flatlist(x): - """Flattens a multi-list x to a single index list.""" - if isinstance(x[0], Sized): - return x[0] - return x - - -def is_singleton(x): +def is_singleton(x: Any) -> bool: """Checks if x is list-like.""" return not isinstance(x, Sized) @@ -191,20 +187,6 @@ def ensure_listlike(x: T | Sequence[T], dups: int = 1) -> tuple[T, ...]: return (x,) * dups -# TODO(Eivind): Remove. -def ensure_listlike_old(x, dups=1): - """Wraps x in a list if it's not list-like.""" - try: - while len(x) < dups: - x = list(x) - x.append(x[-1]) - return x - except TypeError: - return [x] * dups - except IndexError: - return [] - - def normalize_points(*points: Point | Points) -> FloatArray: """Utility function for functions that accept a single sequence of points or multiple points as parameters. @@ -214,7 +196,7 @@ def normalize_points(*points: Point | Points) -> FloatArray: return np.asarray(points) -def rotate_local_x_axis(xaxis=(1, 0, 0), normal=(0, 0, 1)): +def rotate_local_x_axis(xaxis: Point = (1, 0, 0), normal: Point = (0, 0, 1)) -> Scalar: # rotate xaxis vector back to reference domain (r=1, around origin) theta = atan2(normal[1], normal[0]) phi = atan2(sqrt(normal[0] ** 2 + normal[1] ** 2), normal[2]) @@ -222,29 +204,42 @@ def rotate_local_x_axis(xaxis=(1, 0, 0), normal=(0, 0, 1)): R2 = rotation_matrix(-phi, (0, 1, 0)) if len(xaxis) != 3: # typically 2D geometries xaxis = [xaxis[0], xaxis[1], 0] - xaxis = np.array([xaxis]) - xaxis = xaxis.dot(R1).dot(R2) + xaxis_array = np.array([xaxis]) + xaxis_array = xaxis_array.dot(R1).dot(R2) # if xaxis is orthogonal to normal, then xaxis[2]==0 now. If not then # treating it as such is the closest projection, which makes perfect sense - return atan2(xaxis[0, 1], xaxis[0, 0]) + return atan2(xaxis_array[0, 1], xaxis_array[0, 0]) -def flip_and_move_plane_geometry[T: SplineObject](obj: T, center=(0, 0, 0), normal=(0, 0, 1)) -> T: +def flip_and_move_plane_geometry[T: SplineObject]( + obj: T, + center: Point = (0, 0, 0), + normal: Point = (0, 0, 1), +) -> T: """re-orients a planar geometry by moving it to a different location and tilting it""" # don't touch it if not needed. translate or scale operations may force # object into 3D space - if not np.allclose(normal, np.array([0, 0, 1])): - theta = atan2(normal[1], normal[0]) - phi = atan2(sqrt(normal[0] ** 2 + normal[1] ** 2), normal[2]) + normal_arr = np.asarray(normal) + if not np.allclose(normal_arr, np.array([0, 0, 1])): + theta = atan2(normal_arr[1], normal_arr[0]) + phi = atan2(sqrt(normal_arr[0] ** 2 + normal_arr[1] ** 2), normal_arr[2]) obj.rotate(phi, (0, 1, 0)) obj.rotate(theta, (0, 0, 1)) - if not np.allclose(center, 0): - obj.translate(center) + + center_arr = np.asarray(center) + if not np.allclose(center_arr, 0): + obj.translate(center_arr) + return obj -def reshape(cps, newshape, order="C", ncomps=None): +def reshape( + cps: FloatArray, + newshape: tuple[int, ...], + order: Literal["C", "F"] = "C", + ncomps: Int | None = None, +) -> FloatArray: """Like numpy's reshape, but preserves control points of several dimensions that are stored contiguously. @@ -263,10 +258,12 @@ def reshape(cps, newshape, order="C", ncomps=None): except AttributeError: ncomps = len(cps) // npts + shape: Sequence[Int] if order == "C": shape = list(newshape) + [ncomps] elif order == "F": shape = list(newshape[::-1]) + [ncomps] + cps = np.reshape(cps, shape) if order == "F": spec = list(range(len(newshape)))[::-1] + [len(newshape)] @@ -274,7 +271,7 @@ def reshape(cps, newshape, order="C", ncomps=None): return cps -def uniquify(iterator): +def uniquify[T](iterator: Iterable[T]) -> Iterator[T]: """Iterates over all elements in `iterator`, removing duplicates.""" seen = set() for i in iterator: @@ -284,7 +281,14 @@ def uniquify(iterator): yield i -def raise_order_1D(n, k, T, P, m, periodic): +def raise_order_1D( + n: int, + k: int, + T: FloatArray, + P: FloatArray, + m: int, + periodic: int +) -> FloatArray: """Implementation of method in "Efficient Degree Elevation and Knot Insertion for B-spline Curves using Derivatives" by Qi-Xing Huang a Shi-Min Hu, Ralph R Martin. Only the case of open knot vector is fully implemented. @@ -307,7 +311,7 @@ def raise_order_1D(n, k, T, P, m, periodic): # Find multiplicity of the knot vector T b = BSplineBasis(k, T) - z = [k - 1 - b.continuity(t0) for t0 in b.knot_spans()] + z = [k - 1 - b.knot_continuity(t0) for t0 in b.knot_spans()] # Step 1: Find Pt_i^j Pt = np.zeros((d, n + 1, k)) @@ -373,9 +377,7 @@ def raise_order_1D(n, k, T, P, m, periodic): "section_to_index", "check_section", "check_direction", - "ensure_flatlist", "is_singleton", - "ensure_listlike_old", "rotate_local_x_axis", "flip_and_move_plane_geometry", "reshape", From 11a28faeb2d339094fc7e77509e57f296154afbb Mon Sep 17 00:00:00 2001 From: Eivind Fonn Date: Sun, 4 Jan 2026 16:04:17 +0100 Subject: [PATCH 15/22] Lint and format --- src/splipy/basis.py | 1 - src/splipy/curve.py | 4 ++-- src/splipy/io/grdecl.py | 1 - src/splipy/io/master.py | 12 ++++++++---- src/splipy/io/ofoam.py | 8 +++++--- src/splipy/io/stl.py | 4 +++- src/splipy/io/threedm.py | 6 +++--- src/splipy/splinemodel.py | 26 +++++++++++++++----------- src/splipy/state.py | 5 ++++- src/splipy/utils/__init__.py | 18 +++++------------- src/splipy/utils/curve.py | 15 ++++++++------- 11 files changed, 53 insertions(+), 47 deletions(-) diff --git a/src/splipy/basis.py b/src/splipy/basis.py index e1578be..24a61bc 100644 --- a/src/splipy/basis.py +++ b/src/splipy/basis.py @@ -138,7 +138,6 @@ def greville(self, index: Int) -> float: ... @overload def greville(self) -> FloatArray: ... - def greville(self, index: Int | None = None) -> float | FloatArray: """Fetch greville points, also known as knot averages: diff --git a/src/splipy/curve.py b/src/splipy/curve.py index f40f7f5..d12ee80 100644 --- a/src/splipy/curve.py +++ b/src/splipy/curve.py @@ -6,8 +6,6 @@ import numpy as np import scipy.sparse.linalg as splinalg -from splipy.typing import Point - from . import state from .basis import BSplineBasis from .splineobject import SplineObject @@ -16,6 +14,8 @@ if TYPE_CHECKING: from collections.abc import Callable, Sequence + from splipy.typing import Point + from .typing import ArrayLike, Direction, FloatArray, Scalar __all__ = ["Curve"] diff --git a/src/splipy/io/grdecl.py b/src/splipy/io/grdecl.py index 692e55a..c6afb3c 100644 --- a/src/splipy/io/grdecl.py +++ b/src/splipy/io/grdecl.py @@ -19,7 +19,6 @@ from splipy.volume import Volume from .g2 import G2 -from .master import MasterIO if TYPE_CHECKING: from collections.abc import Sequence diff --git a/src/splipy/io/master.py b/src/splipy/io/master.py index bd84219..3a8d018 100644 --- a/src/splipy/io/master.py +++ b/src/splipy/io/master.py @@ -1,9 +1,13 @@ from __future__ import annotations -from types import TracebackType -from typing import Self, Sequence -from splipy.splinemodel import SplineModel -from splipy.splineobject import SplineObject +from typing import TYPE_CHECKING, Self + +if TYPE_CHECKING: + from collections.abc import Sequence + from types import TracebackType + + from splipy.splinemodel import SplineModel + from splipy.splineobject import SplineObject class MasterIO: diff --git a/src/splipy/io/ofoam.py b/src/splipy/io/ofoam.py index 5f83fb2..45f8dcc 100644 --- a/src/splipy/io/ofoam.py +++ b/src/splipy/io/ofoam.py @@ -3,13 +3,15 @@ from itertools import groupby from operator import itemgetter from pathlib import Path -from types import TracebackType -from typing import Self +from typing import TYPE_CHECKING, Self import numpy as np from splipy.splinemodel import FaceArray, SplineModel +if TYPE_CHECKING: + from types import TracebackType + class OpenFOAM: target: Path @@ -58,7 +60,7 @@ def write(self, model: SplineModel) -> None: model.generate_cp_numbers() model.generate_cell_numbers() faces = model.faces() - ninternal = sum(faces.name == None) + ninternal = sum(faces.name is None) note = ( f"nPoints: {model.ncps} nCells: {model.ncells} nFaces: {len(faces)} nInternalFaces: {ninternal}" ) diff --git a/src/splipy/io/stl.py b/src/splipy/io/stl.py index 766eea4..09df540 100644 --- a/src/splipy/io/stl.py +++ b/src/splipy/io/stl.py @@ -146,7 +146,9 @@ def __enter__(self) -> Self: self.writer = ASCII_STL_Writer(Path(self.filename).open("w")) return self - def write(self, obj: SplineObject | Sequence[SplineObject] | SplineModel, n: int | Sequence[int] | None = None) -> None: + def write( + self, obj: SplineObject | Sequence[SplineObject] | SplineModel, n: int | Sequence[int] | None = None + ) -> None: if isinstance(obj, SplineModel): if obj.pardim == 3: # volume model for node in obj.boundary(): diff --git a/src/splipy/io/threedm.py b/src/splipy/io/threedm.py index 0116ff0..570dd26 100644 --- a/src/splipy/io/threedm.py +++ b/src/splipy/io/threedm.py @@ -1,6 +1,6 @@ from __future__ import annotations -from typing import TYPE_CHECKING, Any, Protocol, Self, Sequence, cast +from typing import TYPE_CHECKING, Any, Protocol, Self, cast import numpy as np from rhino3dm import ( @@ -26,14 +26,14 @@ from rhino3dm import Surface as threedmSurface # name conflict with splipy from splipy import BSplineBasis, Curve, Surface, curve_factory -from splipy.splinemodel import SplineModel from .master import MasterIO if TYPE_CHECKING: - from collections.abc import Iterable + from collections.abc import Iterable, Sequence from types import TracebackType + from splipy.splinemodel import SplineModel from splipy.splineobject import SplineObject from splipy.typing import Point diff --git a/src/splipy/splinemodel.py b/src/splipy/splinemodel.py index fa39404..9608b63 100644 --- a/src/splipy/splinemodel.py +++ b/src/splipy/splinemodel.py @@ -1,19 +1,26 @@ from __future__ import annotations import bisect -from collections import Counter, OrderedDict, namedtuple -from collections.abc import Callable, Iterator +from collections import Counter, OrderedDict +from collections.abc import Callable, Iterator, Sequence from dataclasses import dataclass from itertools import chain, islice, permutations, product from operator import itemgetter from pathlib import Path -from re import S -from typing import Any, Sequence, Unpack, cast, overload +from typing import Any, Unpack, cast import numpy as np import numpy.typing as npt -from splipy.typing import ControlPoints, FloatArray, IntArray, Scalar, Section, Int, SectionElement, SectionKwargs +from splipy.typing import ( + FloatArray, + Int, + IntArray, + Scalar, + Section, + SectionElement, + SectionKwargs, +) from . import state from .splineobject import SplineObject @@ -181,6 +188,7 @@ class OrientationError(RuntimeError): :class:`splipy.SplineModel.Orientation` indicating an inability to match two objects. """ + pass @@ -188,6 +196,7 @@ class TwinError(RuntimeError): """A `TwinError` is raised when two objects with identical interfaces are added, but different interiors. """ + pass @@ -426,12 +435,7 @@ class TopologicalNode: cell_numbers: IntArray | None cp_numbers: IntArray | None - def __init__( - self, - obj: SplineObject, - lower_nodes: list[tuple[TopologicalNode, ...]], - index: int - ) -> None: + def __init__(self, obj: SplineObject, lower_nodes: list[tuple[TopologicalNode, ...]], index: int) -> None: """Initialize a `TopologicalNode` object associated with the given `SplineObject` and lower order nodes. diff --git a/src/splipy/state.py b/src/splipy/state.py index 006fb2d..8078405 100644 --- a/src/splipy/state.py +++ b/src/splipy/state.py @@ -4,7 +4,10 @@ import sys from contextlib import contextmanager -from typing import Iterator, TypedDict, Unpack +from typing import TYPE_CHECKING, TypedDict, Unpack + +if TYPE_CHECKING: + from collections.abc import Iterator states = [ "controlpoint_relative_tolerance", diff --git a/src/splipy/utils/__init__.py b/src/splipy/utils/__init__.py index 7eb49d7..a96e5f5 100644 --- a/src/splipy/utils/__init__.py +++ b/src/splipy/utils/__init__.py @@ -1,22 +1,21 @@ from __future__ import annotations -from collections.abc import Iterator, Sequence, Sized +from collections.abc import Iterable, Iterator, Sequence, Sized from itertools import combinations, product, repeat from math import atan2, sqrt -from typing import TYPE_CHECKING, Any, Iterable, Literal, SupportsFloat, TypeVar, Unpack, cast +from typing import TYPE_CHECKING, Any, Literal, SupportsFloat, TypeVar, Unpack, cast import numpy as np -from splipy.typing import Int, Knots, Scalar - if TYPE_CHECKING: from splipy.splineobject import SplineObject from splipy.typing import ( Direction, FloatArray, - IntArray, + Int, Point, Points, + Scalar, Section, SectionElement, SectionKwargs, @@ -281,14 +280,7 @@ def uniquify[T](iterator: Iterable[T]) -> Iterator[T]: yield i -def raise_order_1D( - n: int, - k: int, - T: FloatArray, - P: FloatArray, - m: int, - periodic: int -) -> FloatArray: +def raise_order_1D(n: int, k: int, T: FloatArray, P: FloatArray, m: int, periodic: int) -> FloatArray: """Implementation of method in "Efficient Degree Elevation and Knot Insertion for B-spline Curves using Derivatives" by Qi-Xing Huang a Shi-Min Hu, Ralph R Martin. Only the case of open knot vector is fully implemented. diff --git a/src/splipy/utils/curve.py b/src/splipy/utils/curve.py index 3ae6612..3cb0f83 100644 --- a/src/splipy/utils/curve.py +++ b/src/splipy/utils/curve.py @@ -1,14 +1,15 @@ from __future__ import annotations -from itertools import repeat - -from splipy.curve import Curve -from splipy.typing import FloatArray, Points - __doc__ = "Implementation of various curve utilities" +from typing import TYPE_CHECKING + import numpy as np +if TYPE_CHECKING: + from splipy.curve import Curve + from splipy.typing import FloatArray, Points + def curve_length_parametrization(pts: Points, normalize: bool = False, reps: int = 1) -> FloatArray: """Calculate knots corresponding to a curvelength parametrization of a set of @@ -26,10 +27,10 @@ def curve_length_parametrization(pts: Points, normalize: bool = False, reps: int distances = np.linalg.norm(points[1:, ...] - points[:-1, ...], axis=1) distances = np.cumsum(distances) - knots[reps:reps-1+npts] = distances + knots[reps : reps - 1 + npts] = distances if reps > 1: - knots[-reps+1:] = knots[-reps] + knots[-reps + 1 :] = knots[-reps] if normalize: knots /= knots[-1] From 592cf9aa8f4fb48dfee5d2c581b26017bbffabca Mon Sep 17 00:00:00 2001 From: Eivind Fonn Date: Sun, 4 Jan 2026 16:14:18 +0100 Subject: [PATCH 16/22] Add mypy to test pipeline --- Makefile | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Makefile b/Makefile index 1317b7f..1b3c780 100644 --- a/Makefile +++ b/Makefile @@ -53,7 +53,7 @@ examples: uv run python examples/write.py .PHONY: test -test: pytest ruff examples +test: pytest mypy ruff examples # Build targets (used from CI) From 766726210c6bd0ea71c64b7167f25f6382edffa5 Mon Sep 17 00:00:00 2001 From: Eivind Fonn Date: Sun, 4 Jan 2026 16:18:05 +0100 Subject: [PATCH 17/22] ofoam.py: fix lint-induced regression in --- src/splipy/io/ofoam.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/splipy/io/ofoam.py b/src/splipy/io/ofoam.py index 45f8dcc..a7f0291 100644 --- a/src/splipy/io/ofoam.py +++ b/src/splipy/io/ofoam.py @@ -60,7 +60,7 @@ def write(self, model: SplineModel) -> None: model.generate_cp_numbers() model.generate_cell_numbers() faces = model.faces() - ninternal = sum(faces.name is None) + ninternal = sum(faces.name == None) # noqa: E711 note = ( f"nPoints: {model.ncps} nCells: {model.ncells} nFaces: {len(faces)} nInternalFaces: {ninternal}" ) From 7ceb75315a2c2af57401409c8c63358791563ac0 Mon Sep 17 00:00:00 2001 From: Eivind Fonn Date: Sun, 4 Jan 2026 16:28:14 +0100 Subject: [PATCH 18/22] Ignore nutils, rhino when not installed This adds some mypy ignore rules to make type checking work even when certain libraries are not installed. --- src/splipy/io/threedm.py | 5 ++++- src/splipy/surface_factory.py | 4 ++-- src/splipy/utils/nutils.py | 2 +- 3 files changed, 7 insertions(+), 4 deletions(-) diff --git a/src/splipy/io/threedm.py b/src/splipy/io/threedm.py index 570dd26..c7fb1ba 100644 --- a/src/splipy/io/threedm.py +++ b/src/splipy/io/threedm.py @@ -1,9 +1,12 @@ +# Make type-checking work in cases where rhino3dm is not installed +# mypy: disable-error-code="no-any-unimported" + from __future__ import annotations from typing import TYPE_CHECKING, Any, Protocol, Self, cast import numpy as np -from rhino3dm import ( +from rhino3dm import ( # type: ignore[import-not-found,unused-ignore] Arc, BezierCurve, Brep, diff --git a/src/splipy/surface_factory.py b/src/splipy/surface_factory.py index c96c638..f5464c3 100644 --- a/src/splipy/surface_factory.py +++ b/src/splipy/surface_factory.py @@ -376,8 +376,8 @@ def coons_patch(bottom: Curve, right: Curve, top: Curve, left: Curve) -> Surface def poisson_patch(bottom: Curve, right: Curve, top: Curve, left: Curve) -> Surface: - from nutils import function as fn # type: ignore[import-untyped] - from nutils import mesh + from nutils import function as fn # type: ignore[import-untyped,import-not-found,unused-ignore] + from nutils import mesh # type: ignore[import-untyped,import-not-found,unused-ignore] # error test input if left.rational or right.rational or top.rational or bottom.rational: diff --git a/src/splipy/utils/nutils.py b/src/splipy/utils/nutils.py index 18a9091..f2d2732 100644 --- a/src/splipy/utils/nutils.py +++ b/src/splipy/utils/nutils.py @@ -44,7 +44,7 @@ def degree(spline: SplineObject) -> list[int]: def splipy_to_nutils(spline: SplineObject) -> Any: """Returns nutils domain and geometry object for spline mapping given by the argument""" - from nutils import function, mesh # type: ignore[import-untyped] + from nutils import function, mesh # type: ignore[import-untyped,import-not-found,unused-ignore] domain, geom = mesh.rectilinear(spline.knots()) cp = controlpoints(spline) From e87ea66a7cb233ab936f017df7cef237a028cc85 Mon Sep 17 00:00:00 2001 From: Eivind Fonn Date: Sun, 4 Jan 2026 17:32:54 +0100 Subject: [PATCH 19/22] Clean up some type hints --- src/splipy/basis.py | 18 +++++++++--------- src/splipy/curve.py | 16 ++++++++-------- src/splipy/curve_factory.py | 16 ++++++++-------- src/splipy/splineobject.py | 26 +++++++++++++------------- src/splipy/surface.py | 10 ++++++---- src/splipy/surface_factory.py | 6 +++--- src/splipy/typing.py | 1 + src/splipy/volume.py | 4 +++- 8 files changed, 51 insertions(+), 46 deletions(-) diff --git a/src/splipy/basis.py b/src/splipy/basis.py index 24a61bc..62e0b0e 100644 --- a/src/splipy/basis.py +++ b/src/splipy/basis.py @@ -15,7 +15,7 @@ from .utils import ensure_listlike if TYPE_CHECKING: - from splipy.typing import Knots + from splipy.typing import Knots, Params from .typing import FloatArray, Int, Scalar @@ -90,16 +90,16 @@ def num_functions(self) -> int: .. warning:: This is different from :func:`splipy.BSplineBasis.__len__`.""" return len(self.knots) - self.order - (self.periodic + 1) - def start(self) -> float: + def start(self) -> Scalar: """Start point of parametric domain. For open knot vectors, this is the first knot. :return: Knot number *p*, where *p* is the spline order :rtype: float """ - return float(self.knots.flat[self.order - 1]) + return self.knots.flat[self.order - 1] - def end(self) -> float: + def end(self) -> Scalar: """End point of parametric domain. For open knot vectors, this is the last knot. @@ -107,7 +107,7 @@ def end(self) -> float: the number of knots :rtype: Float """ - return float(self.knots.flat[-self.order]) + return self.knots.flat[-self.order] def greville_all(self) -> FloatArray: """Fetch all greville points, also known as knot averages: @@ -153,7 +153,7 @@ def greville(self, index: Int | None = None) -> float | FloatArray: @overload def evaluate( self, - t: Knots | Scalar, + t: Params | Scalar, d: int = 0, from_right: bool = ..., ) -> npt.NDArray[np.double]: ... @@ -161,7 +161,7 @@ def evaluate( @overload def evaluate( self, - t: Knots | Scalar, + t: Params | Scalar, d: int = 0, from_right: bool = ..., sparse: Literal[False] = ..., @@ -170,7 +170,7 @@ def evaluate( @overload def evaluate( self, - t: Knots | Scalar, + t: Params | Scalar, d: int = 0, from_right: bool = ..., sparse: Literal[True] = ..., @@ -178,7 +178,7 @@ def evaluate( def evaluate( self, - t: Knots | Scalar, + t: Params | Scalar, d: int = 0, from_right: bool = True, sparse: bool = False, diff --git a/src/splipy/curve.py b/src/splipy/curve.py index d12ee80..6842541 100644 --- a/src/splipy/curve.py +++ b/src/splipy/curve.py @@ -14,7 +14,7 @@ if TYPE_CHECKING: from collections.abc import Callable, Sequence - from splipy.typing import Point + from splipy.typing import Params, Point from .typing import ArrayLike, Direction, FloatArray, Scalar @@ -48,7 +48,7 @@ def __init__( """ super().__init__([basis], controlpoints, rational, raw=raw) - def evaluate(self, *params: ArrayLike | Scalar, tensor: bool = True) -> FloatArray: + def evaluate(self, *params: Params | Scalar, tensor: bool = True) -> FloatArray: """Evaluate the object at given parametric values. This function returns an *n1* × *n2* × ... × *dim* array, where *ni* is @@ -88,7 +88,7 @@ def evaluate(self, *params: ArrayLike | Scalar, tensor: bool = True) -> FloatArr def derivative( self, - *params: ArrayLike | Scalar, + *params: Params | Scalar, d: int | Sequence[int] = 1, above: bool | Sequence[bool] = True, tensor: bool = True, @@ -147,7 +147,7 @@ def derivative( return result - def binormal(self, t: ArrayLike | Scalar, above: bool = True) -> FloatArray: + def binormal(self, t: Params | Scalar, above: bool = True) -> FloatArray: """Evaluate the normalized binormal of the curve at the given parametric value(s). This function returns an *n* × 3 array, where *n* is the number of @@ -195,7 +195,7 @@ def binormal(self, t: ArrayLike | Scalar, above: bool = True) -> FloatArray: return result / magnitude - def normal(self, t: ArrayLike | Scalar, above: bool = True) -> FloatArray: + def normal(self, t: Params | Scalar, above: bool = True) -> FloatArray: """Evaluate the normal of the curve at the given parametric value(s). This function returns an *n* × 3 array, where *n* is the number of @@ -220,7 +220,7 @@ def normal(self, t: ArrayLike | Scalar, above: bool = True) -> FloatArray: return np.cross(B, T) - def curvature(self, t: ArrayLike | Scalar, above: bool = True) -> FloatArray | float: + def curvature(self, t: Params | Scalar, above: bool = True) -> FloatArray | float: """Evaluate the curvaure at specified point(s). The curvature is defined as .. math:: \\frac{|\\boldsymbol{v}\\times \\boldsymbol{a}|}{|\\boldsymbol{v}|^3} @@ -244,7 +244,7 @@ def curvature(self, t: ArrayLike | Scalar, above: bool = True) -> FloatArray | f speed: FloatArray = np.linalg.norm(v, axis=-1) return magnitude / speed**3 - def torsion(self, t: ArrayLike | Scalar, above: bool = True) -> FloatArray | float: + def torsion(self, t: Params | Scalar, above: bool = True) -> FloatArray | float: """Evaluate the torsion for a 3D curve at specified point(s). The torsion is defined as .. math:: \\frac{(\\boldsymbol{v}\\times \\boldsymbol{a})\\cdot @@ -541,7 +541,7 @@ def arclength_circle(t): err_inf = max(np.max(np.sqrt(error)), err_inf) return (np.array(err2, dtype=np.float64), err_inf) - def get_antiderivative_curve(self, constant: ArrayLike | None = None) -> Curve: + def get_antiderivative_curve(self, constant: Point | None = None) -> Curve: """Compute the antiderivative (integral) of the curve. The antiderivative is computed by inverting the derivative operator on diff --git a/src/splipy/curve_factory.py b/src/splipy/curve_factory.py index e34a994..8f43192 100644 --- a/src/splipy/curve_factory.py +++ b/src/splipy/curve_factory.py @@ -28,7 +28,7 @@ if TYPE_CHECKING: from collections.abc import Callable - from splipy.typing import FloatArray, Knots, Point, Points + from splipy.typing import FloatArray, Params, Point, Points from .typing import Scalar @@ -94,14 +94,14 @@ def line(a: Point, b: Point, relative: bool = False) -> Curve: @overload -def polygon(*points: Point, t: Knots | None = None, relative: bool = False) -> Curve: ... +def polygon(*points: Point, t: Params | None = None, relative: bool = False) -> Curve: ... @overload -def polygon(points: Points, /, *, t: Knots | None = None, relative: bool = False) -> Curve: ... +def polygon(points: Points, /, *, t: Params | None = None, relative: bool = False) -> Curve: ... -def polygon(*points_in: Point | Points, t: Knots | None = None, relative: bool = False) -> Curve: +def polygon(*points_in: Point | Points, t: Params | None = None, relative: bool = False) -> Curve: """Create a linear interpolation between input points. :param [array-like] points: The points to interpolate @@ -349,7 +349,7 @@ def circle_segment( return flip_and_move_plane_geometry(result, center, normal) -def interpolate(x: Points, basis: BSplineBasis, t: Knots | None = None) -> Curve: +def interpolate(x: Points, basis: BSplineBasis, t: Params | None = None) -> Curve: """Perform general spline interpolation on a provided basis. :param matrix-like x: Matrix *X[i,j]* of interpolation points *xi* with @@ -374,7 +374,7 @@ def interpolate(x: Points, basis: BSplineBasis, t: Knots | None = None) -> Curve return Curve(basis, cp) -def least_square_fit(x: Points, basis: BSplineBasis, t: Knots) -> Curve: +def least_square_fit(x: Points, basis: BSplineBasis, t: Params) -> Curve: """Perform a least-square fit of a point cloud onto a spline basis :param matrix-like x: Matrix *X[i,j]* of interpolation points *xi* with @@ -398,7 +398,7 @@ def least_square_fit(x: Points, basis: BSplineBasis, t: Knots) -> Curve: def cubic_curve( x: Points, boundary: Boundary = Boundary.FREE, - t: Knots | None = None, + t: Params | None = None, tangents: Points | None = None, ) -> Curve: """Perform cubic spline interpolation on a provided basis. @@ -742,7 +742,7 @@ def move_along_tangent(t): return crv -def fit_points(x: Points, t: Knots | None = None, rtol: Scalar = 1e-4, atol: Scalar = 0.0) -> Curve: +def fit_points(x: Points, t: Params | None = None, rtol: Scalar = 1e-4, atol: Scalar = 0.0) -> Curve: """Computes an approximation for a list of points up to a specified tolerance. The method will iteratively refine parts where needed resulting in a non-uniform knot vector with as optimized knot locations as possible. The target curve is the diff --git a/src/splipy/splineobject.py b/src/splipy/splineobject.py index eaa3125..29d4e71 100644 --- a/src/splipy/splineobject.py +++ b/src/splipy/splineobject.py @@ -25,7 +25,7 @@ if TYPE_CHECKING: from splipy.typing import ControlPoints, Point - from .typing import ArrayLike, Direction, FloatArray, Scalar, SectionElement, SectionKwargs + from .typing import ArrayLike, Direction, FloatArray, Params, Scalar, SectionElement, SectionKwargs __all__ = ["SplineObject"] @@ -154,7 +154,7 @@ def _validate_domain_old(self, *params): # type: ignore[no-untyped-def] def evaluate( self, - *params: ArrayLike | Scalar, + *params: Params | Scalar, tensor: bool = True, ) -> FloatArray: """Evaluate the object at given parametric values. @@ -206,7 +206,7 @@ def evaluate( def derivative( self, - *params: ArrayLike | Scalar, + *params: Params | Scalar, d: int | Sequence[int] = 1, above: bool | Sequence[bool] = True, tensor: bool = True, @@ -365,11 +365,11 @@ def get_antiderivative_spline(self) -> list[SplineObject]: ... @overload def get_antiderivative_spline( - self, direction: Direction, constant: ArrayLike | None = None + self, direction: Direction, constant: Point | None = None ) -> SplineObject: ... def get_antiderivative_spline( - self, direction: Direction | None = None, constant: ArrayLike | None = None + self, direction: Direction | None = None, constant: Point | None = None ) -> SplineObject | list[SplineObject]: """Compute the antiderivative (integral) of the spline object in a given parametric direction. @@ -488,7 +488,7 @@ def get_antiderivative_spline( @overload def tangent( self, - *params: ArrayLike | Scalar, + *params: Params | Scalar, direction: Direction, above: bool | Sequence[bool] = True, tensor: bool = True, @@ -497,7 +497,7 @@ def tangent( @overload def tangent( self, - *params: ArrayLike | Scalar, + *params: Params | Scalar, direction: None = None, above: bool | Sequence[bool] = True, tensor: bool = True, @@ -505,7 +505,7 @@ def tangent( def tangent( self, - *params: ArrayLike | Scalar, + *params: Params | Scalar, direction: Direction | None = None, above: bool | Sequence[bool] = True, tensor: bool = True, @@ -783,9 +783,9 @@ def lower_order(self, diff: int, *args: int) -> Self: def start(self) -> tuple[float, ...]: ... @overload - def start(self, direction: Direction) -> float: ... + def start(self, direction: Direction) -> Scalar: ... - def start(self, direction: Direction | None = None) -> float | tuple[float, ...]: + def start(self, direction: Direction | None = None) -> Scalar | tuple[Scalar, ...]: """Return the start of the parametric domain. If `direction` is given, returns the start of that direction, as a @@ -806,7 +806,7 @@ def end(self) -> tuple[float, ...]: ... @overload def end(self, direction: Direction) -> float: ... - def end(self, direction: Direction | None = None) -> float | tuple[float, ...]: + def end(self, direction: Direction | None = None) -> Scalar | tuple[Scalar, ...]: """Return the end of the parametric domain. If `direction` is given, returns the end of that direction, as a float. @@ -919,7 +919,7 @@ def swap(self, dir1: Direction = 0, dir2: Direction = 1) -> Self: return self - def insert_knot(self, knot: Scalar | ArrayLike, direction: Direction = 0) -> Self: + def insert_knot(self, knot: Scalar | Params, direction: Direction = 0) -> Self: """Insert a new knot into the spline. :param int direction: The direction to insert in @@ -988,7 +988,7 @@ def refine(self, n: int, *ns: int, direction: Direction | None = None) -> Self: for n, d in zip(args, directions): knots = self.knots(direction=d) # excluding multiple knots - new_knots: list[FloatArray] = [] + new_knots: list[Scalar] = [] for k0, k1 in zip(knots[:-1], knots[1:]): new_knots.extend(np.linspace(k0, k1, n + 2)[1:-1]) self.insert_knot(new_knots, d) diff --git a/src/splipy/surface.py b/src/splipy/surface.py index 6483a88..460c3df 100644 --- a/src/splipy/surface.py +++ b/src/splipy/surface.py @@ -13,6 +13,8 @@ if TYPE_CHECKING: from collections.abc import Sequence + from splipy.typing import Params, Point + from .typing import ArrayLike, Direction, FloatArray, Scalar __all__ = ["Surface"] @@ -49,8 +51,8 @@ def __init__( def normal( self, - u: ArrayLike | Scalar, - v: ArrayLike | Scalar, + u: Params | Scalar, + v: Params | Scalar, above: bool | Sequence[bool] = True, tensor: bool = True, ) -> FloatArray: @@ -107,7 +109,7 @@ def normal( def derivative( self, - *params: ArrayLike | Scalar, + *params: Params | Scalar, d: int | Sequence[int] = 1, above: bool | Sequence[bool] = True, tensor: bool = True, @@ -368,7 +370,7 @@ def __repr__(self) -> str: result += str(self.controlpoints[i, j, :]) + "\n" return result - def get_antiderivative_surface(self, direction: Direction, constant: ArrayLike | None = None) -> Surface: + def get_antiderivative_surface(self, direction: Direction, constant: Point | None = None) -> Surface: """Compute the antiderivative (integral) of the surface in a given parametric direction. The antiderivative is computed by inverting the derivative operator on diff --git a/src/splipy/surface_factory.py b/src/splipy/surface_factory.py index f5464c3..c4c7909 100644 --- a/src/splipy/surface_factory.py +++ b/src/splipy/surface_factory.py @@ -21,7 +21,7 @@ if TYPE_CHECKING: from collections.abc import Callable, Sequence - from splipy.typing import FloatArray, Knots, Point, Scalar + from splipy.typing import FloatArray, Params, Point, Scalar __all__ = [ "square", @@ -823,7 +823,7 @@ def loft(*in_curves: Curve | Sequence[Curve]) -> Surface: def interpolate( x: FloatArray, bases: Sequence[BSplineBasis], - u: Sequence[Knots] | None = None, + u: Sequence[Params] | None = None, ) -> Surface: """Interpolate a surface on a set of regular gridded interpolation points `x`. @@ -853,7 +853,7 @@ def interpolate( return Surface(bases[0], bases[1], cp.transpose(1, 0, 2).reshape((np.prod(surf_shape), dim))) -def least_square_fit(x: FloatArray, bases: Sequence[BSplineBasis], u: Sequence[Knots]) -> Surface: +def least_square_fit(x: FloatArray, bases: Sequence[BSplineBasis], u: Sequence[Params]) -> Surface: """Perform a least-square fit of a point cloud `x` onto a spline basis. The points can be either a matrix (in which case the first index is diff --git a/src/splipy/typing.py b/src/splipy/typing.py index 6b89108..586f0b0 100644 --- a/src/splipy/typing.py +++ b/src/splipy/typing.py @@ -37,6 +37,7 @@ # Again, these are identical, but we can signal intent to a human reader based # on which we use. type Knots = Sequence[Scalar] | FloatArray | IntArray +type Params = Sequence[Scalar] | FloatArray | IntArray type Point = Sequence[Scalar] | FloatArray | IntArray type Direction = Literal["u", "v", "w", "U", "V", "W"] | int diff --git a/src/splipy/volume.py b/src/splipy/volume.py index 22ddd20..29e65a0 100644 --- a/src/splipy/volume.py +++ b/src/splipy/volume.py @@ -11,6 +11,8 @@ if TYPE_CHECKING: from collections.abc import Sequence + from splipy.typing import Point + from .curve import Curve from .surface import Surface from .typing import ArrayLike, Direction, FloatArray @@ -200,7 +202,7 @@ def __repr__(self) -> str: result += str(self.controlpoints[i, j, k, :]) + "\n" return result - def get_antiderivative_volume(self, direction: Direction, constant: ArrayLike | None = None) -> Volume: + def get_antiderivative_volume(self, direction: Direction, constant: Point | None = None) -> Volume: """Compute the antiderivative (integral) of the volume in a given parametric direction. The antiderivative is computed by inverting the derivative operator on From 3989e93177f732dcbc17a7eb3a3b228258605752 Mon Sep 17 00:00:00 2001 From: Eivind Fonn Date: Mon, 5 Jan 2026 11:27:03 +0100 Subject: [PATCH 20/22] Remove unused typevars --- src/splipy/typing.py | 4 +--- src/splipy/utils/__init__.py | 7 ++----- 2 files changed, 3 insertions(+), 8 deletions(-) diff --git a/src/splipy/typing.py b/src/splipy/typing.py index 586f0b0..ed2536f 100644 --- a/src/splipy/typing.py +++ b/src/splipy/typing.py @@ -1,13 +1,11 @@ from __future__ import annotations from collections.abc import Sequence -from typing import Literal, TypedDict, TypeVar +from typing import Literal, TypedDict import numpy as np import numpy.typing as npt -B = TypeVar("B", bound=npt.NBitBase) - # Anything that can be converted to an array. Use this type for parameters that # are (multidimensional) lists of points, such as controlpoints. Generally, only # use this for public-facing functions. Use np.asarray to convert to an actual diff --git a/src/splipy/utils/__init__.py b/src/splipy/utils/__init__.py index a96e5f5..7a69f5e 100644 --- a/src/splipy/utils/__init__.py +++ b/src/splipy/utils/__init__.py @@ -3,7 +3,7 @@ from collections.abc import Iterable, Iterator, Sequence, Sized from itertools import combinations, product, repeat from math import atan2, sqrt -from typing import TYPE_CHECKING, Any, Literal, SupportsFloat, TypeVar, Unpack, cast +from typing import TYPE_CHECKING, Any, Literal, SupportsFloat, Unpack, cast import numpy as np @@ -174,10 +174,7 @@ def is_singleton(x: Any) -> bool: return not isinstance(x, Sized) -T = TypeVar("T") - - -def ensure_listlike(x: T | Sequence[T], dups: int = 1) -> tuple[T, ...]: +def ensure_listlike[T](x: T | Sequence[T], dups: int = 1) -> tuple[T, ...]: if isinstance(x, Sequence): y = tuple(x) while len(y) < dups: From 2aa11d2f800428c92c6c82570ce16aa32e5d09ed Mon Sep 17 00:00:00 2001 From: Eivind Fonn Date: Mon, 5 Jan 2026 11:30:09 +0100 Subject: [PATCH 21/22] Remove dead code --- src/splipy/splineobject.py | 11 ----------- 1 file changed, 11 deletions(-) diff --git a/src/splipy/splineobject.py b/src/splipy/splineobject.py index 29d4e71..93b8b56 100644 --- a/src/splipy/splineobject.py +++ b/src/splipy/splineobject.py @@ -141,17 +141,6 @@ def _validate_domain(self, *params: FloatArray) -> None: if b.periodic < 0 and (np.min(p) < b.start() or b.end() < np.max(p)): raise ValueError("Evaluation outside parametric domain") - # TODO(Eivind): Remove this method - def _validate_domain_old(self, *params): # type: ignore[no-untyped-def] - """Check whether the given evaluation parameters are valid. - - :raises ValueError: If the parameters are outside the domain - """ - for b, p in zip(self.bases, params): - b.snap(p) - if b.periodic < 0 and (min(p) < b.start() or b.end() < max(p)): - raise ValueError("Evaluation outside parametric domain") - def evaluate( self, *params: Params | Scalar, From a474759fc91df6915689adc3890dffee9efedcd0 Mon Sep 17 00:00:00 2001 From: Eivind Fonn Date: Mon, 5 Jan 2026 14:11:59 +0100 Subject: [PATCH 22/22] Remove some unnecessary float casts --- src/splipy/basis.py | 37 +++++------ src/splipy/curve.py | 18 ++---- src/splipy/curve_factory.py | 2 +- src/splipy/io/g2.py | 109 ++++++++++++++++++--------------- src/splipy/io/spl.py | 2 +- src/splipy/io/svg.py | 2 +- src/splipy/splineobject.py | 16 ++--- src/splipy/surface.py | 8 +-- src/splipy/surface_factory.py | 4 +- src/splipy/utils/__init__.py | 18 +++--- src/splipy/utils/image.py | 8 +-- src/splipy/utils/nutils.py | 4 +- src/splipy/utils/refinement.py | 8 +-- src/splipy/volume.py | 14 ++--- src/splipy/volume_factory.py | 3 +- 15 files changed, 118 insertions(+), 135 deletions(-) diff --git a/src/splipy/basis.py b/src/splipy/basis.py index 62e0b0e..3ac56aa 100644 --- a/src/splipy/basis.py +++ b/src/splipy/basis.py @@ -120,7 +120,7 @@ def greville_all(self) -> FloatArray: operator = np.ones((self.order - 1,), dtype=np.float64) / (self.order - 1) return np.convolve(self.knots[1 : -1 - (self.periodic + 1)], operator, mode="valid") - def greville_single(self, index: Int) -> float: + def greville_single(self, index: Int) -> Scalar: """Fetch a greville point, also known as a knot averages: .. math:: \\sum_{j=i+1}^{i+p-1} \\frac{t_j}{p-1} @@ -130,15 +130,15 @@ def greville_single(self, index: Int) -> float: :return: A Greville point :rtype: float """ - return float(np.sum(self.knots[index + 1 : index + self.order]) / (self.order - 1)) + return np.sum(self.knots[index + 1 : index + self.order]) / (self.order - 1) @overload - def greville(self, index: Int) -> float: ... + def greville(self, index: Int) -> Scalar: ... @overload def greville(self) -> FloatArray: ... - def greville(self, index: Int | None = None) -> float | FloatArray: + def greville(self, index: Int | None = None) -> Scalar | FloatArray: """Fetch greville points, also known as knot averages: .. math:: \\sum_{j=i+1}^{i+p-1} \\frac{t_j}{p-1} @@ -266,10 +266,10 @@ def evaluate_old(self, t, d=0, from_right=True, sparse=False): # type: ignore[n for j in range(p - q - 1, p): k = mu - p + j # 'i'-index in global knot vector (ref Hughes book pg.21) if j != p - q - 1: - M[j] = M[j] * float(evalT - self.knots[k]) / (self.knots[k + q] - self.knots[k]) + M[j] = M[j] * (evalT - self.knots[k]) / (self.knots[k + q] - self.knots[k]) if j != p - 1: - M[j] = M[j] + M[j + 1] * float(self.knots[k + q + 1] - evalT) / ( + M[j] = M[j] + M[j + 1] * (self.knots[k + q + 1] - evalT) / ( self.knots[k + q + 1] - self.knots[k + 1] ) @@ -277,9 +277,9 @@ def evaluate_old(self, t, d=0, from_right=True, sparse=False): # type: ignore[n for j in range(p - q - 1, p): k = mu - p + j # 'i'-index in global knot vector (ref Hughes book pg.21) if j != p - q - 1: - M[j] = M[j] * float(q) / (self.knots[k + q] - self.knots[k]) + M[j] = M[j] * q / (self.knots[k + q] - self.knots[k]) if j != p - 1: - M[j] = M[j] - M[j + 1] * float(q) / (self.knots[k + q + 1] - self.knots[k + 1]) + M[j] = M[j] - M[j + 1] * q / (self.knots[k + q + 1] - self.knots[k + 1]) data[i * p : (i + 1) * p] = M indices[i * p : (i + 1) * p] = np.arange(mu - p, mu) % n @@ -334,9 +334,6 @@ def reparam(self, start: Scalar = 0, end: Scalar = 1) -> None: :raises ValueError: If *end* ≤ *start* """ - start = float(start) - end = float(end) - if end <= start: raise ValueError("end must be larger than start") self.normalize() @@ -359,8 +356,6 @@ def knot_continuity(self, knot: Scalar) -> int: knots. :rtype: int or float """ - knot = float(knot) - if self.periodic >= 0: if knot < self.start() or knot > self.end(): knot = (knot - self.start()) % (self.end() - self.start()) + self.start() @@ -377,7 +372,7 @@ def knot_continuity(self, knot: Scalar) -> int: raise NotAKnotError return self.order - (hi - lo) - 1 - def continuity(self, knot: Scalar) -> int | float: + def continuity(self, knot: Scalar) -> int | Scalar: """Get the continuity of the basis functions at a given point. :return: *p*--*m*--1 at a knot with multiplicity *m*, or ``inf`` @@ -512,8 +507,6 @@ def insert_knot(self, new_knot: Scalar) -> FloatArray: :rtype: numpy.array :raises ValueError: If the new knot is outside the domain """ - new_knot = float(new_knot) - if self.periodic >= 0: if new_knot < self.start() or new_knot > self.end(): new_knot = (new_knot - self.start()) % (self.end() - self.start()) + self.start() @@ -639,24 +632,24 @@ def __len__(self) -> int: """Returns the number of knots in this basis.""" return len(self.knots) - def __getitem__(self, i: int) -> float: + def __getitem__(self, i: int) -> Scalar: """Returns the knot at a given index.""" - return float(self.knots[i]) + return self.knots[i] # type: ignore[no-any-return] def __iadd__(self, a: Scalar) -> Self: - self.knots += float(a) + self.knots += a return self def __isub__(self, a: Scalar) -> Self: - self.knots -= float(a) + self.knots -= a return self def __imul__(self, a: Scalar) -> Self: - self.knots *= float(a) + self.knots *= a return self def __itruediv__(self, a: Scalar) -> Self: - self.knots /= float(a) + self.knots /= a return self __ifloordiv__ = __itruediv__ # integer division (should not distinguish) diff --git a/src/splipy/curve.py b/src/splipy/curve.py index 6842541..2d92b3d 100644 --- a/src/splipy/curve.py +++ b/src/splipy/curve.py @@ -220,7 +220,7 @@ def normal(self, t: Params | Scalar, above: bool = True) -> FloatArray: return np.cross(B, T) - def curvature(self, t: Params | Scalar, above: bool = True) -> FloatArray | float: + def curvature(self, t: Params | Scalar, above: bool = True) -> FloatArray | Scalar: """Evaluate the curvaure at specified point(s). The curvature is defined as .. math:: \\frac{|\\boldsymbol{v}\\times \\boldsymbol{a}|}{|\\boldsymbol{v}|^3} @@ -238,7 +238,7 @@ def curvature(self, t: Params | Scalar, above: bool = True) -> FloatArray | floa w = v[..., 0] * a[..., 1] - v[..., 1] * a[..., 0] if self.dimension == 2 else np.cross(v, a) if len(v.shape) == 1: # single evaluation point - return float(np.linalg.norm(w) / np.linalg.norm(v)) + return np.linalg.norm(w) / np.linalg.norm(v) magnitude: FloatArray = np.abs(w) if self.dimension == 2 else np.linalg.norm(w, axis=-1) speed: FloatArray = np.linalg.norm(v, axis=-1) @@ -269,11 +269,7 @@ def torsion(self, t: Params | Scalar, above: bool = True) -> FloatArray | float: da = self.derivative(t, d=3, above=above) w = np.cross(v, a) - # magnitude: float = float(np.linalg.norm(w)) - # magnitude: FloatArray = np.linalg.norm(w) if v.ndim == 1 else np.linalg.norm(w, axis=-1) - if v.ndim == 1: # single evaluation point - # magnitude = np.linalg.norm(w) dot: FloatArray = np.dot(w, a) return dot / np.linalg.norm(w) ** 2 @@ -368,7 +364,7 @@ def knot_continuity(self, knot: Scalar) -> int: """ return self.bases[0].knot_continuity(knot) - def continuity(self, knot: Scalar) -> int | float: + def continuity(self, knot: Scalar) -> int | Scalar: """Get the parametric continuity of the curve at a given point. Will return p-1-m, where m is the knot multiplicity and inf between knots""" return self.bases[0].continuity(knot) @@ -378,7 +374,7 @@ def get_kinks(self) -> FloatArray: continuity""" return np.array([k for k in self.knots(0) if self.continuity(k) < 1], dtype=np.float64) - def length(self, t0: Scalar | None = None, t1: Scalar | None = None) -> float: + def length(self, t0: Scalar | None = None, t1: Scalar | None = None) -> Scalar: """Computes the euclidian length of the curve in geometric space .. math:: \\int_{t_0}^{t_1} \\left\\| \\frac{d\\boldsymbol{x}}{dt}(t) \\right\\| \\, dt @@ -391,12 +387,10 @@ def length(self, t0: Scalar | None = None, t1: Scalar | None = None) -> float: (x, w) = np.polynomial.legendre.leggauss(quadrature_points) # keep only integration boundaries within given start (t0) and stop (t1) interval if t0 is not None: - t0 = float(t0) i = bisect_left(knots, t0) knots = np.insert(knots, i, t0) knots = knots[i:] if t1 is not None: - t1 = float(t1) i = bisect_right(knots, t1) knots = knots[:i] knots = np.insert(knots, i, t1) @@ -407,7 +401,7 @@ def length(self, t0: Scalar | None = None, t1: Scalar | None = None) -> float: w = np.array([w / 2 * (t1 - t0) for t0, t1 in zip(knots[:-1], knots[1:])], dtype=np.float64).flatten() dx = self.derivative(t) detJ = np.sqrt(np.sum(dx**2, axis=1)) - return float(np.dot(detJ, w)) + return np.dot(detJ, w) # type: ignore[no-any-return] def rebuild(self, p: int, n: int) -> Curve: """Creates an approximation to this curve by resampling it using a @@ -457,7 +451,7 @@ def _closest_point_linear_curve(self, pt: Point) -> tuple[FloatArray, float]: t = t1 return self(t), t - def closest_point(self, pt: Point, t0: Scalar | None = None) -> tuple[FloatArray, float]: + def closest_point(self, pt: Point, t0: Scalar | None = None) -> tuple[FloatArray, Scalar]: """Computes the closest point on this curve to a given point. This is done by newton iteration and is using the state variables `controlpoint_absolute_tolerance` and `controlpoint_relative_tolerance` to determine convergence; but limited to 15 iterations. diff --git a/src/splipy/curve_factory.py b/src/splipy/curve_factory.py index 8f43192..7dd1c02 100644 --- a/src/splipy/curve_factory.py +++ b/src/splipy/curve_factory.py @@ -331,7 +331,7 @@ def circle_segment( knot = knot_vector(0, theta, num_intervals=knot_spans, interior_reps=2, endpoint_reps=3) n = (knot_spans - 1) * 2 + 3 # number of control points needed - dt = float(theta) / knot_spans / 2 # angle step + dt = theta / knot_spans / 2 # angle step cp = np.column_stack( [ r * np.cos(np.arange(n) * dt), diff --git a/src/splipy/io/g2.py b/src/splipy/io/g2.py index bd3e147..a9f54cd 100644 --- a/src/splipy/io/g2.py +++ b/src/splipy/io/g2.py @@ -21,6 +21,7 @@ from types import TracebackType from splipy.curve import Curve + from splipy.typing import FloatArray from splipy.volume import Volume @@ -42,14 +43,23 @@ def read_next_param_range(self) -> tuple[float, float]: start, end = map(float, next(self.fstream).split()) return start, end + def read_next_bool(self) -> bool: + return next(self.fstream).strip() != "0" + + def read_next_float(self) -> float: + return float(next(self.fstream).strip()) + + def read_next_array(self) -> FloatArray: + return np.array(next(self.fstream).split(), dtype=float) + def circle(self) -> Curve: int(self.read_next_non_whitespace().strip()) - r = float(next(self.fstream).strip()) - center = np.array(next(self.fstream).split(), dtype=float) - normal = np.array(next(self.fstream).split(), dtype=float) - xaxis = np.array(next(self.fstream).split(), dtype=float) + r = self.read_next_float() + center = self.read_next_array() + normal = self.read_next_array() + xaxis = self.read_next_array() param = self.read_next_param_range() - reverse = next(self.fstream).strip() != "0" + reverse = self.read_next_bool() result = curve_factory.circle(r=r, center=center, normal=normal, xaxis=xaxis) result.reparam(param) @@ -59,13 +69,13 @@ def circle(self) -> Curve: def ellipse(self) -> Curve: int(self.read_next_non_whitespace().strip()) - r1 = float(next(self.fstream).strip()) - r2 = float(next(self.fstream).strip()) - center = np.array(next(self.fstream).split(), dtype=float) - normal = np.array(next(self.fstream).split(), dtype=float) - xaxis = np.array(next(self.fstream).split(), dtype=float) + r1 = self.read_next_float() + r2 = self.read_next_float() + center = self.read_next_array() + normal = self.read_next_array() + xaxis = self.read_next_array() param = self.read_next_param_range() - reverse = next(self.fstream).strip() != "0" + reverse = self.read_next_bool() result = curve_factory.ellipse(r1=r1, r2=r2, center=center, normal=normal, xaxis=xaxis) result.reparam(param) @@ -75,18 +85,15 @@ def ellipse(self) -> Curve: def line(self) -> Curve: int(self.read_next_non_whitespace().strip()) - start = np.array(next(self.fstream).split(), dtype=float) - direction = np.array(next(self.fstream).split(), dtype=float) - finite = next(self.fstream).strip() != "0" + start = self.read_next_array() + direction = self.read_next_array() + finite = self.read_next_bool() param = self.read_next_param_range() - reverse = next(self.fstream).strip() != "0" - d = np.array(direction) - s = np.array(start) - # d /= np.linalg.norm(d) + reverse = self.read_next_bool() if not finite: param = (-state.unlimited, state.unlimited) - result = curve_factory.line(s + d * param[0], s + d * param[1]) + result = curve_factory.line(start + direction * param[0], start + direction * param[1]) if reverse: result.reverse() return result @@ -105,14 +112,14 @@ def line(self) -> Curve: def cylinder(self) -> Surface: int(self.read_next_non_whitespace().strip()) - r = float(next(self.fstream).strip()) - center = np.array(next(self.fstream).split(), dtype=float) - z_axis = np.array(next(self.fstream).split(), dtype=float) - x_axis = np.array(next(self.fstream).split(), dtype=float) - finite = next(self.fstream).strip() != "0" + r = self.read_next_float() + center = self.read_next_array() + z_axis = self.read_next_array() + x_axis = self.read_next_array() + finite = self.read_next_bool() param_u = self.read_next_param_range() param_v = self.read_next_param_range() if finite else (-state.unlimited, state.unlimited) - swap = next(self.fstream).strip() != "0" + swap = self.read_next_bool() center = center + z_axis * param_v[0] h = param_v[1] - param_v[0] @@ -124,15 +131,15 @@ def cylinder(self) -> Surface: def disc(self) -> Surface: int(self.read_next_non_whitespace().strip()) - center = np.array(next(self.fstream).split(), dtype=float) - r = float(next(self.fstream).strip()) - z_axis = np.array(next(self.fstream).split(), dtype=float) - x_axis = np.array(next(self.fstream).split(), dtype=float) - degen = next(self.fstream).strip() != "0" - angles = [float(next(self.fstream).strip()) for i in range(4)] + center = self.read_next_array() + r = self.read_next_float() + z_axis = self.read_next_array() + x_axis = self.read_next_array() + degen = self.read_next_bool() + angles = [self.read_next_float() for _ in range(4)] param_u = self.read_next_param_range() param_v = self.read_next_param_range() - swap = next(self.fstream).strip() != "0" + swap = self.read_next_bool() if degen: result = surface_factory.disc(r=r, center=center, xaxis=x_axis, normal=z_axis, type="radial") @@ -147,17 +154,17 @@ def disc(self) -> Surface: def plane(self) -> Surface: int(self.read_next_non_whitespace().strip()) - center = np.array(next(self.fstream).split(), dtype=float) - normal = np.array(next(self.fstream).split(), dtype=float) - x_axis = np.array(next(self.fstream).split(), dtype=float) - finite = next(self.fstream).strip() != "0" + center = self.read_next_array() + normal = self.read_next_array() + x_axis = self.read_next_array() + finite = self.read_next_bool() if finite: param_u = self.read_next_param_range() param_v = self.read_next_param_range() else: param_u = (-state.unlimited, +state.unlimited) param_v = (-state.unlimited, +state.unlimited) - swap = next(self.fstream).strip() != "0" + swap = self.read_next_bool() result = Surface() * [param_u[1] - param_u[0], param_v[1] - param_v[0]] + [param_u[0], param_v[0]] result.rotate(rotate_local_x_axis(x_axis, normal)) @@ -169,15 +176,15 @@ def plane(self) -> Surface: def torus(self) -> Surface: int(self.read_next_non_whitespace().strip()) - r2 = float(next(self.fstream).strip()) - r1 = float(next(self.fstream).strip()) - center = np.array(next(self.fstream).split(), dtype=float) - z_axis = np.array(next(self.fstream).split(), dtype=float) - x_axis = np.array(next(self.fstream).split(), dtype=float) - next(self.fstream).strip() != "0" # I have no idea what this does :( + r2 = self.read_next_float() + r1 = self.read_next_float() + center = self.read_next_array() + z_axis = self.read_next_array() + x_axis = self.read_next_array() + self.read_next_bool() # I have no idea what this does :( param_u = self.read_next_param_range() param_v = self.read_next_param_range() - swap = next(self.fstream).strip() != "0" + swap = self.read_next_bool() result = surface_factory.torus(minor_r=r1, major_r=r2, center=center, normal=z_axis, xaxis=x_axis) result.reparam(param_u, param_v) @@ -187,13 +194,13 @@ def torus(self) -> Surface: def sphere(self) -> Surface: int(self.read_next_non_whitespace().strip()) - r = float(next(self.fstream).strip()) - center = np.array(next(self.fstream).split(), dtype=float) - z_axis = np.array(next(self.fstream).split(), dtype=float) - x_axis = np.array(next(self.fstream).split(), dtype=float) + r = self.read_next_float() + center = self.read_next_array() + z_axis = self.read_next_array() + x_axis = self.read_next_array() param_u = self.read_next_param_range() param_v = self.read_next_param_range() - swap = next(self.fstream).strip() != "0" + swap = self.read_next_bool() result = surface_factory.sphere(r=r, center=center, xaxis=x_axis, zaxis=z_axis).swap() if swap: @@ -228,10 +235,10 @@ def surface_of_linear_extrusion(self) -> Surface: int(self.read_next_non_whitespace().strip()) crv = self.splines(1) normal = np.array(self.read_next_non_whitespace().split(), dtype=float) - finite = next(self.fstream).strip() != "0" + finite = self.read_next_bool() param_u = self.read_next_param_range() param_v = self.read_next_param_range() if finite else (-state.unlimited, +state.unlimited) - swap = next(self.fstream).strip() != "0" + swap = self.read_next_bool() result = surface_factory.extrude(crv + normal * param_v[0], normal * (param_v[1] - param_v[0])) result.reparam(param_u, param_v) diff --git a/src/splipy/io/spl.py b/src/splipy/io/spl.py index 6f0988e..4965607 100644 --- a/src/splipy/io/spl.py +++ b/src/splipy/io/spl.py @@ -56,7 +56,7 @@ def read(self) -> list[SplineObject]: knots = [[float(k) for k in islice(lines, nkts)] for nkts in nknots] bases = [BSplineBasis(p, kts, -1) for p, kts in zip(orders, knots)] - cpts = np.array([float(k) for k in islice(lines, totcoeffs * physdim)]) + cpts = np.array(list(islice(lines, totcoeffs * physdim)), dtype=float) cpts = cpts.reshape(physdim, *(ncoeffs[::-1])).transpose() obj = SplineObject.construct_subclass(bases, cpts, rational=False, raw=True) diff --git a/src/splipy/io/svg.py b/src/splipy/io/svg.py index cb02f39..5e55495 100644 --- a/src/splipy/io/svg.py +++ b/src/splipy/io/svg.py @@ -152,7 +152,7 @@ def __exit__( # compute scaling factors by keeping aspect ratio, and never exceed # width or height size (including margins) - geometryRatio = float(boundingbox[3] - boundingbox[1]) / (boundingbox[2] - boundingbox[0]) + geometryRatio = (boundingbox[3] - boundingbox[1]) / (boundingbox[2] - boundingbox[0]) imageRatio = 1.0 * self.height / self.width if geometryRatio > imageRatio: # scale by y-coordinate marginPixels = self.height * self.margin diff --git a/src/splipy/splineobject.py b/src/splipy/splineobject.py index 93b8b56..cd92324 100644 --- a/src/splipy/splineobject.py +++ b/src/splipy/splineobject.py @@ -333,14 +333,14 @@ def get_derivative_spline(self, direction: Direction | None = None) -> SplineObj if self.bases[d].periodic < 0: C = np.zeros((n - 1, n)) for i in range(n - 1): - C[i, i] = -float(p) / (k[i + p + 1] - k[i + 1]) - C[i, i + 1] = float(p) / (k[i + p + 1] - k[i + 1]) + C[i, i] = -p / (k[i + p + 1] - k[i + 1]) + C[i, i + 1] = p / (k[i + p + 1] - k[i + 1]) else: C = np.zeros((n, n)) for i in range(n): ip1 = np.mod(i + 1, n) - C[i, i] = -float(p) / (k[i + p + 1] - k[i + 1]) - C[i, ip1] = float(p) / (k[i + p + 1] - k[i + 1]) + C[i, i] = -p / (k[i + p + 1] - k[i + 1]) + C[i, ip1] = p / (k[i + p + 1] - k[i + 1]) derivative_cps = np.tensordot(C, self.controlpoints, axes=(1, d)) derivative_cps = derivative_cps.transpose(transpose_fix(self.pardim, d)) @@ -1245,7 +1245,7 @@ def project(self, plane: Literal["", "x", "y", "z", "xy", "xz", "yz", "xyz"]) -> return self - def bounding_box(self) -> list[tuple[float, float]]: + def bounding_box(self) -> list[tuple[Scalar, Scalar]]: """Gets the bounding box of a spline object, computed from the control-point values. Could be inaccurate for rational splines. @@ -1257,12 +1257,12 @@ def bounding_box(self) -> list[tuple[float, float]]: """ dim = self.dimension - result: list[tuple[float, float]] = [] + result: list[tuple[Scalar, Scalar]] = [] for i in range(dim): result.append( ( - float(np.min(self.controlpoints[..., i])), - float(np.max(self.controlpoints[..., i])), + np.min(self.controlpoints[..., i]), + np.max(self.controlpoints[..., i]), ) ) return result diff --git a/src/splipy/surface.py b/src/splipy/surface.py index 460c3df..c13cd17 100644 --- a/src/splipy/surface.py +++ b/src/splipy/surface.py @@ -246,12 +246,8 @@ def derivative( return result - def area(self) -> float: + def area(self) -> Scalar: """Computes the area of the surface in geometric space""" - - w1: FloatArray - w2: FloatArray - # fetch integration points (x1, w1) = np.polynomial.legendre.leggauss(self.order(0) + 1) (x2, w2) = np.polynomial.legendre.leggauss(self.order(1) + 1) @@ -276,7 +272,7 @@ def area(self) -> float: J = du[..., 0] * dv[..., 1] - du[..., 1] * dv[..., 0] if self.dimension == 2 else np.cross(du, dv) J = np.sqrt(np.sum(J**2, axis=2)) if self.dimension == 3 else np.abs(J) - return float(w1.dot(J).dot(w2)) + return w1.dot(J).dot(w2) # type: ignore[no-any-return] def edges(self) -> tuple[Curve, Curve, Curve, Curve]: """Return the four edge curves in (parametric) order: umin, umax, vmin, vmax diff --git a/src/splipy/surface_factory.py b/src/splipy/surface_factory.py index c4c7909..12ad5f0 100644 --- a/src/splipy/surface_factory.py +++ b/src/splipy/surface_factory.py @@ -232,7 +232,7 @@ def torus( """ circle = curve_factory.circle(minor_r) circle.rotate(pi / 2, (1, 0, 0)) # flip up into xz-plane - circle.translate((float(major_r), 0, 0)) # move into position to spin around z-axis + circle.translate((major_r, 0, 0)) # move into position to spin around z-axis result = revolve(circle) result.rotate(rotate_local_x_axis(xaxis, normal)) @@ -669,7 +669,7 @@ def thicken(curve: Curve, amount: Scalar | Callable[..., Scalar]) -> Surface: left_points[i, 0] = x[i, 0] + v[i, 1] * dist # x at top left_points[i, 1] = x[i, 1] - v[i, 0] * dist # y at top else: - a = float(cast("Scalar", amount)) + a = cast("Scalar", amount) right_points[:, 0] = x[:, 0] - v[:, 1] * a # x at bottom right_points[:, 1] = x[:, 1] + v[:, 0] * a # y at bottom left_points[:, 0] = x[:, 0] + v[:, 1] * a # x at top diff --git a/src/splipy/utils/__init__.py b/src/splipy/utils/__init__.py index 7a69f5e..f5f7d1a 100644 --- a/src/splipy/utils/__init__.py +++ b/src/splipy/utils/__init__.py @@ -3,7 +3,7 @@ from collections.abc import Iterable, Iterator, Sequence, Sized from itertools import combinations, product, repeat from math import atan2, sqrt -from typing import TYPE_CHECKING, Any, Literal, SupportsFloat, Unpack, cast +from typing import TYPE_CHECKING, Any, Literal, Unpack, cast import numpy as np @@ -23,8 +23,8 @@ def knot_vector( - start: SupportsFloat = 0.0, - end: SupportsFloat | None = None, + start: Scalar = 0.0, + end: Scalar | None = None, num_intervals: int | None = None, interior_reps: int = 1, endpoint_reps: int = 1, @@ -34,17 +34,17 @@ def knot_vector( if end is None: assert num_intervals is not None - end = float(start) + num_intervals + end = start + num_intervals elif num_intervals is None: - num_intervals = int(float(end) - float(start)) + num_intervals = int(end - start) - def iter() -> Iterator[float]: - yield from repeat(float(start), endpoint_reps) + def iter() -> Iterator[Scalar]: + yield from repeat(start, endpoint_reps) for i in range(1, num_intervals): a = i / num_intervals - val = float(start) * (1 - a) + float(end) * a + val = start * (1 - a) + end * a yield from repeat(val, interior_reps) - yield from repeat(float(end), endpoint_reps) + yield from repeat(end, endpoint_reps) count = interior_reps * (num_intervals - 1) + endpoint_reps * 2 return np.fromiter(iter(), dtype=np.float64, count=count) diff --git a/src/splipy/utils/image.py b/src/splipy/utils/image.py index ec9debf..d436456 100644 --- a/src/splipy/utils/image.py +++ b/src/splipy/utils/image.py @@ -55,7 +55,7 @@ def get_corners( if M[0] == 0: dCand = abs(X[index - 1, 0] - X[i - 1, 0]) else: - m = float(M[1]) / M[0] + m = M[1] / M[0] dCand = abs(X[index - 1, 1] - m * X[index - 1, 0] + m * X[i - 1, 0] - X[i - 1, 1]) / sqrt( m**2 + 1 ) @@ -142,9 +142,7 @@ def image_curves(filename: str) -> list[Curve]: corners = get_corners(pts) # recompute corners, since previous sem might be smooth n = len(pts) - parpt: list[float] = list(range(n)) - for i in range(n): - parpt[i] = float(parpt[i]) / (n - 1) + parpt = np.linspace(0, 1, num=n) # the choice of knot vector is a tricky one. We'll go with the following strategy: # - cubic, p=3 curve @@ -180,7 +178,7 @@ def image_curves(filename: str) -> list[Curve]: # make it span [0,1] instead of [0,n-1] for i in range(len(knot)): - knot[i] /= float(n - 1) + knot[i] /= n - 1 # make it periodic since these are all closed curves knot[0] -= knot[-1] - knot[-5] diff --git a/src/splipy/utils/nutils.py b/src/splipy/utils/nutils.py index f2d2732..356749a 100644 --- a/src/splipy/utils/nutils.py +++ b/src/splipy/utils/nutils.py @@ -12,7 +12,7 @@ if TYPE_CHECKING: from splipy.splineobject import SplineObject - from splipy.typing import FloatArray + from splipy.typing import FloatArray, Scalar def controlpoints(spline: SplineObject) -> FloatArray: @@ -28,7 +28,7 @@ def controlpoints(spline: SplineObject) -> FloatArray: raise RuntimeError("Non-spline argument detected") -def multiplicities(spline: SplineObject) -> list[list[float]]: +def multiplicities(spline: SplineObject) -> list[list[Scalar]]: """Returns the multiplicity of the knots at all knot values as a 2D array for all parametric directions, for all knots""" return [ diff --git a/src/splipy/utils/refinement.py b/src/splipy/utils/refinement.py index 1364c7c..3e4f634 100644 --- a/src/splipy/utils/refinement.py +++ b/src/splipy/utils/refinement.py @@ -1,6 +1,6 @@ from __future__ import annotations -from typing import TYPE_CHECKING, SupportsFloat, cast +from typing import TYPE_CHECKING, cast __doc__ = "Implementation of various refinement schemes." @@ -14,14 +14,14 @@ from collections.abc import Sequence from splipy.splineobject import SplineObject - from splipy.typing import Direction, FloatArray + from splipy.typing import Direction, FloatArray, Scalar # TODO(Eivind): put control over these tolerances somewhere. Modstate in splipy # seems to be the place for it, but we can't let splipy.utils influence the # structure of splipy. -def knot_exists(existing_knots: FloatArray, new_knot: SupportsFloat) -> bool: - return bool(np.any(np.isclose(existing_knots, float(new_knot), atol=1e-7, rtol=1e-10))) +def knot_exists(existing_knots: FloatArray, new_knot: Scalar) -> bool: + return bool(np.any(np.isclose(existing_knots, new_knot, atol=1e-7, rtol=1e-10))) def geometric_refine( diff --git a/src/splipy/volume.py b/src/splipy/volume.py index 29e65a0..9998715 100644 --- a/src/splipy/volume.py +++ b/src/splipy/volume.py @@ -11,11 +11,11 @@ if TYPE_CHECKING: from collections.abc import Sequence - from splipy.typing import Point + from splipy.typing import Point, Scalar from .curve import Curve from .surface import Surface - from .typing import ArrayLike, Direction, FloatArray + from .typing import ArrayLike, Direction __all__ = ["Volume"] @@ -102,17 +102,13 @@ def faces( tuple(boundary_faces), ) - def volume(self) -> float: + def volume(self) -> Scalar: """Computes the volume of the object in geometric space""" - - w1: FloatArray - w2: FloatArray - w3: FloatArray - # fetch integration points (x1, w1) = np.polynomial.legendre.leggauss(self.order(0) + 1) (x2, w2) = np.polynomial.legendre.leggauss(self.order(1) + 1) (x3, w3) = np.polynomial.legendre.leggauss(self.order(2) + 1) + # map points to parametric coordinates (and update the weights) (knots1, knots2, knots3) = self.knots() u = np.array([(x1 + 1) / 2 * (t1 - t0) + t0 for t0, t1 in zip(knots1[:-1], knots1[1:])]) @@ -141,7 +137,7 @@ def volume(self) -> float: J = du[:, :, :, 0] * c1 - du[:, :, :, 1] * c2 + du[:, :, :, 2] * c3 - return float(np.abs(J).dot(w3).dot(w2).dot(w1)) + return np.abs(J).dot(w3).dot(w2).dot(w1) # type: ignore[no-any-return] def rebuild(self, p: int | Sequence[int], n: int | Sequence[int]) -> Volume: """Creates an approximation to this volume by resampling it using diff --git a/src/splipy/volume_factory.py b/src/splipy/volume_factory.py index ae2eb4b..05d1ef6 100644 --- a/src/splipy/volume_factory.py +++ b/src/splipy/volume_factory.py @@ -139,7 +139,6 @@ def revolve(surf: Surface, theta: Scalar = 2 * pi, axis: Point = (0, 0, 1)) -> V surf = surf.clone() # clone input surface, throw away old reference surf.set_dimension(3) # add z-components (if not already present) surf.force_rational() # add weight (if not already present) - theta = float(theta) axis_np = np.asarray(axis, dtype=float) @@ -192,7 +191,7 @@ def torus( disc = surface_factory.disc(minor_r, type=type) disc.rotate(pi / 2, (1, 0, 0)) # flip up into xz-plane - disc.translate((float(major_r), 0, 0)) # move into position to spin around z-axis + disc.translate((major_r, 0, 0)) # move into position to spin around z-axis result = revolve(disc) result.rotate(rotate_local_x_axis(xaxis, normal))