pulsesuite.libpulsesuite.nrutils¶

High-performance numerical utilities for scientific computing, ported from Fortran’s nrutils.F90.

1:1 port of the Numerical Recipes nrutils module. Every public Fortran routine has a Python callable with the same camelCase name. Vectorised with NumPy; no Numba or guardrails required.

Fortran interfaces that dispatch on type/rank are replaced by a single Python function that inspects shape and dtype at runtime.

Author: Rahul R. Sah

Attributes¶

`sp`
`dp`

Functions¶

`arrayCopy`(→ Tuple[int, int])	Copy elements from src into dest up to the smaller size.
`swap`(→ Tuple[Any, Any])	Swap two scalars, arrays, or any objects.
`maskedSwap`(→ None)	Swap elements of a and b where mask is True (in-place).
`reallocate`(→ numpy.ndarray)	Reallocate arr to newShape, copying the overlapping region.
`imaxloc`(→ int)	Index of the maximum value (0-based).
`iminloc`(→ int)	Index of the minimum value (0-based).
`ifirstloc`(→ int)	Index of the first True in a boolean array (0-based).
`nrerror`(→ None)	Report a fatal error and stop (Fortran's `nrerror`).
`assertTrue`(→ None)	Assert test is `True`; call `nrerror` on failure.
`assertEq`(→ Any)	Assert all positional arguments are equal; return the common value.
`arth`(→ numpy.ndarray)	Arithmetic progression of length n.
`geop`(→ numpy.ndarray)	Geometric progression of length n.
`cumsum`(→ numpy.ndarray)	Cumulative sum, optionally offset by seed.
`cumprod`(→ numpy.ndarray)	Cumulative product, optionally scaled by seed.
`poly`(→ Union[float, numpy.ndarray])	Evaluate polynomial at x via Horner's method.
`polyTerm`(→ numpy.ndarray)	Recursive polynomial term (Fortran's `poly_term`).
`outerprod`(→ numpy.ndarray)	Outer product: $M_{ij} = a_i cdot b_j$.
`outerdiff`(→ numpy.ndarray)	Outer difference: $M_{ij} = a_i - b_j$.
`outersum`(→ numpy.ndarray)	Outer sum: $M_{ij} = a_i + b_j$.
`outerand`(→ numpy.ndarray)	Outer logical AND: $M_{ij} = a_i wedge b_j$.
`scatterAdd`(→ None)	Scatter-add: `dest[i] += source[j]` at 0-based indices.
`scatterMax`(→ None)	Scatter-max: `dest[i] = max(dest[i], source[j])` at 0-based indices.
`diagAdd`(→ None)	Add scalar or vector diag to the diagonal of mat in-place.
`diagMult`(→ None)	Multiply the diagonal of mat by scalar or vector diag in-place.
`getDiag`(→ numpy.ndarray)	Extract diagonal of mat into a new 1-D array.
`putDiag`(→ None)	Set the diagonal of mat to diagv in-place.
`unitMatrix`(→ None)	Set mat to the identity matrix in-place ($I_{ij} = delta_{ij}$).
`upperTriangle`(→ numpy.ndarray)	Float mask for the upper triangle of a `(j, k)` matrix.
`lowerTriangle`(→ numpy.ndarray)	Float mask for the lower triangle of a `(j, k)` matrix.
`vabs`(→ float)	Euclidean norm $\|v\|_2 = sqrt{sum_i v_i^2}$.
`dummy_jacobian_dp`(...)	Zero-Jacobian placeholder for real (float64) ODE systems.
`dummy_jacobian_dpc`(...)	Zero-Jacobian placeholder for complex (complex128) ODE systems.

Module Contents¶

pulsesuite.libpulsesuite.nrutils.sp¶

pulsesuite.libpulsesuite.nrutils.dp¶

pulsesuite.libpulsesuite.nrutils.arrayCopy(src: numpy.ndarray, dest: numpy.ndarray) → Tuple[int, int]¶

Copy elements from src into dest up to the smaller size.

Mirrors Fortran’s array_copy.

Parameters:

src (ndarray) – Source array.
dest (ndarray) – Destination array (modified in-place).

Returns:

nCopied (int) – Number of elements actually copied.
nNotCopied (int) – src.size - nCopied.

pulsesuite.libpulsesuite.nrutils.swap(a: Any, b: Any) → Tuple[Any, Any]¶

Swap two scalars, arrays, or any objects.

Fortran’s swap_* interfaces (scalar, vector, matrix) are handled by Python’s generic assignment.

Parameters:

a (Any) – Values to swap.
b (Any) – Values to swap.

Returns:

(b, a) – Swapped pair.

Return type:

tuple

pulsesuite.libpulsesuite.nrutils.maskedSwap(a: numpy.ndarray, b: numpy.ndarray, mask: numpy.ndarray) → None¶

Swap elements of a and b where mask is True (in-place).

Fortran’s masked_swap_* interfaces. Uses vectorised boolean indexing — no Python loop.

Parameters:

a (ndarray (same shape, modified in-place))
b (ndarray (same shape, modified in-place))
mask (ndarray of bool) – Boolean mask; swap only where True.

pulsesuite.libpulsesuite.nrutils.reallocate(arr: numpy.ndarray, newShape: int | Tuple[int, Ellipsis]) → numpy.ndarray¶

Reallocate arr to newShape, copying the overlapping region.

Fortran’s reallocate_* family (1-D, 2-D, etc.).

Parameters:

arr (ndarray) – Original array.
newShape (int or tuple of int) – Desired shape.

Returns:

New array (same dtype) with data copied from the overlap.

Return type:

ndarray

pulsesuite.libpulsesuite.nrutils.imaxloc(arr: numpy.ndarray) → int¶: Index of the maximum value (0-based).

pulsesuite.libpulsesuite.nrutils.iminloc(arr: numpy.ndarray) → int¶: Index of the minimum value (0-based).

pulsesuite.libpulsesuite.nrutils.ifirstloc(mask: numpy.ndarray) → int¶

Index of the first True in a boolean array (0-based).

Returns -1 when no True is found.

pulsesuite.libpulsesuite.nrutils.nrerror(msg: str) → None¶

Report a fatal error and stop (Fortran’s nrerror).

Parameters:: msg (str) – Error message printed to stderr.
Raises:: SystemExit –

pulsesuite.libpulsesuite.nrutils.assertTrue(test: bool, msg: str = 'Assertion failed', file: str | None = None, line: int | None = None) → None¶

Assert test is True; call nrerror on failure.

Parameters:

test (bool)
msg (str)
file (optional) – Source location for diagnostics.
line (optional) – Source location for diagnostics.

pulsesuite.libpulsesuite.nrutils.assertEq(*args: Any, msg: str = 'Equality assertion failed') → Any¶

Assert all positional arguments are equal; return the common value.

Replaces Fortran’s assert_eq2, assert_eq3, assert_eq4, and assert_eqn interfaces with a single variadic dispatcher.

Parameters:

*args – Values that must be equal.
msg (str) – Message on failure.

Returns:

The common value.

Return type:

value

pulsesuite.libpulsesuite.nrutils.arth(first: float, increment: float, n: int) → numpy.ndarray¶

Arithmetic progression of length n.

$a_k = text{first} + k cdot text{increment},quad k = 0, dots, n-1$

Fortran’s arth (SP/DP/I4B interfaces).

Parameters:

first (float) – Starting value.
increment (float) – Common difference.
n (int) – Length.

Return type:

ndarray of float64

pulsesuite.libpulsesuite.nrutils.geop(first: float, factor: float, n: int) → numpy.ndarray¶

Geometric progression of length n.

$a_k = text{first} cdot text{factor}^k,quad k = 0, dots, n-1$

Fortran’s geop (SP/DP interfaces).

Parameters:

first (float) – Starting value.
factor (float) – Common ratio.
n (int) – Length.

Return type:

ndarray of float64

pulsesuite.libpulsesuite.nrutils.cumsum(arr: numpy.ndarray, seed: float | int | None = None) → numpy.ndarray¶

Cumulative sum, optionally offset by seed.

$S_j = text{seed} + sum_{i=0}^{j} a_i$

Fortran’s cumsum (SP/DP interfaces).

Parameters:

arr (ndarray)
seed (float or int, optional) – Additive offset applied to the running total.

Return type:

ndarray

pulsesuite.libpulsesuite.nrutils.cumprod(arr: numpy.ndarray, seed: float | int | None = None) → numpy.ndarray¶

Cumulative product, optionally scaled by seed.

$P_j = text{seed} cdot prod_{i=0}^{j} a_i$

Fortran’s cumprod (SP/DP interfaces).

Parameters:

arr (ndarray)
seed (float or int, optional) – Multiplicative scale applied to all products.

Return type:

ndarray

pulsesuite.libpulsesuite.nrutils.poly(x: float | numpy.ndarray, coeffs: numpy.ndarray) → float | numpy.ndarray¶

Evaluate polynomial at x via Horner’s method.

$P(x) = c_0 + c_1 x + c_2 x^2 + cdots$

Coefficients are ordered lowest-degree first (Fortran convention). Handles both real and complex coefficients/arguments.

Fortran’s poly (SP/DP/SPC/DPC interfaces).

Parameters:

x (float or ndarray) – Evaluation point(s).
coeffs (ndarray) – Polynomial coefficients, lowest degree first.

Return type:

float or ndarray

pulsesuite.libpulsesuite.nrutils.polyTerm(a: numpy.ndarray, b: float | int) → numpy.ndarray¶

Recursive polynomial term (Fortran’s poly_term).

$u_0 = a_0,quad u_j = a_j + b cdot u_{j-1}$

Handles both real and complex arrays.

Fortran’s poly_term (SP/DP interfaces).

Parameters:

a (ndarray) – Coefficient array.
b (float or int) – Multiplier.

Return type:

ndarray

pulsesuite.libpulsesuite.nrutils.outerprod(a: numpy.ndarray, b: numpy.ndarray) → numpy.ndarray¶

Outer product: $M_{ij} = a_i cdot b_j$.

Uses np.outer (BLAS-backed for large arrays). Handles real and complex vectors.

Fortran’s outerprod (SP/DP interfaces).

Parameters:

a (ndarray (1-D))
b (ndarray (1-D))

Return type:

ndarray, shape (a.size, b.size)

pulsesuite.libpulsesuite.nrutils.outerdiff(a: numpy.ndarray, b: numpy.ndarray) → numpy.ndarray¶

Outer difference: $M_{ij} = a_i - b_j$.

Vectorised via broadcasting.

Fortran’s outerdiff (SP/DP interfaces).

Parameters:

a (ndarray (1-D))
b (ndarray (1-D))

Return type:

ndarray, shape (a.size, b.size)

pulsesuite.libpulsesuite.nrutils.outersum(a: numpy.ndarray, b: numpy.ndarray) → numpy.ndarray¶

Outer sum: $M_{ij} = a_i + b_j$.

Vectorised via broadcasting.

Fortran’s outersum (SP/DP interfaces).

Parameters:

a (ndarray (1-D))
b (ndarray (1-D))

Return type:

ndarray, shape (a.size, b.size)

pulsesuite.libpulsesuite.nrutils.outerand(a: numpy.ndarray, b: numpy.ndarray) → numpy.ndarray¶

Outer logical AND: $M_{ij} = a_i wedge b_j$.

Fortran’s outerand.

Parameters:

a (ndarray of bool (1-D))
b (ndarray of bool (1-D))

Return type:

ndarray of bool, shape (a.size, b.size)

pulsesuite.libpulsesuite.nrutils.scatterAdd(dest: numpy.ndarray, source: numpy.ndarray, destIndex: numpy.ndarray) → None¶

Scatter-add: dest[i] += source[j] at 0-based indices.

Uses np.add.at for unbuffered accumulation (correct with duplicate indices).

Fortran’s scatter_add (SP/DP interfaces).

Parameters:

dest (ndarray (modified in-place))
source (ndarray)
destIndex (ndarray of int) – 0-based target indices.

pulsesuite.libpulsesuite.nrutils.scatterMax(dest: numpy.ndarray, source: numpy.ndarray, destIndex: numpy.ndarray) → None¶

Scatter-max: dest[i] = max(dest[i], source[j]) at 0-based indices.

Uses np.maximum.at for unbuffered reduction.

Fortran’s scatter_max (SP/DP interfaces).

Parameters:

dest (ndarray (modified in-place))
source (ndarray)
destIndex (ndarray of int) – 0-based target indices.

pulsesuite.libpulsesuite.nrutils.diagAdd(mat: numpy.ndarray, diag: float | numpy.ndarray) → None¶

Add scalar or vector diag to the diagonal of mat in-place.

Dispatcher replacing Fortran’s diagadd interface (scalar / vector).

Parameters:

mat (ndarray, shape (M, N) — modified in-place)
diag (float or ndarray)

pulsesuite.libpulsesuite.nrutils.diagMult(mat: numpy.ndarray, diag: float | numpy.ndarray) → None¶

Multiply the diagonal of mat by scalar or vector diag in-place.

Dispatcher replacing Fortran’s diagmult interface (scalar / vector).

Parameters:

mat (ndarray, shape (M, N) — modified in-place)
diag (float or ndarray)

pulsesuite.libpulsesuite.nrutils.getDiag(mat: numpy.ndarray) → numpy.ndarray¶

Extract diagonal of mat into a new 1-D array.

Fortran’s get_diag (SP/DP interfaces).

Parameters:: mat (ndarray, shape (M, N))
Return type:: ndarray, shape (min(M, N),)

pulsesuite.libpulsesuite.nrutils.putDiag(diagv: numpy.ndarray, mat: numpy.ndarray) → None¶

Set the diagonal of mat to diagv in-place.

Fortran’s put_diag (SP/DP interfaces).

Parameters:

diagv (ndarray)
mat (ndarray, shape (M, N) — modified in-place)

pulsesuite.libpulsesuite.nrutils.unitMatrix(mat: numpy.ndarray) → None¶

Set mat to the identity matrix in-place ($I_{ij} = delta_{ij}$).

Fortran’s unit_matrix.

Parameters:: mat (ndarray, shape (M, N) — modified in-place)

pulsesuite.libpulsesuite.nrutils.upperTriangle(j: int, k: int, extra: int = 0) → numpy.ndarray¶

Float mask for the upper triangle of a (j, k) matrix.

Entry $(r, c)$ is 1.0 where $r < c + text{extra}$, else 0.0.

Fortran’s upper_triangle.

Parameters:

j (int — rows)
k (int — columns)
extra (int, optional) – Diagonal offset (default 0).

Return type:

ndarray of float64, shape (j, k)

pulsesuite.libpulsesuite.nrutils.lowerTriangle(j: int, k: int, extra: int = 0) → numpy.ndarray¶

Float mask for the lower triangle of a (j, k) matrix.

Entry $(r, c)$ is 1.0 where $r > c - text{extra}$, else 0.0.

Fortran’s lower_triangle.

Parameters:

j (int — rows)
k (int — columns)
extra (int, optional) – Diagonal offset (default 0).

Return type:

ndarray of float64, shape (j, k)

pulsesuite.libpulsesuite.nrutils.vabs(v: numpy.ndarray) → float¶

Euclidean norm $|v|_2 = sqrt{sum_i v_i^2}$.

Uses np.linalg.norm (BLAS-backed).

Fortran’s vabs.

Parameters:: v (ndarray)
Return type:: float

pulsesuite.libpulsesuite.nrutils.dummy_jacobian_dp(x: float, y: numpy.typing.NDArray[numpy.float64]) → Tuple[numpy.typing.NDArray[numpy.float64], numpy.typing.NDArray[numpy.float64]]¶

Zero-Jacobian placeholder for real (float64) ODE systems.

Returns $(mathbf{0},; mathbf{0})$ with dimensions matching y.

Parameters:

x (float) – Independent variable (unused).
y (ndarray of float64, shape (n,)) – State vector.

Returns:

dfdx (ndarray of float64, shape (n,))
dfdy (ndarray of float64, shape (n, n))

pulsesuite.libpulsesuite.nrutils.dummy_jacobian_dpc(x: float, y: numpy.typing.NDArray[numpy.complex128]) → Tuple[numpy.typing.NDArray[numpy.complex128], numpy.typing.NDArray[numpy.complex128]]¶

Zero-Jacobian placeholder for complex (complex128) ODE systems.

Returns $(mathbf{0},; mathbf{0})$ with dimensions matching y.

Parameters:

x (float) – Independent variable (unused).
y (ndarray of complex128, shape (n,)) – State vector.

Returns:

dfdx (ndarray of complex128, shape (n,))
dfdy (ndarray of complex128, shape (n, n))