rkstiff.grids

Grid construction utilities for spectral and transform methods

This module provides functions for constructing spatial and spectral grids used in Fourier, Chebyshev, and Hankel spectral methods. These grids form the basis of spatial discretizations for PDE solvers such as those implemented in rkstiff.

Overview

Each function constructs either:

  • A uniform grid for FFT-based spectral differentiation.

  • A Chebyshev–Gauss–Lobatto grid for non-periodic domains.

  • A radial grid for axisymmetric problems using the Hankel transform.

All grids and wavenumber sets are compatible with NumPy and SciPy FFTs and can be directly used with the derivative utilities in rkstiff.derivatives.

Contents

Functions

construct_r_kr_hankel(nr, rmax)

Construct Hankel transform grids for axisymmetric domains.

construct_x_cheb(n[, a, b])

Construct a 1D grid with Chebyshev–Gauss–Lobatto points.

construct_x_dx_cheb(n[, a, b])

Construct Chebyshev–Gauss–Lobatto grid and its differentiation matrix.

construct_x_kx_fft(n[, a, b])

Construct a uniform 1D spatial grid and FFT-compatible wavenumber grid.

construct_x_kx_rfft(n[, a, b])

Construct a uniform 1D spatial grid and rFFT-compatible wavenumber grid.

mirror_grid(r[, u, axis])

Mirror a radial grid (and optional function) to produce a symmetric domain.

rkstiff.grids.construct_x_kx_rfft(n, a=0.0, b=2 * np.pi)[source]

Construct a uniform 1D spatial grid and rFFT-compatible wavenumber grid.

Parameters:
  • n (int) – Number of grid points. Must be an even integer greater than 2.

  • a (float, optional) – Left endpoint of the spatial domain. Default is 0.0.

  • b (float, optional) – Right endpoint of the spatial domain. Default is .

Return type:

Tuple[ndarray, ndarray]

Returns:

  • x (np.ndarray) – Uniform spatial grid with n points in \([a, b)\).

  • kx (np.ndarray) – Spectral wavenumber grid with n/2 + 1 points (for use with rFFT).

Notes

The grid spacing and wavenumbers are given by:

\[\Delta x = \frac{b - a}{n}, \qquad k_x = 2\pi \, \mathrm{rfftfreq}(n, \Delta x).\]

Examples

>>> x, kx = construct_x_kx_rfft(128, a=0, b=2*np.pi)
>>> x.shape
(128,)
>>> kx.shape
(65,)
rkstiff.grids.construct_x_kx_fft(n, a=0.0, b=2 * np.pi)[source]

Construct a uniform 1D spatial grid and FFT-compatible wavenumber grid.

Parameters:
  • n (int) – Number of grid points. Must be an even integer greater than 2.

  • a (float, optional) – Left endpoint of the spatial domain. Default is 0.0.

  • b (float, optional) – Right endpoint of the spatial domain. Default is .

Return type:

Tuple[ndarray, ndarray]

Returns:

  • x (np.ndarray) – Uniform spatial grid with n points in \([a, b)\).

  • kx (np.ndarray) – Spectral wavenumber grid with n points (for use with FFT).

Notes

The uniform grid and Fourier frequencies are defined by:

\[\Delta x = \frac{b - a}{n}, \qquad k_x = 2\pi \, \mathrm{fftfreq}(n, \Delta x).\]

Examples

>>> x, kx = construct_x_kx_fft(128)
>>> x.shape, kx.shape
((128,), (128,))
rkstiff.grids.construct_x_cheb(n, a=-1.0, b=1.0)[source]

Construct a 1D grid with Chebyshev–Gauss–Lobatto points.

Parameters:
  • n (int) – Polynomial order (number of subintervals). The grid has n + 1 points.

  • a (float, optional) – Interval endpoints. Defaults are a = -1, b = 1.

  • b (float, optional) – Interval endpoints. Defaults are a = -1, b = 1.

Returns:

Grid of n + 1 Chebyshev points mapped to the interval \([a, b]\).

Return type:

np.ndarray

Notes

The Chebyshev points on \([-1, 1]\) are given by:

\[x_j = \cos\!\left(\frac{j\pi}{n}\right), \quad j = 0, 1, \dots, n,\]

which are then linearly mapped to \([a, b]\).

Examples

>>> x = construct_x_cheb(10, a=-1, b=1)
>>> len(x)
11
>>> x[0], x[-1]
(1.0, -1.0)
rkstiff.grids.construct_x_dx_cheb(n, a=-1, b=1)[source]

Construct Chebyshev–Gauss–Lobatto grid and its differentiation matrix.

Parameters:
  • n (int) – Polynomial order (number of subintervals). Produces n + 1 points.

  • a (float, optional) – Interval endpoints. Defaults are a = -1, b = 1.

  • b (float, optional) – Interval endpoints. Defaults are a = -1, b = 1.

Return type:

Tuple[ndarray, ndarray]

Returns:

  • x (np.ndarray) – Chebyshev grid points on \([a, b]\).

  • d_cheb_matrix (np.ndarray, shape (n+1, n+1)) – Differentiation matrix \(D\) such that D @ f df/dx.

Notes

The entries of the differentiation matrix are:

\[\begin{split}D_{ij} = \begin{cases} \dfrac{c_i}{c_j (x_i - x_j)}, & i \neq j, \\ -\sum_{k \neq i} D_{ik}, & i = j, \end{cases}\end{split}\]

where \(c_i = 2(-1)^i\) for endpoints and \(c_i = (-1)^i\) otherwise.

The matrix achieves spectral accuracy for smooth functions and satisfies the property \(D\mathbf{1} = 0\).

rkstiff.grids.construct_r_kr_hankel(nr, rmax)[source]

Construct Hankel transform grids for axisymmetric domains.

Parameters:
  • nr (int) – Number of radial grid points (≥ 4).

  • rmax (float) – Maximum radius of the domain.

Return type:

Tuple[ndarray, ndarray, ndarray, float]

Returns:

  • r (np.ndarray) – Radial grid points in \((0, r_\max]\).

  • kr (np.ndarray) – Spectral grid points in wavenumber space.

  • bessel_zeros (np.ndarray) – First nr zeros of \(J_0\).

  • jN (float) – The (nr+1)-th zero of \(J_0\) used for normalization.

Notes

The grid is defined from the zeros of the Bessel function \(J_0\):

\[r_i = \frac{j_i}{j_{N+1}} \, r_\max, \qquad k_i = \frac{j_i}{r_\max},\]

where \(j_i\) are the zeros of \(J_0\).

rkstiff.grids.mirror_grid(r, u=None, axis=-1)[source]

Mirror a radial grid (and optional function) to produce a symmetric domain.

Parameters:
  • r (np.ndarray) – Radial grid on \([0, r_\max]\).

  • u (np.ndarray, optional) – Function values at r. If provided, mirrored values are returned.

  • axis (int, optional) –

    Axis along which to mirror u. Default is -1.

    • -1 — Stack horizontally (for 1D)

    • 0 — Stack vertically (rows)

    • 1 — Stack horizontally (columns)

Return type:

Tuple[ndarray, Optional[ndarray]]

Returns:

  • rnew (np.ndarray) – Mirrored grid on \([-r_\max, r_\max]\).

  • unew (np.ndarray, optional) – Mirrored function values (if u was provided).

Notes

Useful for visualizing radially symmetric solutions or enforcing symmetric boundary conditions.

Examples

>>> r = np.array([0, 1, 2, 3])
>>> u = np.array([1, 2, 3, 4])
>>> rnew, unew = mirror_grid(r, u)
>>> rnew
array([-3, -2, -1, 0, 0, 1, 2, 3])
>>> unew
array([4, 3, 2, 1, 1, 2, 3, 4])