Numpy support for arrays with dimension¶

A Quantity object can have any numerical-like object as its value attribute, including numpy's ndarray.

Physipy support numpy for many functionnalties :

• common creation routines
• mathematical operations
• numpy's functions and universal functions
• comparison
• indexing and fancy indexing
• iterators

Creation¶

Basic creation of dimension-full arrays :

In [1]:

import matplotlib.pyplot as plt
import numpy as np

import physipy
from physipy import m, s, Quantity, Dimension, rad

import matplotlib.pyplot as plt import numpy as np import physipy from physipy import m, s, Quantity, Dimension, rad

In [2]:

x_samples = np.array([1, 2, 3, 4]) * m
y_samples = Quantity(np.array([1, 2, 3, 4]), Dimension("T"))
print(x_samples)
print(y_samples)
print(m*np.array([1, 2, 3, 4]) == x_samples) # multiplication is commutativ

x_samples = np.array([1, 2, 3, 4]) * m y_samples = Quantity(np.array([1, 2, 3, 4]), Dimension("T")) print(x_samples) print(y_samples) print(m*np.array([1, 2, 3, 4]) == x_samples) # multiplication is commutativ

[1 2 3 4] m
[1 2 3 4] s
[ True  True  True  True]

Operation¶

Basic array operation are handled the 'expected' way : note that the resulting dimension are consistent with the operation applied :

In [3]:

print(x_samples + 1*m)
print(x_samples * 2)
print(x_samples**2)
print(1/x_samples)

print(x_samples + 1*m) print(x_samples * 2) print(x_samples**2) print(1/x_samples)

[2 3 4 5] m
[2 4 6 8] m
[ 1  4  9 16] m**2
[1.         0.5        0.33333333 0.25      ] 1/m

Comparison¶

Comparison is allowed only for quantities that have the same units :

In [4]:

# allowed
print(x_samples > 1.5*m)

try: 
    # not allowed
    x_samples > 1.5*s
except Exception as e:
    print(e)

# allowed print(x_samples > 1.5*m) try: # not allowed x_samples > 1.5*s except Exception as e: print(e)

[False  True  True  True]
Dimension error : dimensions of operands are L and T, and are differents (length vs time).

Numpy ufuncs¶

Most numpy ufuncs are handled the expected way, but still check for dimension correctness :

In [5]:

q = 3*m
q_arr = np.arange(3)*m

print(np.add(q, q_arr))
print(np.multiply(q, q_arr))
print(np.sign(q_arr))
print(np.greater_equal(q_arr, 2*m))
print(np.sqrt(q_arr))
print(np.cbrt(q_arr))

print(np.cos(np.pi*rad))
print(np.tan(np.pi/4*rad))

print(np.ceil(q_arr**1.6))
print(np.negative(q_arr))

q = 3*m q_arr = np.arange(3)*m print(np.add(q, q_arr)) print(np.multiply(q, q_arr)) print(np.sign(q_arr)) print(np.greater_equal(q_arr, 2*m)) print(np.sqrt(q_arr)) print(np.cbrt(q_arr)) print(np.cos(np.pi*rad)) print(np.tan(np.pi/4*rad)) print(np.ceil(q_arr**1.6)) print(np.negative(q_arr))

[3 4 5] m
[0 3 6] m**2
[0 1 1]
[False False  True]
[0.         1.         1.41421356] m**0.5
[0.         1.         1.25992105] m**0.333333333333333
-1.0
0.9999999999999999
[0. 1. 4.] m**1.6
[ 0 -1 -2] m

Trigonometric functions expect dimensionless quantities, and regular dimension correctness is expected :

In [6]:

try:
    np.cos(3*m)
except Exception as e:
    print(e)

try:
    np.add(3*s, q_arr)
except Exception as e:
    print(e)

try: np.cos(3*m) except Exception as e: print(e) try: np.add(3*s, q_arr) except Exception as e: print(e)

Dimension error : dimensions of operands are L and no-dimension, and are differents (length vs dimensionless).
Dimension error : dimensions of operands are T and L, and are differents (time vs length).

Numpy's functions¶

Most classic numpy's functions are also handled :

In [7]:

print(np.linspace(3*m, 10*m, 5))
print(np.argmax(q_arr))
print(np.around(q_arr*2.3))
print(np.cross(q_arr, q_arr[::-1]))
print(np.dstack((q_arr, q_arr)))
print(np.mean(q_arr))
print(np.var(q_arr))
print(np.trapz(q_arr))
print(np.meshgrid(q_arr, q_arr))
print(np.fft.fft(q_arr))
print(np.convolve(q_arr, q_arr))
print(np.ravel(q_arr))
print(np.std(q_arr))
print(np.median(np.abs(q_arr-np.median(q_arr))))

print(np.linspace(3*m, 10*m, 5)) print(np.argmax(q_arr)) print(np.around(q_arr*2.3)) print(np.cross(q_arr, q_arr[::-1])) print(np.dstack((q_arr, q_arr))) print(np.mean(q_arr)) print(np.var(q_arr)) print(np.trapz(q_arr)) print(np.meshgrid(q_arr, q_arr)) print(np.fft.fft(q_arr)) print(np.convolve(q_arr, q_arr)) print(np.ravel(q_arr)) print(np.std(q_arr)) print(np.median(np.abs(q_arr-np.median(q_arr))))

[ 3.    4.75  6.5   8.25 10.  ] m
2 m
[0. 2. 5.] m
[-2  4 -2] m**2
[[[0 0]
  [1 1]
  [2 2]]] m
1.0 m
0.6666666666666666 m**2
2.0 m
(<Quantity : [[0 1 2]
 [0 1 2]
 [0 1 2]] m>, <Quantity : [[0 0 0]
 [1 1 1]
 [2 2 2]] m>)
[ 3. +0.j        -1.5+0.8660254j -1.5-0.8660254j] m
[0 0 1 4 4] m**2
[0 1 2] m
0.816496580927726 m
1.0 m

Reduce with ufuncs :

In [8]:

import numpy as np
from physipy import m
q = np.arange(10)*m

import numpy as np from physipy import m q = np.arange(10)*m

In [9]:

q = np.arange(10)*m
print(np.add.reduce(q))
print(np.multiply.reduce(q))

q = np.arange(10)*m print(np.add.reduce(q)) print(np.multiply.reduce(q))

45 m
0 m**10

Indexing¶

Indexing works just like with regular numpy arrays :

In [10]:

big_arr = np.arange(20).reshape(4,5)*s

print(big_arr)
print(big_arr[0])
print(big_arr[:, 2])

big_arr = np.arange(20).reshape(4,5)*s print(big_arr) print(big_arr[0]) print(big_arr[:, 2])

[[ 0  1  2  3  4]
 [ 5  6  7  8  9]
 [10 11 12 13 14]
 [15 16 17 18 19]] s
[0 1 2 3 4] s
[ 2  7 12 17] s

Fancy indexing¶

In [11]:

print(big_arr)
print(np.greater_equal(big_arr, 12*s))
print(big_arr[np.greater_equal(big_arr, 12*s)])

print(big_arr) print(np.greater_equal(big_arr, 12*s)) print(big_arr[np.greater_equal(big_arr, 12*s)])

[[ 0  1  2  3  4]
 [ 5  6  7  8  9]
 [10 11 12 13 14]
 [15 16 17 18 19]] s
[[False False False False False]
 [False False False False False]
 [False False  True  True  True]
 [ True  True  True  True  True]]
[12 13 14 15 16 17 18 19] s

Common array methods¶

flat iterator¶

In [12]:

print(big_arr.flat)

for q in q_arr.flat:
    print(q)

print(big_arr.flat) for q in q_arr.flat: print(q)

<physipy.quantity.quantity.FlatQuantityIterator object at 0x0000023DCF3F98B0>
0 m
1 m
2 m

Known issues¶

logical fucntions¶

The expected behavior of logical functions is not trivial :

• logical_and
• logical_or
• logical_xor
• logical_not

Hence they are not implemented.

np.arange¶

The commonly used np.arange cannot be overriden the same way the ufuncs or classic numpy function can be. Hence, a wrapped version is provided

In [13]:

from physipy.quantity.utils import qarange

from physipy.quantity.utils import qarange

In [14]:

try:
    np.arange(10*m)
except Exception as e:
    print(e)

try: np.arange(10*m) except Exception as e: print(e)

Dimension error : dimensions of operands are L and no-dimension, and are differents (length vs dimensionless).

In [15]:

# using range
print(np.array(range(10))*m)
# using np.arange
print(np.arange(10)*m)
# using physipy's qarange : note that the "m" quantity is inside the function call
print(qarange(10*m))

# using range print(np.array(range(10))*m) # using np.arange print(np.arange(10)*m) # using physipy's qarange : note that the "m" quantity is inside the function call print(qarange(10*m))

[0 1 2 3 4 5 6 7 8 9] m
[0 1 2 3 4 5 6 7 8 9] m
[0 1 2 3 4 5 6 7 8 9] m

With this wrapper, you can then do the following :

In [16]:

print(np.arange(2.5, 12)*m)
print(qarange(2.5*m, 12*m))

print(np.arange(2.5, 12)*m) print(qarange(2.5*m, 12*m))

[ 2.5  3.5  4.5  5.5  6.5  7.5  8.5  9.5 10.5 11.5] m
[ 2.5  3.5  4.5  5.5  6.5  7.5  8.5  9.5 10.5 11.5] m

The qarange wrapper still cares about dimension correctness :

In [17]:

try:
    print(qarange(2*m, 10*s))
except Exception as e:
    print(e)

try: print(qarange(2*m, 10*s)) except Exception as e: print(e)

Dimension error : dimensions of operands are L and T, and are differents (length vs time).

In [18]:

np.reshape(q_arr, (1, len(q_arr)))

np.reshape(q_arr, (1, len(q_arr)))

Out[18]:

[[0. 1. 2.]]$\,m$

numpy coverage¶

physipy makes Quantity objects interoperate with numpy through two dispatch mechanisms:

• ufuncs (np.add, np.sin, np.sqrt, ...) via __array_ufunc__;
• array functions (np.concatenate, np.unique, np.linalg.norm, the np.fft family, ...) via __array_function__.

What is supported is introspectable at runtime with physipy.numpy_coverage(), which compares the running numpy against what physipy implements. Everything below is generated from it, so it always reflects the installed versions.

In [1]:

import physipy

coverage = physipy.numpy_coverage()
print(coverage)

import physipy coverage = physipy.numpy_coverage() print(coverage)

physipy numpy coverage (numpy 2.2.5)
  ufuncs          :  63/ 86 implemented ( 73.3%), 4 n/a
  array functions : 169/326 implemented ( 51.8%), 9 n/a

Functions that cannot be implemented¶

A Quantity is a single array of magnitudes sharing one Dimension — every element carries the same physical dimension. Any numpy function whose natural output would be a single array with per-element heterogeneous dimensions therefore has no faithful representation, and is reported as not applicable (rather than merely missing).

The clearest case is the polynomial-coefficient family: a 1-D coefficient array is inherently heterogeneous (c_n carries y, c_{n-1} carries y/x, ...). That is exactly why np.polyfit returns a tuple of separate Quantity coefficients rather than a single array. For the same reason np.vander (each column is a different power of x) and np.roots (whose input is such a coefficient array) are not applicable.

When the result is representable, physipy uses one of two escapes: restrict to dimensionless input and raise DimensionError otherwise (np.cumprod), or return a tuple of homogeneous Quantity objects (np.polyfit).

In [2]:

from IPython.display import Markdown

# Full per-function breakdown, generated from the running numpy.
Markdown(coverage.to_markdown())

from IPython.display import Markdown # Full per-function breakdown, generated from the running numpy. Markdown(coverage.to_markdown())

Out[2]:

physipy numpy coverage (numpy 2.2.5)¶

Array functions¶

169/326 implemented (52%), 9 not applicable

• Implemented (169): allclose, amax, amin, angle, append, apply_along_axis, arange, argmax, argmin, argsort, argwhere, around, array_equal, array_split, asanyarray, atleast_1d, atleast_2d, atleast_3d, average, bincount, block, broadcast_arrays, broadcast_to, choose, clip, column_stack, compress, concat, concatenate, convolve, copy, copyto, corrcoef, count_nonzero, cov, cross, cumprod, cumsum, cumulative_prod, cumulative_sum, delete, diag, diagflat, diagonal, diff, dot, dsplit, dstack, ediff1d, empty_like, expand_dims, extract, fft.fft, fft.fft2, fft.fftn, fft.fftshift, fft.hfft, fft.ifft, fft.ifft2, fft.ifftn, fft.ifftshift, fft.ihfft, fft.irfft, fft.irfft2, fft.irfftn, fft.rfft, fft.rfft2, fft.rfftn, fix, flatnonzero, flip, fliplr, flipud, full, full_like, gradient, histogram, histogram2d, hsplit, hstack, imag, insert, interp, isclose, iscomplex, iscomplexobj, isreal, isrealobj, linalg.eig, linalg.inv, linalg.lstsq, linalg.norm, linspace, max, may_share_memory, mean, median, meshgrid, min, moveaxis, nan_to_num, nanargmax, nanargmin, nancumsum, nanmax, nanmean, nanmedian, nanmin, nanpercentile, nanprod, nanquantile, nanstd, nansum, nanvar, ndim, nonzero, ones_like, outer, pad, percentile, permute_dims, piecewise, place, polyfit, polyval, prod, ptp, put, put_along_axis, putmask, quantile, ravel, real, real_if_close, repeat, reshape, resize, roll, rollaxis, rot90, round, searchsorted, shape, sinc, sort, sort_complex, split, squeeze, stack, std, sum, swapaxes, take, take_along_axis, tile, transpose, trapezoid, trapz, tril, trim_zeros, triu, unique, unstack, var, vsplit, vstack, where, zeros, zeros_like
• Missing (157): all, any, apply_over_axes, argpartition, array, array2string, array_equiv, array_repr, array_str, asarray, asarray_chkfinite, ascontiguousarray, asfortranarray, asmatrix, astype, bartlett, base_repr, binary_repr, blackman, bmat, broadcast_shapes, busday_count, busday_offset, can_cast, common_type, correlate, datetime_as_string, datetime_data, diag_indices, diag_indices_from, digitize, einsum, einsum_path, empty, eye, fft.fftfreq, fft.rfftfreq, fill_diagonal, format_float_positional, format_float_scientific, from_dlpack, frombuffer, fromfile, fromfunction, fromiter, frompyfunc, fromregex, fromstring, genfromtxt, geomspace, get_include, get_printoptions, getbufsize, geterr, geterrcall, hamming, hanning, histogram_bin_edges, histogramdd, i0, identity, in1d, indices, info, inner, intersect1d, is_busday, isdtype, isfortran, isin, isneginf, isposinf, isscalar, issubdtype, iterable, ix_, kaiser, kron, lexsort, linalg.cholesky, linalg.cond, linalg.cross, linalg.det, linalg.diagonal, linalg.eigh, linalg.eigvals, linalg.eigvalsh, linalg.matmul, linalg.matrix_norm, linalg.matrix_power, linalg.matrix_rank, linalg.matrix_transpose, linalg.multi_dot, linalg.outer, linalg.pinv, linalg.qr, linalg.slogdet, linalg.solve, linalg.svd, linalg.svdvals, linalg.tensordot, linalg.tensorinv, linalg.tensorsolve, linalg.trace, linalg.vecdot, linalg.vector_norm, load, loadtxt, logspace, mask_indices, matrix_transpose, min_scalar_type, mintypecode, nancumprod, nested_iters, ones, packbits, partition, printoptions, promote_types, ravel_multi_index, require, result_type, row_stack, save, savetxt, savez, savez_compressed, select, set_printoptions, setbufsize, setdiff1d, seterr, seterrcall, setxor1d, shares_memory, show_config, show_runtime, size, tensordot, test, trace, tri, tril_indices, tril_indices_from, triu_indices, triu_indices_from, typename, union1d, unique_all, unique_counts, unique_inverse, unique_values, unpackbits, unravel_index, unwrap, vdot
• Not applicable (9): poly, polyadd, polyder, polydiv, polyint, polymul, polysub, roots, vander

Ufuncs¶

63/86 implemented (73%), 4 not applicable

• Implemented (63): absolute, add, arccos, arccosh, arcsin, arcsinh, arctan, arctan2, arctanh, cbrt, ceil, conjugate, copysign, cos, cosh, deg2rad, divide, equal, exp, exp2, expm1, fabs, floor, floor_divide, fmax, fmin, fmod, greater, greater_equal, hypot, isfinite, isinf, isnan, less, less_equal, log, log10, log1p, log2, logaddexp, logaddexp2, matmul, maximum, minimum, modf, multiply, negative, nextafter, not_equal, power, rad2deg, reciprocal, remainder, rint, sign, sin, sinh, sqrt, square, subtract, tan, tanh, trunc
• Missing (23): bitwise_and, bitwise_count, bitwise_or, bitwise_xor, degrees, divmod, float_power, frexp, gcd, heaviside, invert, isnat, lcm, ldexp, left_shift, matvec, positive, radians, right_shift, signbit, spacing, vecdot, vecmat
• Not applicable (4): logical_and, logical_not, logical_or, logical_xor

Proxy support for numpy.random functions¶

In [24]:

from physipy import calculus, s
import numpy as np

from physipy import calculus, s import numpy as np

For now you have to manually create random vectors since numpy's random functions do not support interface :

In [25]:

np.random.normal(1, 2, 10000)*s

np.random.normal(1, 2, 10000)*s

Out[25]:

[ 4.08696168 0.92381521 4.5491531 ... 3.45169288 -0.19269588 2.48078724]$\,s$