I’m new to Numba
and am attempting to vectorize a function. I’ve gone through the Generalized Universal Function API documentation but haven’t succeeded yet.
My function processes a 2D data
matrix, another 2D matrix of (x, y) coordinates, and a 1D array for a coordinate pair (x0, y0). It creates a submatrix X
using the local neighborhood of (x0, y0), decomposes it, and returns the first n_out
singular values. My goal is to modify this function to accept an xy
array of shape (N x 2).
import numpy as np
import numba
np.random.seed(4)
data = np.random.normal(size=(20, 5))
XY = np.random.normal(size=(20, 2))
@numba.njit(fastmath=True)
def euclidean_distance(XY, xy):
"""Compute the Euclidean distance between two arrays."""
dist = np.zeros(XY.shape[0])
for r in range(XY.shape[0]):
d = 0
for c in range(XY.shape[1]):
d += (xy[c] - XY[r, c]) ** 2
dist[r] = d
return np.sqrt(dist)
@numba.njit(fastmath=True)
def vectorize_me(data, XY, xy, n_out: int = 3):
distance = euclidean_distance(XY, xy)
X = data[distance < 0.5]
U, s, Vt = np.linalg.svd(X)
return s[:n_out]
xy_sample = XY[0]
vectorize_me(data, XY, xy_sample, n_out=3)
While a simple loop with @numba.njit(parallel=True)
could work for 2D coordinate pairs, I’m aiming to use xarray
’s apply_ufunc
, which I believe requires a guvectorize
approach.
My current attempt looks like this:
@numba.guvectorize(
[
(
numba.float64[:, :],
numba.float64[:, :],
numba.float64[:, :],
numba.int64,
numba.float64[:, :],
)
],
"(n,m),(n,o),(o),()->(n_out)",
)
def i_am_vectorized(data, XY, xy, n_out, result):
distance = euclidean_distance(XY, xy)
X = A[distance < 0.5]
s = np.linalg.svd(X)
for i in range(n_out):
result[i] = s[i]
i_am_vectorized(A, XY, XY, n_out)
which results in a NameError: undefined output symbols: n_out
.
However, citing from the GUF API mentioned above:
Typically, the size of all core dimensions in an output will be determined by the size of a core dimension with the same label in an input array. This is not a requirement, and it is possible to define a signature where a label comes up for the first time in an output, although some precautions must be taken when calling such a function. (emphasis is mine)
So I don’t really understand what I’m doing wrong