Python interface (GPU)¶
Quick-start examples¶
As mentioned in the Python GPU installation instructions, the easiest way to install the Python interface for cuFINUFFT is to run pip install cufinufft
.
Note: The interface to cuFINUFFT has changed between versions 1.3 and 2.2. Please see Migration to cuFINUFFT v2.2 for details.
Assuming cuFINUFFT has been installed, we will now consider how to calculate a 1D type 1 transform.
To manage the GPU and transfer to and from host and device, we will use the cupy
library.
Consequently, we start with a few import statements.
import cupy as cp
import cufinufft
We then proceed to setting up a few parameters.
# number of nonuniform points
M = 100000
# grid size
N = 200000
# generate positions for the nonuniform points and the coefficients
x_gpu = 2 * cp.pi * cp.random.uniform(size=M)
c_gpu = (cp.random.standard_normal(size=M)
+ 1J * cp.random.standard_normal(size=M))
Now that the data is prepared, we can call cuFINUFFT to compute the transform from the M
source points onto a grid of N
modes
f_gpu = cufinufft.nufft1d1(x_gpu, c_gpu, (N,))
# move results off the GPU
f = f_gpu.get()
The last line is needed here to copy the data from the device (GPU) to the host (CPU).
Now f
is a size-N
array containing the NUDFT of c
at the points x
.
This example, along with similar examples for other frameworks, can be found in python/cufinufft/examples/getting_started_*.py
.
Other possible calculations are possible by supplying different options during plan creation. See the full API documentation below for more information.
Full documentation¶
- cufinufft.nufft1d1(x, data, n_modes=None, out=None, eps=1e-06, isign=1, **kwargs)¶
1D type-1 (nonuniform to uniform) complex NUFFT
M-1 f[k1] = SUM c[j] exp(+/-i k1 x(j)) j=0 for -N1/2 <= k1 <= (N1-1)/2
- Parameters:
x (float[M]) – nonuniform points, in the interval [-pi, pi), values outside will be folded
c (complex[M] or complex[n_tr, M]) – source strengths.
n_modes (integer or integer tuple of length 1, optional) – number of uniform Fourier modes requested (N1, ). May be even or odd; in either case, modes k1 are integers satisfying -N1/2 <= k1 <= (N1-1)/2. Must be specified if
out
is not given.out (complex[N1] or complex[n_tr, N1], optional) – output array for Fourier mode values. If
n_modes
is specifed, the shape must match, otherwisen_modes
is inferred fromout
.eps (float, optional) – precision requested (>1e-16).
isign (int, optional) – if non-negative, uses positive sign in exponential, otherwise negative sign.
**kwargs (optional) – other options may be specified, see the
Plan
constructor for details.
Note
The output is written into the
out
array if supplied.- Returns:
The resulting array.
- Return type:
complex[N1] or complex[n_tr, N1]
Example (CuPy):
import cupy as cp import cufinufft # number of nonuniform points M = 100 # the nonuniform points x = 2 * cp.pi * cp.random.uniform(size=M) # their complex strengths c = (cp.random.standard_normal(size=M) + 1J * cp.random.standard_normal(size=M)) # desired number of Fourier modes N1 = 50 # calculate the type-1 NUFFT f = cufinufft.nufft1d1(x, c, (N1, ))
- cufinufft.nufft1d2(x, data, out=None, eps=1e-06, isign=-1, **kwargs)¶
1D type-2 (uniform to nonuniform) complex NUFFT
c[j] = SUM f[k1] exp(+/-i k1 x(j)) k1 for j = 0, ..., M-1, where the sum is over -N1/2 <= k1 <= (N1-1)/2
- Parameters:
x (float[M]) – nonuniform points, in the interval [-pi, pi), values outside will be folded
f (complex[N1] or complex[n_tr, N1]) – Fourier mode coefficients, where N1 may be even or odd. In either case the mode indices k1 satisfy -N1/2 <= k1 <= (N1-1)/2.
out (complex[M] or complex[n_tr, M], optional) – output array at targets.
eps (float, optional) – precision requested (>1e-16).
isign (int, optional) – if non-negative, uses positive sign in exponential, otherwise negative sign.
**kwargs (optional) – other options may be specified, see the
Plan
constructor for details.
Note
The output is written into the
out
array if supplied.- Returns:
The resulting array.
- Return type:
complex[M] or complex[n_tr, M]
Example (CuPy):
import cupy as cp import cufinufft # number of nonuniform points M = 100 # the nonuniform points x = 2 * cp.pi * cp.random.uniform(size=M) # number of Fourier modes N1 = 50 # the Fourier mode coefficients f = (cp.random.standard_normal(size=(N1, )) + 1J * cp.random.standard_normal(size=(N1, ))) # calculate the type-2 NUFFT c = cufinufft.nufft1d2(x, f)
- cufinufft.nufft2d1(x, y, data, n_modes=None, out=None, eps=1e-06, isign=1, **kwargs)¶
2D type-1 (nonuniform to uniform) complex NUFFT
M-1 f[k1, k2] = SUM c[j] exp(+/-i (k1 x(j) + k2 y(j))) j=0 for -N1/2 <= k1 <= (N1-1)/2, -N2/2 <= k2 <= (N2-1)/2
- Parameters:
x (float[M]) – nonuniform points, in the interval [-pi, pi), values outside will be folded
y (float[M]) – nonuniform points, in the interval [-pi, pi), values outside will be folded
c (complex[M] or complex[n_tr, M]) – source strengths.
n_modes (integer or integer tuple of length 2, optional) – number of uniform Fourier modes requested (N1, N2). May be even or odd; in either case, modes k1, k2 are integers satisfying -N1/2 <= k1 <= (N1-1)/2, -N2/2 <= k2 <= (N2-1)/2. Must be specified if
out
is not given.out (complex[N1, N2] or complex[n_tr, N1, N2], optional) – output array for Fourier mode values. If
n_modes
is specifed, the shape must match, otherwisen_modes
is inferred fromout
.eps (float, optional) – precision requested (>1e-16).
isign (int, optional) – if non-negative, uses positive sign in exponential, otherwise negative sign.
**kwargs (optional) – other options may be specified, see the
Plan
constructor for details.
Note
The output is written into the
out
array if supplied.- Returns:
The resulting array.
- Return type:
complex[N1, N2] or complex[n_tr, N1, N2]
Example (CuPy):
import cupy as cp import cufinufft # number of nonuniform points M = 100 # the nonuniform points x = 2 * cp.pi * cp.random.uniform(size=M) y = 2 * cp.pi * cp.random.uniform(size=M) # their complex strengths c = (cp.random.standard_normal(size=M) + 1J * cp.random.standard_normal(size=M)) # desired number of Fourier modes N1, N2 = 50, 75 # calculate the type-1 NUFFT f = cufinufft.nufft2d1(x, y, c, (N1, N2))
- cufinufft.nufft2d2(x, y, data, out=None, eps=1e-06, isign=-1, **kwargs)¶
2D type-2 (uniform to nonuniform) complex NUFFT
c[j] = SUM f[k1, k2] exp(+/-i (k1 x(j) + k2 y(j))) k1, k2 for j = 0, ..., M-1, where the sum is over -N1/2 <= k1 <= (N1-1)/2, -N2/2 <= k2 <= (N2-1)/2
- Parameters:
x (float[M]) – nonuniform points, in the interval [-pi, pi), values outside will be folded
y (float[M]) – nonuniform points, in the interval [-pi, pi), values outside will be folded
f (complex[N1, N2] or complex[n_tr, N1, N2]) – Fourier mode coefficients, where N1, N2 may be even or odd. In either case the mode indices k1, k2 satisfy -N1/2 <= k1 <= (N1-1)/2, -N2/2 <= k2 <= (N2-1)/2.
out (complex[M] or complex[n_tr, M], optional) – output array at targets.
eps (float, optional) – precision requested (>1e-16).
isign (int, optional) – if non-negative, uses positive sign in exponential, otherwise negative sign.
**kwargs (optional) – other options may be specified, see the
Plan
constructor for details.
Note
The output is written into the
out
array if supplied.- Returns:
The resulting array.
- Return type:
complex[M] or complex[n_tr, M]
Example (CuPy):
import cupy as cp import cufinufft # number of nonuniform points M = 100 # the nonuniform points x = 2 * cp.pi * cp.random.uniform(size=M) y = 2 * cp.pi * cp.random.uniform(size=M) # number of Fourier modes N1, N2 = 50, 75 # the Fourier mode coefficients f = (cp.random.standard_normal(size=(N1, N2)) + 1J * cp.random.standard_normal(size=(N1, N2))) # calculate the type-2 NUFFT c = cufinufft.nufft2d2(x, y, f)
- cufinufft.nufft3d1(x, y, z, data, n_modes=None, out=None, eps=1e-06, isign=1, **kwargs)¶
3D type-1 (nonuniform to uniform) complex NUFFT
M-1 f[k1, k2, k3] = SUM c[j] exp(+/-i (k1 x(j) + k2 y(j) + k3 z(j))) j=0 for -N1/2 <= k1 <= (N1-1)/2, -N2/2 <= k2 <= (N2-1)/2, -N3/2 <= k3 <= (N3-1)/2
- Parameters:
x (float[M]) – nonuniform points, in the interval [-pi, pi), values outside will be folded
y (float[M]) – nonuniform points, in the interval [-pi, pi), values outside will be folded
z (float[M]) – nonuniform points, in the interval [-pi, pi), values outside will be folded
c (complex[M] or complex[n_tr, M]) – source strengths.
n_modes (integer or integer tuple of length 3, optional) – number of uniform Fourier modes requested (N1, N2, N3). May be even or odd; in either case, modes k1, k2, k3 are integers satisfying -N1/2 <= k1 <= (N1-1)/2, -N2/2 <= k2 <= (N2-1)/2, -N3/2 <= k3 <= (N3-1)/2. Must be specified if
out
is not given.out (complex[N1, N2, N3] or complex[n_tr, N1, N2, N3], optional) – output array for Fourier mode values. If
n_modes
is specifed, the shape must match, otherwisen_modes
is inferred fromout
.eps (float, optional) – precision requested (>1e-16).
isign (int, optional) – if non-negative, uses positive sign in exponential, otherwise negative sign.
**kwargs (optional) – other options may be specified, see the
Plan
constructor for details.
Note
The output is written into the
out
array if supplied.- Returns:
The resulting array.
- Return type:
complex[N1, N2, N3] or complex[n_tr, N1, N2, N3]
Example (CuPy):
import cupy as cp import cufinufft # number of nonuniform points M = 100 # the nonuniform points x = 2 * cp.pi * cp.random.uniform(size=M) y = 2 * cp.pi * cp.random.uniform(size=M) z = 2 * cp.pi * cp.random.uniform(size=M) # their complex strengths c = (cp.random.standard_normal(size=M) + 1J * cp.random.standard_normal(size=M)) # desired number of Fourier modes N1, N2, N3 = 50, 75, 100 # calculate the type-1 NUFFT f = cufinufft.nufft3d1(x, y, z, c, (N1, N2, N3))
- cufinufft.nufft3d2(x, y, z, data, out=None, eps=1e-06, isign=-1, **kwargs)¶
3D type-2 (uniform to nonuniform) complex NUFFT
c[j] = SUM f[k1, k2, k3] exp(+/-i (k1 x(j) + k2 y(j) + k3 z(j))) k1, k2, k3 for j = 0, ..., M-1, where the sum is over -N1/2 <= k1 <= (N1-1)/2, -N2/2 <= k2 <= (N2-1)/2, -N3/2 <= k3 <= (N3-1)/2
- Parameters:
x (float[M]) – nonuniform points, in the interval [-pi, pi), values outside will be folded
y (float[M]) – nonuniform points, in the interval [-pi, pi), values outside will be folded
z (float[M]) – nonuniform points, in the interval [-pi, pi), values outside will be folded
f (complex[N1, N2, N3] or complex[n_tr, N1, N2, N3]) – Fourier mode coefficients, where N1, N2, N3 may be even or odd. In either case the mode indices k1, k2, k3 satisfy -N1/2 <= k1 <= (N1-1)/2, -N2/2 <= k2 <= (N2-1)/2, -N3/2 <= k3 <= (N3-1)/2.
out (complex[M] or complex[n_tr, M], optional) – output array at targets.
eps (float, optional) – precision requested (>1e-16).
isign (int, optional) – if non-negative, uses positive sign in exponential, otherwise negative sign.
**kwargs (optional) – other options may be specified, see the
Plan
constructor for details.
Note
The output is written into the
out
array if supplied.- Returns:
The resulting array.
- Return type:
complex[M] or complex[n_tr, M]
Example (CuPy):
import cupy as cp import cufinufft # number of nonuniform points M = 100 # the nonuniform points x = 2 * cp.pi * cp.random.uniform(size=M) y = 2 * cp.pi * cp.random.uniform(size=M) z = 2 * cp.pi * cp.random.uniform(size=M) # number of Fourier modes N1, N2, N3 = 50, 75, 100 # the Fourier mode coefficients f = (cp.random.standard_normal(size=(N1, N2, N3)) + 1J * cp.random.standard_normal(size=(N1, N2, N3))) # calculate the type-2 NUFFT c = cufinufft.nufft3d2(x, y, z, f)
- class cufinufft.Plan(nufft_type, n_modes, n_trans=1, eps=1e-06, isign=None, dtype='complex64', **kwargs)¶
A non-uniform fast Fourier transform (NUFFT) plan
The
Plan
class lets the user exercise more fine-grained control over the execution of an NUFFT. First, the plan is created with a certain set of parameters (type, mode configuration, precision, sign, number of simultaneous transforms, and so on). Then the nonuniform points are set (source or target depending on the type). Finally, the plan is executed on some data, yielding the desired output.- Parameters:
nufft_type (int) – type of NUFFT (1 or 2).
n_modes (tuple of ints) – the number of modes in each dimension (for example (50, 100)).
n_trans (int, optional) – number of transforms to compute.
eps (float, optional) – precision requested (>1e-16).
isign (int, optional) – if +1, uses the positive sign exponential, otherwise the negative sign; defaults to +1 for type 1 and to -1 for type 2.
dtype (string, optional) – the precision of the transform, ‘complex64’ or ‘complex128’.
**kwargs (optional) – additional options corresponding may be specified as keyword arguments. These include upsampling ratio
upsampfac
(default2.0
),gpu_method
(1: nonuniform points-driven, 2: shared memory),gpu_sort
(forgpu_method == 1
, 0: no sort, 1: sort),gpu_kerevalmeth
(0: direct exp(sqrt), Horner evaluation),gpu_device_id
(GPU ID),gpu_stream
(CUDA stream pointer) andmodeord
(0: CMCL-compatible mode ordering, 1: FFT-style mode ordering).
- property type¶
- property dtype¶
- property dim¶
- property n_modes¶
- property n_trans¶
- setpts(x, y=None, z=None, s=None, t=None, u=None)¶
Set the nonuniform points
For type 1, this sets the coordinates of the
M
nonuniform source points and for type 2, it sets the coordinates of theM
target points.The dimension of the plan determines the number of arguments supplied. For example, if
dim == 2
, we providex
andy
.- Parameters:
x (float[M]) – first coordinate of the nonuniform points (source for type 1, target for type 2).
y (float[M], optional) – second coordinate of the nonuniform points (source for type 1, target for type 2).
z (float[M], optional) – third coordinate of the nonuniform points (source for type 1, target for type 2).
- execute(data, out=None)¶
Execute the plan
Performs the NUFFT specified at plan instantiation with the points set by
setpts
. For type-1 and type-3 transforms, the input is a set of source strengths, while for a type-2 transform, it consists of an array of sizen_modes
. Ifn_trans
is greater than one,n_trans
inputs are expected, stacked along the first axis.- Parameters:
data (complex[M], complex[n_tr, M], complex[n_modes], or complex[n_tr, n_modes]) – The input source strengths (type 1) or source modes (type 2).
out (complex[n_modes], complex[n_tr, n_modes], complex[M], or complex[n_tr, M], optional) – The array where the output is stored. Must be of the correct size.
- Returns:
The output array of the transform(s).
- Return type:
complex[n_modes], complex[n_tr, n_modes], complex[M], or complex[n_tr, M]