RFFT

Onesided FFT of a real signal, on the same array-agnostic surface as the spectral transforms. Any length runs on any backend: powers of two, smooth sizes, primes, and lengths far past a GPU block’s shared memory (a fused two-kernel factor pipeline takes over there, transparently).

specux.rfft(x, backend="auto", out=None)    # (..., n) real -> (..., n//2+1)

x = np.random.randn(8, 4096).astype(np.float32)
F = specux.rfft(x)                    # (8, 2049) complex64, unscaled

xt = torch.randn(8, 262144, device="cuda", requires_grad=True)
F = specux.rfft(xt)                   # resident, any length, differentiable
F.abs().sum().backward()              # analytic adjoint, stays on the GPU

xc = cupy.random.standard_normal((8, 4096), dtype=cupy.float32)
F = specux.rfft(xc)                   # on the GPU, on cupy's own stream

The forward transform is unscaled; irfft carries the 1 / n.

out= buffers

out= fills a caller-provided array of the result’s shape and dtype, which keeps steady-state loops allocation-free:

blocks = np.random.randn(100, 8, 4096).astype(np.float32)
F = np.empty((8, 2049), np.complex64)
for block in blocks:                  # one warm buffer for the whole stream
    specux.rfft(block, out=F)

On the autograd / torch.compile path out= is ignored (a graph node owns its output).

Large lengths

Past a single GPU block’s shared memory the transform runs as exactly two fused kernels: a factor stage over n = n1 * n2, then a mirror-tile stage that computes the sub-FFTs and stores the onesided bins in one pass (each block holds a column tile and its Hermitian mirror, so the split costs no extra trip through memory). Nothing changes at the call site:

x = torch.randn(4, 1_048_576, device="cuda")
F = specux.rfft(x)                    # one million points, on the GPU

Parameters

x: real rows shaped (..., n); leading dimensions are flattened for the kernels and restored on the result.
backend: "auto" follows the input’s device.
out: optional preallocated result buffer (copy semantics).

The configured form

x = np.random.randn(8, 4096).astype(np.float32)
t = specux.RFFT()                     # composes with the other Transform
F = t(x)                              # objects in serializable pipelines

For repeated same-length transforms, fft_plan binds the length once.

References

Cooley-Tukey [cooley1965], Rader prime-length [rader1968], and Bluestein chirp-z [bluestein1970]; see References.