pygeoinf.symmetric_space package¶

Submodules¶

pygeoinf.symmetric_space.circle module¶

class pygeoinf.symmetric_space.circle.Lebesgue(kmax: int, /, *, radius: float = 1.0)[source]¶

Bases: AbstractSymmetricLebesgueSpace

Implementation of the Lebesgue space L² on a circle.

This class represents square-integrable functions on a circle. A function is represented by its values on an evenly spaced grid. The co-ordinate basis for the space is through Fourier expansions.

degree_multiplicity(degree: int) → int[source]¶: Returns the number of eigenfunctions associated with a given degree.

degree_transfer_operator(target_degree: int) → LinearOperator[source]¶

Returns the transfer operator from this space to one with a different degree.

This operator upsamples (by zero-padding Fourier coefficients) or downsamples (by truncating Fourier coefficients) the function grid.

property fft_factor: float¶: The factor by which the Fourier coefficients are scaled in forward transformations.

from_coefficients(coeff: ndarray) → ndarray[source]¶: Maps complex Fourier coefficients to a function vector.

from_components(c: ndarray) → ndarray[source]¶

Maps a NumPy component array back to a vector in the space.

Parameters:: c – The components of the vector as a NumPy array.
Returns:: The corresponding vector in the space.

classmethod from_covariance(covariance_function: Callable[[float], float], /, *, radius: float = 1.0, rtol: float = 1e-06, min_degree: int = 0, max_degree: int | None = None, power_of_two: bool = False) → Lebesgue[source]¶

Creates an instance of the L² space with a Fourier truncation degree (kmax) automatically chosen to capture the expected energy of functions drawn from a specified prior measure.

This factory method calculates the expected squared norm (energy) of a random field whose spectral variances are defined by the provided covariance_function. It iteratively adds higher frequency modes until the relative contribution of the next degree drops below the specified relative tolerance.

Parameters:

covariance_function – A callable mapping a Laplacian eigenvalue to its spectral variance.
radius – The radius of the circle. Defaults to 1.0.
rtol – The relative tolerance for the energy truncation. Defaults to 1e-6.
min_degree – The absolute minimum truncation degree to return. Defaults to 0.
max_degree – An optional safety ceiling for the truncation degree. If convergence is not reached by this degree, the search terminates and returns this value.
power_of_two – If True, rounds the resulting kmax up to the nearest power of two (useful for FFT optimizations). Defaults to False.

Returns:

A fully instantiated L² space on the circle with the optimal kmax.

Return type:

Lebesgue

Raises:

RuntimeError – If the energy sequence fails to converge within 100,000 iterations and no max_degree is specified.

classmethod from_heat_kernel_prior(kernel_scale: float, /, *, radius: float = 1.0, rtol: float = 1e-06, min_degree: int = 0, max_degree: int | None = None, power_of_two: bool = False) → Lebesgue[source]¶

Creates an instance of the L² space on the circle, tuned to the expected energy of a Heat Kernel prior measure.

Parameters:

kernel_scale – The length-scale parameter of the heat kernel covariance.
radius – The radius of the circle. Defaults to 1.0.
rtol – The relative tolerance for the energy truncation. Defaults to 1e-6.
min_degree – The absolute minimum truncation degree to return. Defaults to 0.
max_degree – An optional safety ceiling for the truncation degree.
power_of_two – If True, rounds the resulting kmax up to the nearest power of two.

classmethod from_sobolev_kernel_prior(kernel_order: float, kernel_scale: float, /, *, radius: float = 1.0, rtol: float = 1e-06, min_degree: int = 0, max_degree: int | None = None, power_of_two: bool = False) → Lebesgue[source]¶

Creates an instance of the L² space on the circle, tuned to the expected energy of a Sobolev-type prior measure.

Parameters:

kernel_order – The smoothness order of the Sobolev prior measure.
kernel_scale – The length-scale parameter of the Sobolev prior measure.
radius – The radius of the circle. Defaults to 1.0.
rtol – The relative tolerance for the energy truncation. Defaults to 1e-6.
min_degree – The absolute minimum truncation degree to return. Defaults to 0.
max_degree – An optional safety ceiling for the truncation degree.
power_of_two – If True, rounds the resulting kmax up to the nearest power of two.

property gaussian_curvature: float¶: Returns the Gaussian curvature (K) of the symmetric manifold.

geodesic_distance(p1: float, p2: float) → float[source]¶

Returns the shortest distance along the manifold between two points.

Parameters:

p1 – The starting point.
p2 – The end point.

geodesic_quadrature(p1: float, p2: float, n_points: int) → Tuple[List[float], ndarray][source]¶

Returns quadrature points and weights for a geodesic between p1 and p2.

Returns:: List of manifold coordinates. weights: Integration weights scaled by the line element.
Return type:: points

index_to_integer(k: int) → int[source]¶: Maps an eigenvalue index to an integer.

integer_to_index(i: int) → int[source]¶: Maps an integer to the eigenvalue index

invariant_covariance_function(spectral_variances: ndarray) → Callable[[ndarray], ndarray][source]¶: Returns a vectorized function that evaluates the two-point covariance of an invariant measure as a function of geodesic distance.

property inverse_fft_factor: float¶: The factor by which the Fourier coefficients are scaled in inverse transformations.

is_element(x: Any) → bool[source]¶

Checks if an object is a valid element of the space.

Note: The default implementation checks the object’s type against the type of the zero vector. This may not be robust for all vector representations and can be overridden if needed.

Parameters:: x – The object to check.
Returns:: True if the object is an element of the space, False otherwise.

property kmax: int¶: The maximum Fourier degree represented in this space.

laplacian_eigenvalue(k: int) → float[source]¶

Returns the eigenvalue of the Laplacian for a given index.

The index k can be a single integer (e.g., for a circle) or a tuple of integers (e.g., for a sphere or torus), depending on the geometry of the space.

Parameters:: k – The index of the eigenvalue to return.

laplacian_eigenvector_squared_norm(k: int) → float[source]¶

Returns the squared norm of the eigenvalue of the Laplacian for a given index.

The index k can be a single integer (e.g., for a circle) or a tuple of integers (e.g., for a sphere or torus), depending on the geometry of the space.

Parameters:: k – The index of the eigenvalue to return.

laplacian_eigenvectors_at_point(theta: float) → ndarray[source]¶

Returns a list of the values of the eigenvectors at a given point

Parameters:: x – The evaluation point.

property point_spacing: float¶: The angular spacing between grid points.

points() → ndarray[source]¶: Returns a numpy array of the grid point angles.

project_function(f: Callable[[float], float]) → ndarray[source]¶

Returns an element of the space by projecting a given function.

The function f is evaluated at the grid points of the space.

Parameters:: f – A function that takes an angle (float) and returns a value.

property radius: float¶: The radius of the circle.

random_point() → float[source]¶: Returns a single random point from the underlying symmetric space.

representative_index(degree: int) → int[source]¶: Returns a valid eigenfunction index for the given degree.

to_coefficients(u: ndarray) → ndarray[source]¶: Maps a function vector to its complex Fourier coefficients.

to_components(x: ndarray) → ndarray[source]¶

Maps a vector to its representation as a NumPy component array.

Parameters:: x – A vector in the space.
Returns:: The components of the vector as a NumPy array.

vector_multiply(x1: ndarray, x2: ndarray) → ndarray[source]¶: Computes the pointwise product of two vectors.

vector_sqrt(x: ndarray) → ndarray[source]¶: Returns the pointwise square root of a vector.

with_degree(degree: int) → Lebesgue[source]¶: Returns a new instance of the space with a modified truncation degree.

class pygeoinf.symmetric_space.circle.Sobolev(kmax: int, order: float, scale: float, /, *, radius: float = 1.0, safe: bool = True)[source]¶

Bases: SymmetricSobolevSpace

Implementation of the Sobolev space Hˢ on a circle.

angles() → ndarray[source]¶: Returns a numpy array of the grid point angles.

property derivative_operator: LinearOperator¶: Returns the derivative operator from the space to one with a lower order.

property fft_factor: float¶: The factor by which the Fourier coefficients are scaled in forward transformations.

from_coefficients(coeff: ndarray) → ndarray[source]¶: Maps complex Fourier coefficients to a function vector.

classmethod from_covariance(covariance_function: Callable[[float], float], order: float, scale: float, /, *, radius: float = 1.0, rtol: float = 1e-06, min_degree: int = 0, max_degree: int | None = None, power_of_two: bool = False, safe: bool = True) → Sobolev[source]¶

Creates an instance of the Sobolev space with a Fourier truncation degree (kmax) automatically chosen to capture the expected energy of functions drawn from a specified prior measure.

This factory method calculates the expected squared norm (energy) of a random field whose spectral variances are defined by the provided covariance_function, accounting for the Sobolev mass-weighting factor. It iteratively adds higher frequency modes until the relative contribution of the next degree drops below the specified tolerance.

Parameters:

covariance_function – A callable mapping a Laplacian eigenvalue to its spectral variance.
order – The Sobolev order, controlling the smoothness of functions.
scale – The Sobolev length-scale.
radius – The radius of the circle. Defaults to 1.0.
rtol – The relative tolerance for the energy truncation. Defaults to 1e-6.
min_degree – The absolute minimum truncation degree to return. Defaults to 0.
max_degree – An optional safety ceiling for the truncation degree. If convergence is not reached by this degree, the search terminates and returns this value.
power_of_two – If True, rounds the resulting kmax up to the nearest power of two (useful for FFT optimizations). Defaults to False.
safe – If True, enables mathematical correctness checks during operations.

Returns:

A fully instantiated Sobolev space on the circle with the optimal kmax.

Return type:

Sobolev

Raises:

RuntimeError – If the energy sequence fails to converge within 100,000 iterations and no max_degree is specified.

classmethod from_heat_kernel_prior(kernel_scale: float, order: float, scale: float, /, *, radius: float = 1.0, rtol: float = 1e-06, min_degree: int = 0, max_degree: int | None = None, power_of_two: bool = False, safe: bool = True) → Sobolev[source]¶

Creates an instance of the Sobolev space on the circle, tuned to the expected energy of a Heat Kernel prior measure.

Parameters:

kernel_scale – The length-scale parameter of the heat kernel covariance.
order – The Sobolev order defining the function space.
scale – The Sobolev length-scale defining the function space.
radius – The radius of the circle. Defaults to 1.0.
rtol – The relative tolerance for the energy truncation. Defaults to 1e-6.
min_degree – The absolute minimum truncation degree to return. Defaults to 0.
max_degree – An optional safety ceiling for the truncation degree.
power_of_two – If True, rounds the resulting kmax up to the nearest power of two.
safe – If True, enables mathematical correctness checks during operations.

classmethod from_sobolev_kernel_prior(kernel_order: float, kernel_scale: float, order: float, scale: float, /, *, radius: float = 1.0, rtol: float = 1e-06, min_degree: int = 0, max_degree: int | None = None, power_of_two: bool = False, safe: bool = True) → Sobolev[source]¶

Creates an instance of the Sobolev space on the circle, tuned to the expected energy of a Sobolev-type prior measure.

Parameters:

kernel_order – The smoothness order of the Sobolev prior measure.
kernel_scale – The length-scale parameter of the Sobolev prior measure.
order – The Sobolev order defining the function space.
scale – The Sobolev length-scale defining the function space.
radius – The radius of the circle. Defaults to 1.0.
rtol – The relative tolerance for the energy truncation. Defaults to 1e-6.
min_degree – The absolute minimum truncation degree to return. Defaults to 0.
max_degree – An optional safety ceiling for the truncation degree.
power_of_two – If True, rounds the resulting kmax up to the nearest power of two.
safe – If True, enables mathematical correctness checks during operations.

static from_sobolev_parameters(order: float, scale: float, /, *, radius: float = 1.0, rtol: float = 1e-06, power_of_two: bool = False, safe: bool = True) → Sobolev[source]¶

Creates an instance with kmax chosen based on Sobolev parameters.

The method estimates the truncation error for the Dirac measure and is only applicable for spaces with order > 0.5.

Parameters:

order – The Sobolev order. Must be > 0.5.
scale – The Sobolev length-scale.
radius – The radius of the circle. Defaults to 1.0.
rtol – Relative tolerance used in assessing truncation error. Defaults to 1e-8.
power_of_two – If True, kmax is set to the next power of two.
safe – If true, the class checks for mathematical correctness of operations where possible.

Returns:

An instance of the Sobolev class with an appropriate kmax.

Raises:

ValueError – If order is <= 0.5.

property inverse_fft_factor: float¶: The factor by which the Fourier coefficients are scaled in inverse transformations.

property kmax: int¶: The maximum Fourier degree represented in this space.

property point_spacing: float¶: The angular spacing between grid points.

points() → ndarray[source]¶: Returns a numpy array of the grid point angles.

project_function(f: Callable[[float], float]) → ndarray[source]¶

Returns an element of the space by projecting a given function.

The function f is evaluated at the grid points of the space.

Parameters:: f – A function that takes an angle (float) and returns a value.

property radius: float¶: The radius of the circle.

to_coefficients(u: ndarray) → ndarray[source]¶: Maps a function vector to its complex Fourier coefficients.

with_degree(degree: int) → Sobolev[source]¶: Returns a new instance of the space with a modified truncation degree.

with_order(order: float) → Sobolev[source]¶: Returns a new instance of the exact same space but with a modified Sobolev order.

pygeoinf.symmetric_space.circle.plot(space: Lebesgue | Sobolev, u: ndarray, /, *, ax: Axes | None = None, **kwargs) → Axes[source]¶

Creates a simple line plot of a function on the circle.

Parameters:

space – The function space.
u – A 1D numpy array representing the function values (the y-axis).
ax – An existing Matplotlib Axes object. If None, plots to the current active axes.
**kwargs – Additional keyword arguments forwarded directly to ax.plot().

Returns:

The Matplotlib Axes object.

pygeoinf.symmetric_space.circle.plot_error_bounds(space: Lebesgue | Sobolev, u: ndarray, u_bound: ndarray, /, *, ax: Axes | None = None, **kwargs) → Axes[source]¶

Plots a function on the circle along with its pointwise error bounds.

This is particularly useful for visualizing Gaussian measures or Bayesian posterior uncertainties over the circular domain.

Parameters:

space – The function space.
u – A 1D numpy array representing the mean function values.
u_bound – A 1D numpy array giving the pointwise standard deviations or bounds.
ax – An existing Matplotlib Axes object. If None, plots to the current active axes.
**kwargs – Additional keyword arguments forwarded directly to ax.fill_between().

Returns:

The Matplotlib Axes object.

pygeoinf.symmetric_space.line module¶

class pygeoinf.symmetric_space.line.Lebesgue(kmax: int, /, *, a: float = 0.0, b: float = 1.0, c: float = 0.1)[source]¶

Bases: AbstractSymmetricLebesgueSpace

Implementation of the Lebesgue space L² on a line of interval. In the case of the line, it is assumed that the function has support limited to the chosen interval.

A function is represented by its values on an evenly spaced grid. The co-ordinate basis for the space is through Fourier expansions. Details of the implementation are handled by mapping the function to one defined on a circle.

property a: float¶: Returns the left end point.

angle_to_point(th: float) → float[source]¶: Maps an angle to the corresponding point.

property b: float¶: Returns the right end point.

property c: float¶: Returns the padding distance.

property circle_space¶: Returns the isomorphic space of functions on a circle.

degree_multiplicity(degree: int) → int[source]¶: Returns the number of eigenfunctions associated with a given degree.

degree_transfer_operator(target_degree: int) → LinearOperator[source]¶: Returns the transfer operator from this space to one with a different degree.

from_coefficients(coeff: ndarray) → ndarray[source]¶: Maps complex Fourier coefficients to a function vector.

from_components(c: ndarray) → ndarray[source]¶

Maps a NumPy component array back to a vector in the space.

Parameters:: c – The components of the vector as a NumPy array.
Returns:: The corresponding vector in the space.

classmethod from_covariance(covariance_function: Callable[[float], float], /, *, a: float = 0.0, b: float = 1.0, c: float = 0.1, rtol: float = 1e-06, min_degree: int = 0, max_degree: int | None = None, power_of_two: bool = False) → Lebesgue[source]¶

Creates an instance of the L² space with a Fourier truncation degree (kmax) automatically chosen to capture the expected energy of functions drawn from a specified prior measure.

Parameters:

covariance_function – A callable mapping a Laplacian eigenvalue to its spectral variance.
a – The left endpoint of the interval. Defaults to 0.0.
b – The right endpoint of the interval. Defaults to 1.0.
c – The padding distance for the computational domain. Defaults to 0.1.
rtol – The relative tolerance for the energy truncation. Defaults to 1e-6.
min_degree – The absolute minimum truncation degree to return. Defaults to 0.
max_degree – An optional safety ceiling for the truncation degree. If convergence is not reached by this degree, the search terminates and returns this value.
power_of_two – If True, rounds the resulting kmax up to the nearest power of two (useful for FFT optimizations). Defaults to False.

Returns:

A fully instantiated L² space on the line with the optimal kmax.

Return type:

Lebesgue

Raises:

RuntimeError – If the energy sequence fails to converge within 100,000 iterations and no max_degree is specified.

classmethod from_heat_kernel_prior(kernel_scale: float, /, *, a: float = 0.0, b: float = 1.0, c: float = 0.1, rtol: float = 1e-06, min_degree: int = 0, max_degree: int | None = None, power_of_two: bool = False) → Lebesgue[source]¶

Creates an instance of the L² space on the line, tuned to the expected energy of a Heat Kernel prior measure.

Parameters:

kernel_scale – The length-scale parameter of the heat kernel covariance.
a – The left endpoint of the interval. Defaults to 0.0.
b – The right endpoint of the interval. Defaults to 1.0.
c – The padding distance for the computational domain. Defaults to 0.1.
rtol – The relative tolerance for the energy truncation. Defaults to 1e-6.
min_degree – The absolute minimum truncation degree to return. Defaults to 0.
max_degree – An optional safety ceiling for the truncation degree.
power_of_two – If True, rounds the resulting kmax up to the nearest power of two.

classmethod from_sobolev_kernel_prior(kernel_order: float, kernel_scale: float, /, *, a: float = 0.0, b: float = 1.0, c: float = 0.1, rtol: float = 1e-06, min_degree: int = 0, max_degree: int | None = None, power_of_two: bool = False) → Lebesgue[source]¶

Creates an instance of the L² space on the line, tuned to the expected energy of a Sobolev-type prior measure.

Parameters:

kernel_order – The smoothness order of the Sobolev prior measure.
kernel_scale – The length-scale parameter of the Sobolev prior measure.
a – The left endpoint of the interval. Defaults to 0.0.
b – The right endpoint of the interval. Defaults to 1.0.
c – The padding distance for the computational domain. Defaults to 0.1.
rtol – The relative tolerance for the energy truncation. Defaults to 1e-6.
min_degree – The absolute minimum truncation degree to return. Defaults to 0.
max_degree – An optional safety ceiling for the truncation degree.
power_of_two – If True, rounds the resulting kmax up to the nearest power of two.

property gaussian_curvature: float¶: Returns the Gaussian curvature (K) of the symmetric manifold.

geodesic_distance(p1: float, p2: float) → float[source]¶

Returns the shortest distance along the manifold between two points.

Parameters:

p1 – The starting point.
p2 – The end point.

geodesic_quadrature(p1: float, p2: float, n_points: int) → Tuple[List[float], ndarray][source]¶

Returns quadrature points and weights for a geodesic between p1 and p2.

Returns:: List of manifold coordinates. weights: Integration weights scaled by the line element.
Return type:: points

index_to_integer(k: int) → int[source]¶: Maps an eigenvalue index to an integer.

integer_to_index(i: int) → int[source]¶: Maps an integer to the eigenvalue index

invariant_covariance_function(spectral_variances: ndarray) → Callable[[ndarray], ndarray][source]¶: Returns a vectorized function that evaluates the two-point covariance of an invariant measure as a function of geodesic distance.

is_element(x: Any) → bool[source]¶

Checks if an object is a valid element of the space.

Note: The default implementation checks the object’s type against the type of the zero vector. This may not be robust for all vector representations and can be overridden if needed.

Parameters:: x – The object to check.
Returns:: True if the object is an element of the space, False otherwise.

property kmax: int¶: Returns the maximum Fourier degree.

laplacian_eigenvalue(k: int) → float[source]¶

Returns the eigenvalue of the Laplacian for a given index.

The index k can be a single integer (e.g., for a circle) or a tuple of integers (e.g., for a sphere or torus), depending on the geometry of the space.

Parameters:: k – The index of the eigenvalue to return.

laplacian_eigenvector_squared_norm(k: int) → float[source]¶

Returns the squared norm of the eigenvalue of the Laplacian for a given index.

The index k can be a single integer (e.g., for a circle) or a tuple of integers (e.g., for a sphere or torus), depending on the geometry of the space.

Parameters:: k – The index of the eigenvalue to return.

laplacian_eigenvectors_at_point(x: float) → ndarray[source]¶

Returns a list of the values of the eigenvectors at a given point

Parameters:: x – The evaluation point.

point_to_angle(x: float) → float[source]¶: Maps a point to the corresponding angle.

points() → ndarray[source]¶: Returns a numpy array of the grid points.

project_function(f: Callable[[float], float]) → ndarray[source]¶

Returns an element of the space by projecting a given function.

The function f is evaluated at the grid points of the space. To prevent spectral ringing (Gibbs phenomenon) from non-periodic functions, the function is smoothly tapered to zero within the padding regions [a-c, a] and [b, b+c] using a raised cosine window.

Parameters:: f – A function that takes a point (float) and returns a value.

random_point() → float[source]¶: Returns a single random point from the underlying symmetric space.

representative_index(degree: int) → int[source]¶: Returns a valid eigenfunction index for the given degree.

to_coefficients(u: ndarray) → ndarray[source]¶: Maps a function vector to its complex Fourier coefficients.

to_components(x: ndarray) → ndarray[source]¶

Maps a vector to its representation as a NumPy component array.

Parameters:: x – A vector in the space.
Returns:: The components of the vector as a NumPy array.

vector_multiply(x1: ndarray, x2: ndarray) → ndarray[source]¶: Computes the pointwise product of two vectors.

vector_sqrt(x: ndarray) → ndarray[source]¶: Returns the pointwise square root of a vector.

with_degree(degree: int) → Lebesgue[source]¶: Returns a new instance of the space with a modified truncation degree.

class pygeoinf.symmetric_space.line.Sobolev(kmax: int, order: float, scale: float, /, *, a: float = 0.0, b: float = 1.0, c: float = None, safe: bool = True)[source]¶

Bases: SymmetricSobolevSpace

Implementation of the Sobolev space Hˢ on a circle.

property a: float¶: Returns the left end point.

angle_to_point(th: float) → float[source]¶: Maps an angle to the corresponding point.

property b: float¶: Returns the right end point.

property c: float¶: Returns the padding distance.

property circle_space¶: Returns the isomorphic space of functions on a circle.

property derivative_operator¶: Returns the derivative operator from the space to one with a lower order.

from_coefficients(coeff: ndarray) → ndarray[source]¶: Maps complex Fourier coefficients to a function vector.

classmethod from_covariance(covariance_function: Callable[[float], float], order: float, scale: float, /, *, a: float = 0.0, b: float = 1.0, c: float | None = None, rtol: float = 1e-06, min_degree: int = 0, max_degree: int | None = None, power_of_two: bool = False, safe: bool = True) → Sobolev[source]¶

Creates an instance of the Sobolev space with a Fourier truncation degree (kmax) automatically chosen to capture the expected energy of functions drawn from a specified prior measure.

Parameters:

covariance_function – A callable mapping a Laplacian eigenvalue to its spectral variance.
order – The Sobolev order, controlling the smoothness of functions.
scale – The Sobolev length-scale.
a – The left endpoint of the interval. Defaults to 0.0.
b – The right endpoint of the interval. Defaults to 1.0.
c – The padding distance for the computational domain. Defaults to 6 * scale.
rtol – The relative tolerance for the energy truncation. Defaults to 1e-6.
min_degree – The absolute minimum truncation degree to return. Defaults to 0.
max_degree – An optional safety ceiling for the truncation degree. If convergence is not reached by this degree, the search terminates and returns this value.
power_of_two – If True, rounds the resulting kmax up to the nearest power of two (useful for FFT optimizations). Defaults to False.
safe – If True, enables mathematical correctness checks during operations.

Returns:

A fully instantiated Sobolev space on the line with the optimal kmax.

Return type:

Sobolev

Raises:

RuntimeError – If the energy sequence fails to converge within 100,000 iterations and no max_degree is specified.

classmethod from_heat_kernel_prior(kernel_scale: float, order: float, scale: float, /, *, a: float = 0.0, b: float = 1.0, c: float | None = None, rtol: float = 1e-06, min_degree: int = 0, max_degree: int | None = None, power_of_two: bool = False, safe: bool = True) → Sobolev[source]¶

Creates an instance of the Sobolev space on the line, tuned to the expected energy of a Heat Kernel prior measure.

Parameters:

kernel_scale – The length-scale parameter of the heat kernel covariance.
order – The Sobolev order defining the function space.
scale – The Sobolev length-scale defining the function space.
a – The left endpoint of the interval. Defaults to 0.0.
b – The right endpoint of the interval. Defaults to 1.0.
c – The padding distance for the computational domain. Defaults to 6 * scale.
rtol – The relative tolerance for the energy truncation. Defaults to 1e-6.
min_degree – The absolute minimum truncation degree to return. Defaults to 0.
max_degree – An optional safety ceiling for the truncation degree.
power_of_two – If True, rounds the resulting kmax up to the nearest power of two.
safe – If True, enables mathematical correctness checks during operations.

classmethod from_sobolev_kernel_prior(kernel_order: float, kernel_scale: float, order: float, scale: float, /, *, a: float = 0.0, b: float = 1.0, c: float | None = None, rtol: float = 1e-06, min_degree: int = 0, max_degree: int | None = None, power_of_two: bool = False, safe: bool = True) → Sobolev[source]¶

Creates an instance of the Sobolev space on the line, tuned to the expected energy of a Sobolev-type prior measure.

Parameters:

kernel_order – The smoothness order of the Sobolev prior measure.
kernel_scale – The length-scale parameter of the Sobolev prior measure.
order – The Sobolev order defining the function space.
scale – The Sobolev length-scale defining the function space.
a – The left endpoint of the interval. Defaults to 0.0.
b – The right endpoint of the interval. Defaults to 1.0.
c – The padding distance for the computational domain. Defaults to 6 * scale.
rtol – The relative tolerance for the energy truncation. Defaults to 1e-6.
min_degree – The absolute minimum truncation degree to return. Defaults to 0.
max_degree – An optional safety ceiling for the truncation degree.
power_of_two – If True, rounds the resulting kmax up to the nearest power of two.
safe – If True, enables mathematical correctness checks during operations.

static from_sobolev_parameters(order: float, scale: float, /, *, a: float = 0.0, b: float = 1.0, c: float = None, rtol: float = 1e-06, power_of_two: bool = False, safe: bool = True) → Sobolev[source]¶

Creates an instance with kmax chosen based on Sobolev parameters.

The method estimates the truncation error for the Dirac measure and is only applicable for spaces with order > 0.5.

Parameters:

order – The Sobolev order. Must be > 0.5.
scale – The Sobolev length-scale.
a – The left endpoint of the interval. Defaults to 0.0
b – The right endpoint of the interval. Dafaults to 0.1
c – The padding distance for the computational domain. Defaults to 6*scale
rtol – Relative tolerance used in assessing truncation error. Defaults to 1e-8.
power_of_two – If True, kmax is set to the next power of two.
safe – If true, the class checks for mathematical correctness of operations where possible.

Returns:

An instance of the Sobolev class with an appropriate kmax.

Raises:

ValueError – If order is <= 0.5.

property kmax: int¶: Returns the maximum Fourier degree.

point_to_angle(x: float) → float[source]¶: Maps a point to the corresponding angle.

points() → ndarray[source]¶: Returns a numpy array of the grid points.

project_function(f: Callable[[float], float]) → ndarray[source]¶

Returns an element of the space by projecting a given function.

The function f is evaluated at the grid points of the space.

Parameters:: f – A function that takes an angle (float) and returns a value.

to_coefficients(u: ndarray) → ndarray[source]¶: Maps a function vector to its complex Fourier coefficients.

with_degree(degree: int) → Sobolev[source]¶: Returns a new instance of the space with a modified truncation degree.

with_order(order: float) → Sobolev[source]¶: Returns a new instance of the exact same space but with a modified Sobolev order.

pygeoinf.symmetric_space.line.plot(space: Lebesgue | Sobolev, u: ndarray, /, *, ax: Axes | None = None, full: bool | None = False, **kwargs) → Axes[source]¶

Makes a simple plot of a function on the line/interval.

Parameters:

space – The function space.
u – The vector representing the function to be plotted.
ax – An existing Matplotlib Axes object. If None, plots to the current active axes.
full – If False, limits the x-axis to the interval [a, b]. Defaults to False.
**kwargs – Keyword arguments forwarded to ax.plot().

Returns:

The Matplotlib Axes object.

pygeoinf.symmetric_space.line.plot_error_bounds(space: Lebesgue | Sobolev, u: ndarray, u_bound: ndarray, /, *, ax: Axes | None = None, full: bool | None = False, **kwargs) → Axes[source]¶

Plots a function on the line along with its pointwise error bounds.

Parameters:

space – The function space.
u – A 1D numpy array representing the mean function values.
u_bound – A 1D numpy array giving the pointwise standard deviations or bounds.
ax – An existing Matplotlib Axes object. If None, plots to the current active axes.
full – If False, limits the x-axis to the interval [a, b]. Defaults to False.
**kwargs – Additional keyword arguments forwarded directly to ax.fill_between().

Returns:

The Matplotlib Axes object.

pygeoinf.symmetric_space.plane module¶

class pygeoinf.symmetric_space.plane.Lebesgue(kmax: int, /, *, ax: float = 0.0, bx: float = 1.0, cx: float = 0.1, ay: float = 0.0, by: float = 1.0, cy: float = 0.1)[source]¶

Bases: AbstractSymmetricLebesgueSpace

Implementation of the Lebesgue space L² on a compact rectangle in the 2D plane.

Functions are assumed to have support limited to the chosen rectangular interval [ax, bx] x [ay, by]. The coordinate basis for the space is through 2D Fourier expansions. Details of the implementation are handled by smoothly embedding the function into a larger 2D periodic Torus.

angle_to_point(th: Tuple[float, float]) → Tuple[float, float][source]¶

angle_to_point_x(th_x: float) → float[source]¶

angle_to_point_y(th_y: float) → float[source]¶

property bounds_x: Tuple[float, float, float]¶: Returns (ax, bx, cx).

property bounds_y: Tuple[float, float, float]¶: Returns (ay, by, cy).

degree_multiplicity(degree: int) → int[source]¶: Returns the number of eigenfunctions associated with a given degree.

degree_transfer_operator(target_degree: int) → LinearOperator[source]¶: Returns the transfer operator from this space to one with a different degree.

estimate_truncation_degree(*args, **kwargs) → int[source]¶: Estimates the truncation degree required to capture the expected energy of a function drawn from a prior measure.

from_coefficient_operator(kmax: int) → LinearOperator[source]¶

from_coefficients(coeff: ndarray) → ndarray[source]¶

from_components(c: ndarray) → ndarray[source]¶

Maps a NumPy component array back to a vector in the space.

Parameters:: c – The components of the vector as a NumPy array.
Returns:: The corresponding vector in the space.

classmethod from_covariance(covariance_function: Callable[[float], float], /, *, ax: float = 0.0, bx: float = 1.0, cx: float = 0.1, ay: float = 0.0, by: float = 1.0, cy: float = 0.1, rtol: float = 1e-06, min_degree: int = 1, max_degree: int | None = None, power_of_two: bool = False) → Lebesgue[source]¶

classmethod from_heat_kernel_prior(kernel_scale: float, /, **kwargs) → Lebesgue[source]¶

classmethod from_sobolev_kernel_prior(kernel_order: float, kernel_scale: float, /, **kwargs) → Lebesgue[source]¶

property gaussian_curvature: float¶: Returns the Gaussian curvature (K) of the symmetric manifold.

geodesic_distance(p1: Tuple[float, float], p2: Tuple[float, float]) → float[source]¶: Euclidean distance on the plane (ignores Torus wrapping).

geodesic_quadrature(p1: Tuple[float, float], p2: Tuple[float, float], n_points: int) → Tuple[List[Tuple[float, float]], ndarray][source]¶: Standard straight-line Euclidean quadrature.

index_to_integer(k: int) → int[source]¶: Maps an eigenvalue index to an integer.

integer_to_index(i: int) → int[source]¶: Maps an integer to the eigenvalue index

invariant_covariance_function(spectral_variances: ndarray) → Callable[[ndarray], ndarray][source]¶: Returns a vectorized function that evaluates the two-point covariance of an invariant measure as a function of geodesic distance.

is_element(x: Any) → bool[source]¶

Checks if an object is a valid element of the space.

Note: The default implementation checks the object’s type against the type of the zero vector. This may not be robust for all vector representations and can be overridden if needed.

Parameters:: x – The object to check.
Returns:: True if the object is an element of the space, False otherwise.

property kmax: int¶

laplacian_eigenvalue(k: int) → float[source]¶

Returns the eigenvalue of the Laplacian for a given index.

The index k can be a single integer (e.g., for a circle) or a tuple of integers (e.g., for a sphere or torus), depending on the geometry of the space.

Parameters:: k – The index of the eigenvalue to return.

laplacian_eigenvector_squared_norm(k: int) → float[source]¶

Returns the squared norm of the eigenvalue of the Laplacian for a given index.

The index k can be a single integer (e.g., for a circle) or a tuple of integers (e.g., for a sphere or torus), depending on the geometry of the space.

Parameters:: k – The index of the eigenvalue to return.

laplacian_eigenvectors_at_point(p: Tuple[float, float]) → ndarray[source]¶

Returns a list of the values of the eigenvectors at a given point

Parameters:: x – The evaluation point.

point_to_angle(p: Tuple[float, float]) → Tuple[float, float][source]¶

points() → Tuple[ndarray, ndarray][source]¶: Returns two 1D arrays of the flattened X and Y coordinates.

project_function(f: Callable[[Tuple[float, float]], float]) → ndarray[source]¶

Returns an element of the space by projecting a given 2D function.

Applies a separable 2D raised-cosine window to smoothly taper the function to zero within the padding regions.

random_point() → Tuple[float, float][source]¶: Returns a random point strictly within the unpadded bounds.

representative_index(degree: int) → int[source]¶: Returns a valid eigenfunction index for the given degree.

spectral_projection_operator(modes: List[Tuple[int, int]]) → LinearOperator[source]¶

to_coefficient_operator(kmax: int) → LinearOperator[source]¶

to_coefficients(u: ndarray) → ndarray[source]¶

to_components(x: ndarray) → ndarray[source]¶

Maps a vector to its representation as a NumPy component array.

Parameters:: x – A vector in the space.
Returns:: The components of the vector as a NumPy array.

property torus_space: Lebesgue¶: Returns the isomorphic space of functions on the extended Torus.

vector_multiply(x1: ndarray, x2: ndarray) → ndarray[source]¶: Computes the pointwise product of two vectors.

vector_sqrt(x: ndarray) → ndarray[source]¶: Returns the pointwise square root of a vector.

wavevector_indices(kx: int, ky: int) → List[int][source]¶

with_degree(degree: int) → Lebesgue[source]¶: Returns a new instance of the space with a modified truncation degree.

class pygeoinf.symmetric_space.plane.Sobolev(kmax: int, order: float, scale: float, /, *, ax: float = 0.0, bx: float = 1.0, cx: float | None = None, ay: float = 0.0, by: float = 1.0, cy: float | None = None, safe: bool = True)[source]¶

Bases: SymmetricSobolevSpace

Implementation of the Sobolev space Hˢ on a compact 2D plane.

property bounds_x: Tuple[float, float, float]¶

property bounds_y: Tuple[float, float, float]¶

estimate_truncation_degree(covariance_function: Callable[[float], float], /, *, rtol: float = 1e-06, min_degree: int = 1, max_degree: int | None = None) → int[source]¶: Delegates the energy truncation search to the underlying Plane Lebesgue space, ensuring it loops geometrically over (kx, ky) shells rather than flat 1D indices.

from_coefficient_operator(kmax: int) → LinearOperator[source]¶

classmethod from_covariance(covariance_function: Callable[[float], float], order: float, scale: float, /, **kwargs) → Sobolev[source]¶

classmethod from_heat_kernel_prior(kernel_scale: float, order: float, scale: float, /, **kwargs) → Sobolev[source]¶

classmethod from_sobolev_kernel_prior(kernel_order: float, kernel_scale: float, order: float, scale: float, /, **kwargs) → Sobolev[source]¶

static from_sobolev_parameters(order: float, scale: float, /, **kwargs) → Sobolev[source]¶

property kmax: int¶

points() → Tuple[ndarray, ndarray][source]¶

project_function(f: Callable[[Tuple[float, float]], float]) → ndarray[source]¶

Returns an element of the space by projecting a given function.

Parameters:: f – A function that takes a point and returns a value.

spectral_projection_operator(modes: List[Tuple[int, int]]) → LinearOperator[source]¶

to_coefficient_operator(kmax: int) → LinearOperator[source]¶

property torus_space: Sobolev¶

with_degree(degree: int) → Sobolev[source]¶: Returns a new instance of the space with a modified truncation degree.

with_order(order: float) → Sobolev[source]¶: Returns a new instance of the exact same space but with a modified Sobolev order.

pygeoinf.symmetric_space.plane.plot(space: Lebesgue | Sobolev, u: ndarray, /, *, ax: Axes | Tuple[Axes, Any] | None = None, contour: bool = False, cmap: str = 'RdBu', symmetric: bool | float = False, contour_lines: bool = False, contour_lines_kwargs: dict | None = None, num_levels: int = 10, colorbar: bool = False, colorbar_kwargs: dict | None = None, full: bool = False, **kwargs) → Tuple[Axes, Any][source]¶

Creates a high-quality map plot of a function on the 2D Plane.

By default, this crops the view to the primary computational domain [ax, bx] x [ay, by]. Set full=True to view the tapered padding regions.

Parameters:

space – The Plane function space (Lebesgue or Sobolev).
u – The scalar field array to be plotted.
ax – An existing Matplotlib Axes object or (Axes, Artist) tuple.
contour – If True, renders the field as a filled contour plot (contourf).
cmap – The Matplotlib colormap string or object to use.
symmetric – If True, dynamically centers the color scale around zero.
contour_lines – If True, overlays solid contour lines on top of the base plot.
contour_lines_kwargs – A dictionary of keyword arguments passed to ax.contour.
num_levels – The number of color levels to generate automatically. Defaults to 10.
colorbar – If True, attaches a colorbar to the plot.
colorbar_kwargs – Formatting options for the colorbar.
full – If False, crops the plot to the unpadded domain. If True, shows the padding.
**kwargs – Additional keyword arguments forwarded to ax.contourf or ax.pcolormesh.

Returns:

A tuple (Axes, Any) containing the Axes and the rendered image artist.

pygeoinf.symmetric_space.plane.plot_geodesic(p1: Tuple[float, float], p2: Tuple[float, float], /, *, ax: Axes | Tuple[Axes, Any] | None = None, **kwargs) → Tuple[Axes, Any][source]¶

Plots a straight-line geodesic between two points on the 2D plane.

Parameters:

p1 – The starting coordinate as a tuple (x, y).
p2 – The ending coordinate as a tuple (x, y).
ax – An existing Matplotlib Axes object or (Axes, Artist) tuple.
**kwargs – Keyword arguments passed directly to ax.plot.

pygeoinf.symmetric_space.plane.plot_geodesic_network(paths: List[Tuple[Tuple[float, float], Tuple[float, float]]], /, *, ax: Axes | Tuple[Axes, Any] | None = None, plot_sources: bool = True, plot_receivers: bool = True, source_kwargs: dict | None = None, receiver_kwargs: dict | None = None, **kwargs) → Tuple[Axes, List[Any]][source]¶

Plots a network of intersecting straight-line paths on the 2D plane.

Parameters:

paths – A list of point pairs defining the network. Each pair is ((x1, y1), (x2, y2)).
ax – An existing Matplotlib Axes object or (Axes, Artist) tuple.
plot_sources – If True, plots the unique source points.
plot_receivers – If True, plots the unique receiver points.
source_kwargs – Styling arguments for the source scatter points.
receiver_kwargs – Styling arguments for the receiver scatter points.
**kwargs – Default styling arguments applied to the geodesic lines.

pygeoinf.symmetric_space.plane.plot_points(points: List[Tuple[float, float]], /, *, data: List[float] | ndarray | None = None, ax: Axes | Tuple[Axes, Any] | None = None, cmap: str = 'RdBu', color: str = 'red', s: float = 20, marker: str = 'o', edgecolors: str = 'none', zorder: int = 5, symmetric: bool | float = False, colorbar: bool = False, colorbar_kwargs: dict | None = None, **kwargs) → Tuple[Axes, Any][source]¶

Plots discrete observation points on the 2D plane.

Parameters:

points – A list of (x, y) tuples.
data – Optional array of values to color the points by.
ax – An existing Matplotlib Axes object or (Axes, Artist) tuple.
cmap – The colormap to use if data is provided.
color – The fixed color to use if data is NOT provided.
s – The marker size.
marker – The marker style (e.g., ‘o’, ‘*’, ‘^’).
edgecolors – The color of the marker edges.
zorder – The drawing order.
symmetric – If True, dynamically centers the color scale symmetrically around zero.
colorbar – If True and data is provided, attaches a colorbar.
colorbar_kwargs – Formatting options for the colorbar.
**kwargs – Additional keyword arguments passed to ax.scatter.

pygeoinf.symmetric_space.sphere module¶

Provides concrete implementations of function spaces on the two-sphere (S²).

This module uses the abstract framework from the symmetric space module to create fully-featured Lebesgue (L²) and Sobolev (Hˢ) Hilbert spaces for functions defined on the surface of a sphere.

It utilizes the pyshtools library for highly efficient and accurate spherical harmonic transforms. Following a compositional design, this module first defines a base Lebesgue space and then constructs the Sobolev space as a MassWeightedHilbertSpace over it. The module also includes powerful plotting utilities built on cartopy for professional-quality geospatial visualization.

class pygeoinf.symmetric_space.sphere.Lebesgue(lmax: int, /, *, radius: float = 1, grid: str = 'DH', extend: bool = True)[source]¶

Bases: AbstractSymmetricLebesgueSpace

Implementation of the Lebesgue space L² on a sphere.

This class represents square-integrable functions on a sphere. A function is represented by its values on an evenly spaced grid. The co-ordinate basis for the space is through spherical harmonic expansions.

ax(a: float, x: SHGrid) → None[source]¶: Custom in-place ax implementation for pyshtools objects. x := a*x

axpy(a: float, x: SHGrid, y: SHGrid) → None[source]¶: Custom in-place axpy implementation for pyshtools objects. y := a*x + y

property csphase: int¶: The Condon-Shortley phase convention used (1).

degree_multiplicity(degree: int) → int[source]¶: Returns the number of eigenfunctions associated with a given degree.

degree_transfer_operator(target_degree: int) → LinearOperator[source]¶

Returns the transfer operator from this space to one with a different degree.

This operator leverages the hierarchical nature of the 1D SH vector to efficiently truncate or zero-pad the coefficients.

domain_mask(*, keep_ocean: bool = False, resolution: str = '110m') → SHGrid[source]¶

Generates a spatial mask that equals 1.0 in the specified domain and 0.0 outside of it.

Parameters:

keep_ocean – If True, the mask is 1.0 on the ocean and 0.0 on land. If False, the mask is 1.0 on land and 0.0 on the ocean.
resolution – The Cartopy shapefile resolution (‘110m’, ‘50m’, or ‘10m’).

Returns:

An SHGrid object representing the mask.

property extend: bool¶: True if the spatial grid includes both 0 and 360-degree longitudes.

from_coefficient_operator(lmax: int, /, *, lmin: int = 0)[source]¶

Returns a LinearOperator mapping a vector of coefficients to a function.

The operator maps a vector in (mathbb{R}^k) to an element of the Hilbert space. The input vector must follow the native pyshtools ordering: ordered by degree $l$ (major), and for each degree, the order $m$ is sorted as $0, 1, dots, l, -1, dots, -l$.

Parameters:

lmax – The maximum spherical harmonic degree expected in the input.
lmin – The minimum spherical harmonic degree expected. Defaults to 0.

Returns:

A LinearOperator mapping numpy.ndarray -> SHGrid.

from_coefficients(ulm: SHCoeffs) → SHGrid[source]¶: Maps spherical harmonic coefficients to a function vector.

from_components(c: ndarray) → SHGrid[source]¶

Maps a NumPy component array back to a vector in the space.

Parameters:: c – The components of the vector as a NumPy array.
Returns:: The corresponding vector in the space.

classmethod from_covariance(covariance_function: Callable[[float], float], /, *, radius: float = 1.0, grid: str = 'DH', extend: bool = True, rtol: float = 1e-06, min_degree: int = 0, max_degree: int | None = None, power_of_two: bool = False) → Lebesgue[source]¶

Creates an instance of the L² space with a spherical harmonic truncation degree (lmax) automatically chosen to capture the expected energy of functions drawn from a specified prior measure.

This factory method calculates the expected squared norm (energy) of a random field whose spectral variances are defined by the provided covariance_function. It iteratively adds higher spherical harmonic degrees until the relative contribution of the next degree drops below the specified relative tolerance.

Parameters:

covariance_function – A callable mapping a Laplacian eigenvalue to its spectral variance.
radius – The radius of the sphere. Defaults to 1.0.
grid – The pyshtools spatial grid format to use. Defaults to “DH”.
extend – If True, the spatial grid includes both 0 and 360-degree longitudes.
rtol – The relative tolerance for the energy truncation. Defaults to 1e-6.
min_degree – The absolute minimum truncation degree to return. Defaults to 0.
max_degree – An optional safety ceiling for the truncation degree. If convergence is not reached by this degree, the search terminates and returns this value.
power_of_two – If True, rounds the resulting lmax up to the nearest power of two.

Returns:

A fully instantiated L² space on the sphere with the optimal lmax.

Return type:

Lebesgue

Raises:

RuntimeError – If the energy sequence fails to converge within 100,000 iterations and no max_degree is specified.

classmethod from_heat_kernel_prior(kernel_scale: float, /, *, radius: float = 1.0, grid: str = 'DH', extend: bool = True, rtol: float = 1e-06, min_degree: int = 0, max_degree: int | None = None, power_of_two: bool = False) → Lebesgue[source]¶

Creates an instance of the L² space on the sphere, tuned to the expected energy of a Heat Kernel prior measure.

Parameters:

kernel_scale – The length-scale parameter of the heat kernel covariance.
radius – The radius of the sphere. Defaults to 1.0.
grid – The pyshtools spatial grid format to use. Defaults to “DH”.
extend – If True, the spatial grid includes both 0 and 360-degree longitudes.
rtol – The relative tolerance for the energy truncation. Defaults to 1e-6.
min_degree – The absolute minimum truncation degree to return. Defaults to 0.
max_degree – An optional safety ceiling for the truncation degree.
power_of_two – If True, rounds the resulting lmax up to the nearest power of two.

classmethod from_sobolev_kernel_prior(kernel_order: float, kernel_scale: float, /, *, radius: float = 1.0, grid: str = 'DH', extend: bool = True, rtol: float = 1e-06, min_degree: int = 0, max_degree: int | None = None, power_of_two: bool = False) → Lebesgue[source]¶

Creates an instance of the L² space on the sphere, tuned to the expected energy of a Sobolev-type prior measure.

Parameters:

kernel_order – The smoothness order of the Sobolev prior measure.
kernel_scale – The length-scale parameter of the Sobolev prior measure.
radius – The radius of the sphere. Defaults to 1.0.
grid – The pyshtools spatial grid format to use. Defaults to “DH”.
extend – If True, the spatial grid includes both 0 and 360-degree longitudes.
rtol – The relative tolerance for the energy truncation. Defaults to 1e-6.
min_degree – The absolute minimum truncation degree to return. Defaults to 0.
max_degree – An optional safety ceiling for the truncation degree.
power_of_two – If True, rounds the resulting lmax up to the nearest power of two.

property gaussian_curvature: float¶: Returns the Gaussian curvature (K) of the symmetric manifold.

geodesic_ball_integral(center: Tuple[float, float], radius: float, /, *, n_points: int | None = None, normalize: bool = False) → LinearForm[source]¶

Return the exact integral over a geodesic ball on the sphere.

The geodesic radius is given in the same physical units as geodesic_distance. On a sphere this corresponds to a spherical cap with angular radius radius / self.radius.

geodesic_ball_quadrature(center: Tuple[float, float], radius: float, n_points: int) → Tuple[List[Tuple[float, float]], ndarray][source]¶

Return a deterministic quadrature rule for a spherical cap.

The rule uses Gauss-Legendre nodes in the polar variable $u = cos(gamma)$ and a uniform trapezoidal rule in azimuth on each radial ring. The weights include the sphere area element, so their sum equals the physical cap area.

geodesic_distance(p1: Tuple[float, float], p2: Tuple[float, float]) → float[source]¶: Returns the great-circle distance between two points on the sphere.

geodesic_quadrature(p1: Tuple[float, float], p2: Tuple[float, float], n_points: int) → Tuple[List[Tuple[float, float]], ndarray][source]¶: Generates Gauss-Legendre quadrature points and weights along a great-circle arc.

property grid: str¶: The pyshtools grid type used for spatial representations.

property grid_type: str¶: Returns the pyshtools grid type.

index_to_integer(k: Tuple[int, int]) → int[source]¶: Maps an eigenvalue index to an integer.

integer_to_index(i: int) → Tuple[int, int][source]¶: Maps an integer to the eigenvalue index

invariant_covariance_function(spectral_variances: ndarray) → Callable[[ndarray], ndarray][source]¶: Returns a vectorized function that evaluates the two-point covariance of an invariant measure as a function of geodesic distance.

static iris_stations(*, n_stations: int = None, include_names: bool = False)[source]¶

Convenience method to retrieve a globally distributed set of Global Seismograph Network (GSN) stations from IRIS.

Parameters:

n_stations – The number of random stations to return. If None, returns the entire available network (~177 stations).
include_names – If True, returns (Name, Latitude, Longitude). If False, returns (Latitude, Longitude) tuples.

Returns:

A list of tuples representing the station coordinates in degrees.

is_element(x: Any) → bool[source]¶: Checks if an object is a valid element of the space.

laplacian_eigenvalue(k: Tuple[int, int]) → float[source]¶

Returns the (l.m)-th eigenvalue of the Laplacian.

Parameters:: = (k) – The index of the eigenvalue to return.

laplacian_eigenvector_squared_norm(k: Tuple[int, int]) → float[source]¶

Returns the squared norm of the eigenvalue of the Laplacian for a given index.

The index k can be a single integer (e.g., for a circle) or a tuple of integers (e.g., for a sphere or torus), depending on the geometry of the space.

Parameters:: k – The index of the eigenvalue to return.

laplacian_eigenvectors_at_point(point: Tuple[float, float]) → ndarray[source]¶

Returns a list of the values of the eigenvectors at a given point

Parameters:: x – The evaluation point.

property lmax: int¶: The maximum spherical harmonic truncation degree.

property normalization: str¶: The spherical harmonic normalization convention used (‘ortho’).

pairs_within_distance(points: List[Tuple[float, float]], max_distance: float) → Tuple[ndarray, ndarray, ndarray][source]¶: Finds all pairs of points within a given geodesic distance. Returns arrays of row indices, column indices, and their geodesic distances.

project_function(f: ~typing.Callable[[(<class 'float'>, <class 'float'>)], float]) → ndarray[source]¶

Returns an element of the space by projecting a given function.

Parameters:: f – A function that takes a point (lat, lon) and returns a value.

property radius: float¶: The radius of the sphere.

random_domain_points(n: int, /, *, keep_ocean: bool = True, resolution: str = '110m') → List[Tuple[float, float]][source]¶

Generates a list of random points strictly bounded to either the land or the ocean.

This uses Cartopy’s Natural Earth geometries and Shapely’s spatial indexing to perform ultra-fast rejection sampling.

Parameters:

n – The exact number of points to generate.
keep_ocean – If True, returns only ocean points. If False, returns only land points.
resolution – The Cartopy shapefile resolution (‘110m’, ‘50m’, or ‘10m’).

Returns:

A list of (latitude, longitude) tuples in degrees.

static random_earthquakes(n_points: int, min_magnitude: float = 5.0)[source]¶

Returns a random sample of real, sensible earthquake locations.

This provides a highly clustered, non-uniform spatial distribution along tectonic boundaries, perfect for testing robust spatial operators.

Parameters:

n_points – The exact number of locations to return.
min_magnitude – The minimum magnitude threshold for the catalog.

Returns:

A list of (Latitude, Longitude) tuples in degrees.

random_point() → Tuple[float, float][source]¶: Returns a random point as [latitude, longitude].

representative_index(degree: int) → Tuple[int, int][source]¶: Returns a valid eigenfunction index for the given degree.

sample_power_measure(measure, n_samples, /, *, lmin=None, lmax=None, parallel: bool = False, n_jobs: int = -1)[source]¶: Takes in a Gaussian measure on the space, draws n_samples from and returns samples for the spherical harmonic power at degrees in the indicated range.

property sampling: int¶: The sampling factor used for spatial representations.

spherical_cap_average(center: Tuple[float, float], angular_radius: float, /) → LinearForm[source]¶: Return the exact average functional over a spherical cap.

spherical_cap_integral(center: Tuple[float, float], angular_radius: float, /, *, normalize: bool = False) → LinearForm[source]¶

Return the exact spherical-cap integral functional.

The returned form represents $int_C u(x),dS(x)$, where C is the spherical cap centered at center with angular radius angular_radius. If normalize is true, the returned form represents the cap average $|C|^{-1}int_C u(x),dS(x)$.

The implementation uses pyshtools.SHCoeffs.from_cap to compute the spherical-harmonic coefficients of the cap indicator exactly in the truncated basis. The pyshtools cap is normalized to global average one, so it is rescaled to either the physical cap indicator or normalized cap average.

Parameters:

center – (latitude, longitude) cap center in degrees.
angular_radius – Cap radius in radians on the unit sphere.
normalize – If true, return the cap average instead of the raw surface integral.

to_coefficient_operator(lmax: int, /, *, lmin: int = 0)[source]¶

Returns a LinearOperator mapping a function to its spherical harmonic coefficients.

The operator maps an element of the Hilbert space to a vector in (mathbb{R}^k). The coefficients in the output vector follow the native pyshtools ordering: ordered by degree $l$ (major), and for each degree, the order $m$ is sorted as $0, 1, dots, l, -1, dots, -l$.

Parameters:

lmax – The maximum spherical harmonic degree to include in the output.
lmin – The minimum spherical harmonic degree to include. Defaults to 0.

Returns:

A LinearOperator mapping SHGrid -> numpy.ndarray.

to_coefficients(u: SHGrid) → SHCoeffs[source]¶: Maps a function vector to its spherical harmonic coefficients.

to_components(u: SHGrid) → ndarray[source]¶

Maps a vector to its representation as a NumPy component array.

Parameters:: x – A vector in the space.
Returns:: The components of the vector as a NumPy array.

vector_multiply(x1: SHGrid, x2: SHGrid) → SHGrid[source]¶: Computes the pointwise product of two functions.

vector_sqrt(x: SHGrid) → SHGrid[source]¶: Returns the pointwise square root of a function.

with_degree(degree: int) → Lebesgue[source]¶: Returns a new instance of the space with a modified truncation degree.

class pygeoinf.symmetric_space.sphere.Sobolev(lmax: int, order: float, scale: float, /, radius: float = 1, grid: str = 'DH', extend: bool = True, safe: bool = True)[source]¶

Bases: SymmetricSobolevSpace

Implementation of the Sobolev space Hˢ on a circle.

ax(a: float, x: SHGrid) → None[source]¶: Performs in-place scaling x := a*x. Defaults to x *= a.

axpy(a: float, x: SHGrid, y: SHGrid) → None[source]¶: Performs in-place operation y := y + a*x. Defaults to y += a*x.

property csphase: int¶: The Condon-Shortley phase convention used (1).

domain_mask(*, keep_ocean: bool = False, resolution: str = '110m') → SHGrid[source]¶

Generates a spatial mask that equals 1.0 in the specified domain and 0.0 outside of it.

Parameters:

keep_ocean – If True, the mask is 1.0 on the ocean and 0.0 on land. If False, the mask is 1.0 on land and 0.0 on the ocean.
resolution – The Cartopy shapefile resolution (‘110m’, ‘50m’, or ‘10m’).

Returns:

An SHGrid object representing the mask.

property extend: bool¶: True if the spatial grid includes both 0 and 360-degree longitudes.

from_coefficient_operator(lmax: int, /, *, lmin: int = 0)[source]¶

Returns a LinearOperator mapping a vector of coefficients to a function.

The operator maps a vector in $mathbb{R}^k$ to an element of the Hilbert space. The input vector must follow the standard $l$-major, $m$-minor ordering.

Ordering:

\[v = [u_{0,0}, \quad u_{1,-1}, u_{1,0}, u_{1,1}, \quad u_{2,-2}, \dots, u_{2,2}, \quad \dots]\]

(assuming lmin=0).

Parameters:

lmax – The maximum spherical harmonic degree expected in the input.
lmin – The minimum spherical harmonic degree expected. Defaults to 0.

Returns:

A LinearOperator mapping numpy.ndarray -> SHGrid.

from_coefficients(ulm: SHCoeffs) → SHGrid[source]¶: Maps spherical harmonic coefficients to a function vector.

classmethod from_covariance(covariance_function: Callable[[float], float], order: float, scale: float, /, *, radius: float = 1.0, grid: str = 'DH', extend: bool = True, rtol: float = 1e-06, min_degree: int = 0, max_degree: int | None = None, power_of_two: bool = False, safe: bool = True) → Sobolev[source]¶

Creates an instance of the Sobolev space with a spherical harmonic truncation degree (lmax) automatically chosen to capture the expected energy of functions drawn from a specified prior measure.

This factory method calculates the expected squared norm (energy) of a random field whose spectral variances are defined by the provided covariance_function, accounting for the Sobolev mass-weighting factor. It iteratively adds higher spherical harmonic degrees until the relative contribution of the next degree drops below the tolerance.

Parameters:

covariance_function – A callable mapping a Laplacian eigenvalue to its spectral variance.
order – The Sobolev order, controlling the smoothness of functions.
scale – The Sobolev length-scale.
radius – The radius of the sphere. Defaults to 1.0.
grid – The pyshtools spatial grid format to use. Defaults to “DH”.
extend – If True, the spatial grid includes both 0 and 360-degree longitudes.
rtol – The relative tolerance for the energy truncation. Defaults to 1e-6.
min_degree – The absolute minimum truncation degree to return. Defaults to 0.
max_degree – An optional safety ceiling for the truncation degree. If convergence is not reached by this degree, the search terminates and returns this value.
power_of_two – If True, rounds the resulting lmax up to the nearest power of two.
safe – If True, enables mathematical correctness checks during operations.

Returns:

A fully instantiated Sobolev space on the sphere with the optimal lmax.

Return type:

Sobolev

Raises:

RuntimeError – If the energy sequence fails to converge within 100,000 iterations and no max_degree is specified.

classmethod from_heat_kernel_prior(kernel_scale: float, order: float, scale: float, /, *, radius: float = 1.0, grid: str = 'DH', extend: bool = True, rtol: float = 1e-06, min_degree: int = 0, max_degree: int | None = None, power_of_two: bool = False, safe: bool = True) → Sobolev[source]¶

Creates an instance of the Sobolev space on the sphere, tuned to the expected energy of a Heat Kernel prior measure.

Parameters:

kernel_scale – The length-scale parameter of the heat kernel covariance.
order – The Sobolev order defining the function space.
scale – The Sobolev length-scale defining the function space.
radius – The radius of the sphere. Defaults to 1.0.
grid – The pyshtools spatial grid format to use. Defaults to “DH”.
extend – If True, the spatial grid includes both 0 and 360-degree longitudes.
rtol – The relative tolerance for the energy truncation. Defaults to 1e-6.
min_degree – The absolute minimum truncation degree to return. Defaults to 0.
max_degree – An optional safety ceiling for the truncation degree.
power_of_two – If True, rounds the resulting lmax up to the nearest power of two.
safe – If True, enables mathematical correctness checks during operations.

classmethod from_sobolev_kernel_prior(kernel_order: float, kernel_scale: float, order: float, scale: float, /, *, radius: float = 1.0, grid: str = 'DH', extend: bool = True, rtol: float = 1e-06, min_degree: int = 0, max_degree: int | None = None, power_of_two: bool = False, safe: bool = True) → Sobolev[source]¶

Creates an instance of the Sobolev space on the sphere, tuned to the expected energy of a Sobolev-type prior measure.

Parameters:

kernel_order – The smoothness order of the Sobolev prior measure.
kernel_scale – The length-scale parameter of the Sobolev prior measure.
order – The Sobolev order defining the function space.
scale – The Sobolev length-scale defining the function space.
radius – The radius of the sphere. Defaults to 1.0.
grid – The pyshtools spatial grid format to use. Defaults to “DH”.
extend – If True, the spatial grid includes both 0 and 360-degree longitudes.
rtol – The relative tolerance for the energy truncation. Defaults to 1e-6.
min_degree – The absolute minimum truncation degree to return. Defaults to 0.
max_degree – An optional safety ceiling for the truncation degree.
power_of_two – If True, rounds the resulting lmax up to the nearest power of two.
safe – If True, enables mathematical correctness checks during operations.

static from_sobolev_parameters(order: float, scale: float, /, *, radius: float = 1.0, grid: str = 'DH', rtol: float = 1e-08, power_of_two: bool = False, safe: bool = True) → Sobolev[source]¶

Creates an instance with lmax chosen based on the Sobolev parameters.

This factory method estimates the spherical harmonic truncation degree (lmax) required to represent the space while meeting a specified relative tolerance for the truncation error. This is useful when the required lmax is not known a priori.

Parameters:

order – The order of the Sobolev space, controlling smoothness.
scale – The non-dimensional length-scale for the space.
radius – The radius of the sphere. Defaults to 1.0.
grid – The pyshtools grid type (e.g., ‘DH’). Defaults to ‘DH’.
rtol – The relative tolerance used to determine the lmax.
power_of_two – If True, lmax is set to the next power of two.
safe – If true, the class checks for mathematical correctness of operations where possible.

Returns:

An instance of the Sobolev class with a calculated lmax.

property grid: str¶: The pyshtools grid type used for spatial representations.

property grid_type: str¶: Returns the pyshtools grid type.

static iris_stations(n_stations: int = None, include_names: bool = False)[source]¶

Convenience method to retrieve a globally distributed set of Global Seismograph Network (GSN) stations from IRIS.

Parameters:

n_stations – The number of random stations to return. If None, returns the entire available network (~177 stations).
include_names – If True, returns (Name, Latitude, Longitude). If False, returns (Latitude, Longitude) tuples.

Returns:

A list of tuples representing the station coordinates in degrees.

property lmax: int¶: The maximum spherical harmonic truncation degree.

property normalization: str¶: The spherical harmonic normalization convention used (‘ortho’).

path_average_operator(paths: List[Tuple[Tuple[float, float], Tuple[float, float]]], /, *, n_points: int | None = None, matrix_free: bool = False, parallel: bool = False, n_jobs: int = -1, lazy_quadrature: bool = False) → LinearOperator[source]¶

Constructs a tomographic operator mapping a function field to its line integrals along a set of geodesic paths.

Parameters:

paths – A list of start and end point pairs defining the geodesics.
n_points – The number of quadrature points per path.
matrix_free – If True, evaluates the integrals on-the-fly to save memory.
parallel – If True, distributes the matrix-free evaluations across CPU cores.
n_jobs – Number of parallel jobs. -1 uses all available cores.
lazy_quadrature – If True, computes quadrature points on-the-fly per chunk to prevent memory blowouts when evaluating massive numbers of paths.

point_evaluation_operator(points: List[Tuple[float, float]], /, *, matrix_free: bool = False, parallel: bool = False, n_jobs: int = -1) → LinearOperator[source]¶

Returns a linear operator that evaluates a function at a list of points.

Parameters:

points – A list of (latitude, longitude) tuples.
matrix_free – If True, evaluates the operator on-the-fly to save memory.
parallel – If True, distributes the matrix-free evaluations across CPU cores.
n_jobs – Number of parallel jobs. -1 uses all available cores.

project_function(f: ~typing.Callable[[(<class 'float'>, <class 'float'>)], float]) → ndarray[source]¶

Returns an element of the space by projecting a given function.

Parameters:: f – A function that takes a point (lat, lon) and returns a value.

property radius: float¶: The radius of the sphere.

random_domain_points(n: int, /, *, keep_ocean: bool = True, resolution: str = '110m') → List[Tuple[float, float]][source]¶: Generates a list of random points strictly bounded to either the land or the ocean.

static random_earthquakes(n_points: int, min_magnitude: float = 5.0)[source]¶

Returns a random sample of real, sensible earthquake locations.

This provides a highly clustered, non-uniform spatial distribution along tectonic boundaries, perfect for testing robust spatial operators.

Parameters:

n_points – The exact number of locations to return.
min_magnitude – The minimum magnitude threshold for the catalog.

Returns:

A list of (Latitude, Longitude) tuples in degrees.

sample_power_measure(measure, n_samples, /, *, lmin=None, lmax=None, parallel: bool = False, n_jobs: int = -1)[source]¶: Takes in a Gaussian measure on the space, draws n_samples from and returns samples for the spherical harmonic power at degrees in the indicated range.

property sampling: int¶: The sampling factor used for spatial representations.

to_coefficient_operator(lmax: int, /, *, lmin: int = 0)[source]¶

Returns a LinearOperator mapping a function to its spherical harmonic coefficients.

The operator maps an element of the Hilbert space to a vector in $mathbb{R}^k$. The coefficients in the output vector are ordered by degree $l$ (major) and order $m$ (minor), from $-l$ to $+l$.

Ordering:

\[u = [u_{0,0}, \quad u_{1,-1}, u_{1,0}, u_{1,1}, \quad u_{2,-2}, \dots, u_{2,2}, \quad \dots]\]

(assuming lmin=0).

Parameters:

lmax – The maximum spherical harmonic degree to include in the output.
lmin – The minimum spherical harmonic degree to include. Defaults to 0.

Returns:

A LinearOperator mapping SHGrid -> numpy.ndarray.

to_coefficients(u: SHGrid) → SHCoeffs[source]¶: Maps a function vector to its spherical harmonic coefficients.

with_degree(degree: int) → Sobolev[source]¶: Returns a new instance of the space with a modified truncation degree.

with_order(order: float) → Sobolev[source]¶: Returns a new instance of the exact same space but with a modified Sobolev order.

pygeoinf.symmetric_space.sphere.create_map_figure(figsize: Tuple[float, float] | None = None, projection: Projection | None = None, **kwargs) → Tuple[plt.Figure, GeoAxes][source]¶: Creates a modern Matplotlib Figure and Cartopy GeoAxes optimized for maps.

pygeoinf.symmetric_space.sphere.plot(u: SHGrid, /, *, ax: GeoAxes | Tuple[GeoAxes, Any] | None = None, projection: Projection | None = None, contour: bool = False, cmap: str = 'RdBu', coasts: bool = False, rivers: bool = False, borders: bool = False, map_extent: List[float] | None = None, gridlines: bool = True, gridlines_kwargs: dict | None = None, symmetric: bool | float = False, contour_lines: bool = False, contour_lines_kwargs: dict | None = None, num_levels: int = 10, colorbar: bool = False, colorbar_kwargs: dict | None = None, **kwargs) → Tuple[GeoAxes, Any][source]¶: Creates a high-quality map plot of a spherical harmonic function using Cartopy.

pygeoinf.symmetric_space.sphere.plot_geodesic(p1: Tuple[float, float], p2: Tuple[float, float], /, *, ax: GeoAxes | Tuple[GeoAxes, Any] | None = None, **kwargs) → Tuple[GeoAxes, Any][source]¶: Plots a geodesic (great-circle) curve between two points on the sphere.

pygeoinf.symmetric_space.sphere.plot_geodesic_network(paths: List[Tuple[Tuple[float, float], Tuple[float, float]]], /, *, ax: GeoAxes | Tuple[GeoAxes, Any] | None = None, plot_sources: bool = True, plot_receivers: bool = True, source_kwargs: dict | None = None, receiver_kwargs: dict | None = None, **kwargs) → Tuple[GeoAxes, List[Any]][source]¶: Plots a network of intersecting geodesic paths onto a Cartopy map.

pygeoinf.symmetric_space.sphere.plot_points(points: List[Tuple[float, float]], /, *, data: List[float] | np.ndarray | None = None, ax: GeoAxes | Tuple[GeoAxes, Any] | None = None, projection: Projection | None = None, cmap: str = 'RdBu', color: str = 'red', s: float = 20, marker: str = 'o', edgecolors: str = 'none', zorder: int = 5, coasts: bool = False, rivers: bool = False, borders: bool = False, map_extent: List[float] | None = None, gridlines: bool = True, gridlines_kwargs: dict | None = None, symmetric: bool | float = False, colorbar: bool = False, colorbar_kwargs: dict | None = None, **kwargs) → Tuple[GeoAxes, Any][source]¶: Plots discrete observation points on a map.

pygeoinf.symmetric_space.symmetric_space module¶

Provides an abstract framework for function spaces on symmetric manifolds.

This module offers an abstract framework for defining Hilbert spaces of functions on symmetric spaces (like spheres or tori). The core design leverages the spectral properties of the Laplace-Beltrami operator (Δ).

By inheriting from these base classes and implementing a small number of abstract methods a concrete class can automatically gain a rich set of tools for defining invariant operators and probability measures.

class pygeoinf.symmetric_space.symmetric_space.AbstractSymmetricLebesgueSpace(spatial_dim: int, degree: int | tuple[int, ...], dim: int)[source]¶

Bases: SymmetricHilbertSpace, HilbertModuleMixin, ABC

A specialisation for scalar-valued L² function spaces on symmetric manifolds.

To be instantiated, such a class must provide the following additional methods:

vector_multiply vector_sqrt

Constructs the exact 4th-order variable-coefficient flexure operator for a floating elastic plate mapped to a Riemannian manifold.

Mathematical Formulation:¶

This operator evaluates the covariant tensor divergence of the non-Euclidean moment equations. Under the geometer’s Laplacian convention, the exact equation is:

Op(w) = Delta(D Delta w)

Delta(D_{eff}) Delta w

tr(Hess(D_{eff}) Hess(w)) + 2 K nabla(D_{eff}) cdot nabla w

K D_{eff} Delta w

rho_g w

Notation:¶

$D_{eff} = D(1 - nu)$ is the effective flexural rigidity.
$K$ is the Gaussian curvature of the manifold.

Physical & Geometric Breakdown:¶

Principal Bending: Delta(D Delta w)
The standard biharmonic flexure term.
Scalar Property Interaction: - Delta(D_{eff}) Delta w
The scalar component of the Föppl-von Kármán bilinear bracket, accounting for gradients in the effective rigidity.
Tensor Gradients & Commutator: tr(Hess(D_{eff}) Hess(w)) + 2 K nabla(D_{eff}) cdot nabla w
Accounts for the coupling of spatial rigidity gradients with surface twist. The second half of this term (the gradient-curvature coupling) arises analytically from the non-commutativity of covariant derivatives. In this implementation, both terms are elegantly generated together in a coordinate-free manner by subtracting the Bochner identity block.
Geometric Covariant Softening: - K D_{eff} Delta w
Because covariant derivatives do not commute, pure bending on a curved surface inherently penalizes the manifold’s stiffness.
Foundation Restoring Force: rho_g w
The local hydrostatic restoring force (e.g., fluid mantle buoyancy).

param flexural_rigidity:: The flexural rigidity of the plate ($D$). Can be a spatially varying field (Vector) or a constant scalar.
type flexural_rigidity:: Union[Vector, float, int]
param poisson_ratio:: The Poisson’s ratio of the material ($nu$). Can be a spatially varying field (Vector) or a constant scalar. (Used to compute the effective rigidity $D_{eff}$).
type poisson_ratio:: Union[Vector, float, int]
param buoyancy_factor:: The hydrostatic foundation restoring force coefficient ($rho_g$). Can be a spatially varying field (Vector) or a constant scalar.
type buoyancy_factor:: Union[Vector, float, int]
returns:: A formally self-adjoint linear operator representing the exact 4th-order variable-coefficient flexural equations mapped to the manifold.
rtype:: LinearOperator

gradient_dot_product(f: Vector, g: Vector) → Vector[source]¶: Computes the pointwise dot product of the surface gradients of two fields.

Constructs the inverse of the exact variable-coefficient flexure operator.

For manifolds with constant properties, this method computes the exact inverse in a single shot via a spectral invariant automorphism. Under the geometer’s Laplacian convention (Delta), the geometric polynomial spectrum for eigenvalues k evaluates exactly to:

Eigenvalue(k) = D k^2 - K D_{eff} k + rho_g

Where $D_{eff} = D(1 - nu)$ and $K$ is the Gaussian curvature.

For varying properties, it utilizes a preconditioned Conjugate Gradient solver, using the spatially averaged properties to build the exact geometric spectral inverse as an ideal preconditioner.

Parameters:

flexural_rigidity (Union[Vector, float, int]) – The flexural rigidity of the plate ($D$).
poisson_ratio (Union[Vector, float, int]) – The Poisson’s ratio of the material ($nu$). (Used to compute the effective rigidity $D_{eff}$).
buoyancy_factor (Union[Vector, float, int]) – The hydrostatic foundation restoring force coefficient ($rho_g$).
baseline_rigidity (Optional[float], optional) – A constant scalar rigidity used to construct the geometric spectral preconditioner. If None, the spatial average of flexural_rigidity is computed dynamically.
baseline_buoyancy (Optional[float], optional) – A constant scalar foundation restoring force used to construct the geometric spectral preconditioner. If None, the spatial average of buoyancy_factor is computed dynamically.
solver (Optional[IterativeLinearSolver], optional) – An instance of an iterative linear solver used to invert the varying-coefficient operator. If None, a standard Conjugate Gradient solver (CGSolver) is instantiated.

Returns:

A linear operator representing the exact geometric spectral inverse (for constant properties) or the preconditioned iterative inverse (for varying properties).

Return type:

LinearOperator

spatial_multiplication_operator(f: Vector) → LinearOperator[source]¶: Returns a LinearOperator corresponding to the mapping u -> f * u with f a given function that is suitably regular.

class pygeoinf.symmetric_space.symmetric_space.InvariantGaussianMeasure(domain: SymmetricHilbertSpace, spectral_variances: ndarray, /, *, expectation: Vector | None = None)[source]¶

Bases: GaussianMeasure

A class for GaussianMeasures on SymmetricHilbertSpaces whose covariances are invariant under the symmetry group. The covariances can be expressed as a function of the Laplace-Beltrami operator and are diagonal within the eigenfunction basis.

affine_mapping(*, operator: LinearOperator = None, translation: Vector = None, affine_operator: AffineOperator = None, inverse_solver: LinearSolver = None, inverse_preconditioner: LinearOperator = None) → GaussianMeasure[source]¶

Transforms the measure under an affine map y = A(x) + b.

If the operator A is an InvariantLinearAutomorphism (or None), the resulting measure remains invariant. If it is a generic operator and an inverse solve is requested, a spectrally-regularized preconditioner is automatically built.

static from_function(domain: SymmetricHilbertSpace, f: Callable[[float], float], /, *, expectation=None) → InvariantGaussianMeasure[source]¶

Returns an invariant Gaussian measure on a symmetric Hilbert space whose covariance is of the form f(Δ) with f a function that is well-defined on the spectrum of the Laplacian, Δ.

Parameters:

domain – The domain of the operator
f – The function expressed in terms of the eigenvalue index
expectation – The expected value for the measure. Defaults to None which means zero.

static from_index_function(domain: SymmetricHilbertSpace, g: Callable[[int | tuple[int, ...]], float], /, *, expectation=None) → InvariantGaussianMeasure[source]¶

Returns an invariant Gaussian measure on a symmetric Hilbert space whose covariance is of the form f(Δ) with f a function that is well-defined on the spectrum of the Laplacian, Δ.

Here the function, f, is expressed implicitly as a function, g, of the eigenvalue index.

Parameters:

domain – The domain of the operator
g – The function expressed in terms of the eigenvalue index
expectation – The expected value for the measure. Defaults to None which means zero.

hilbert_schmidt_norm(*, exact: bool = False, size_estimate: int = 10, method: Literal['variable', 'fixed'] = 'variable', max_samples: int | None = None, rtol: float = 0.01, block_size: int = 10, parallel: bool = False, n_jobs: int = -1) → float[source]¶

Computes the Hilbert-Schmidt norm of the covariance operator.

By definition, ||C||_HS = sqrt( Tr(C* C) ). Because the covariance is self-adjoint, this evaluates to sqrt( Tr(C^2) ).

Parameters:

exact – If True, computes the exact norm by sweeping the basis.
size_estimate – Initial number of stochastic probe vectors.
method – ‘variable’ to sample dynamically until ‘rtol’ is met; ‘fixed’ otherwise.
max_samples – Hard limit on stochastic samples.
rtol – Relative tolerance for the stochastic estimator.
block_size – Number of probes evaluated between convergence checks.
parallel – If True, evaluates the vectors concurrently.
n_jobs – Number of CPU cores to utilize.

Returns:

The Hilbert-Schmidt norm of the covariance operator.

Return type:

float

Notes

Optimisations within the invariant class mean that the various optional arguments are not used in this instance.

kl_divergence(other: GaussianMeasure, /, *, method: Literal['auto', 'dense', 'randomized'] = 'auto', **kwargs) → float[source]¶

Computes the exact or approximate Kullback-Leibler (KL) divergence D_KL(self || other).

If both measures are InvariantGaussianMeasures on the same domain, and method=”auto”, this calculates the divergence in ultra-fast O(N) time directly from the spectral variances.

nuclear_norm(*, exact: bool = False, size_estimate: int = 10, method: Literal['variable', 'fixed'] = 'variable', max_samples: int | None = None, rtol: float = 0.01, block_size: int = 10, parallel: bool = False, n_jobs: int = -1) → float[source]¶

Computes the Nuclear norm (Trace norm) of the covariance operator.

Because a covariance operator is positive semi-definite, its Nuclear norm is mathematically identical to its trace: ||C||_* = Tr(C).

Parameters:

exact – If True, computes the exact norm by extracting the full diagonal of the component matrix. This is fast for diagonal/dense matrices but O(N^2) for abstract operators. If False, uses stochastic Hutchinson trace estimation.
size_estimate – Initial number of stochastic probe vectors.
method – ‘variable’ to sample dynamically until ‘rtol’ is met; ‘fixed’ otherwise.
max_samples – Hard limit on stochastic samples.
rtol – Relative tolerance for the stochastic estimator.
block_size – Number of probes evaluated between convergence checks.
parallel – If True, evaluates the vectors concurrently.
n_jobs – Number of CPU cores to utilize.

Returns:

The Nuclear norm of the covariance operator.

Return type:

float

Notes

Optimisations within the invariant class mean that the various optional arguments are not used in this instance.

rescale_norm_variance(std: float) → InvariantGaussianMeasure[source]¶

Returns a new measure whose covariance is scaled such that

E[||x - E[x]||^2] = std^2.

The expectation of the measure is unchanged.

property spectral_variances: ndarray¶: Provides instant access to the exact eigenvalues of the covariance operator, useful for log-determinant or KL-divergence calculations.

zero_expectation() → InvariantGaussianMeasure[source]¶: Returns a new invariant measure with the same covariance, but with expectation set to zero.

class pygeoinf.symmetric_space.symmetric_space.InvariantLinearAutomorphism(domain: SymmetricHilbertSpace, eigenvalues: ndarray)[source]¶

Bases: DiagonalSparseMatrixLinearOperator

A class for LinearOperators on SymmetricHilbertSpaces that are invariant under the symmetry group. Such operators can be expressed as a function of the Laplace-Beltrami operator and are diagonal within the eigenfunction basis.

apply_function(func: Callable[[ndarray], ndarray]) → InvariantLinearAutomorphism[source]¶

Applies a mathematical function to the operator via functional calculus.

Because this operator is invariant and diagonal in the spectral basis, the function is simply applied element-wise to the eigenvalues.

Parameters:: func – A NumPy ufunc or a callable that can be vectorized.
Returns:: A new InvariantLinearAutomorphism representing f(A).

property eigenvalues: ndarray¶: Returns the operator’s eigenvalues (the diagonal multipliers).

property exp: InvariantLinearAutomorphism¶: Returns the matrix exponential of the operator.

static from_function(domain: SymmetricHilbertSpace, f: Callable[[float], float]) → InvariantLinearAutomorphism[source]¶

Returns an invariant linear automorphism on a symmetric Hilbert space of the form f(Δ) with f a function that is well-defined on the spectrum of the Laplacian, Δ.

Parameters:

domain – The domain of the operator
f – The function

static from_index_function(domain: SymmetricHilbertSpace, g: Callable[[int | tuple[int, ...]], float]) → InvariantLinearAutomorphism[source]¶

Returns an invariant linear automorphism on a symmetric Hilbert space of the form f(Δ) with f a function that is well-defined on the spectrum of the Laplacian, Δ.

Here the function, f, is expressed implicitly as a function, g, of the eigenvalue index.

Parameters:

domain – The domain of the operator
g – The function expressed in terms of the eigenvalue index

property inverse: InvariantLinearAutomorphism¶

Returns the inverse of the invariant automorphism.

Since the operator is diagonal in the spectral basis, the inverse is formed by taking the reciprocal of the eigenvalues.

property log: InvariantLinearAutomorphism¶: Returns the matrix logarithm of the operator.

property sqrt: InvariantLinearAutomorphism¶: Returns the principal square root of the operator.

property trace: float¶: Returns the trace of the automorphism.

class pygeoinf.symmetric_space.symmetric_space.SymmetricHilbertSpace(spatial_dim: int, degree: int | tuple[int, ...], dim: int)[source]¶

Bases: OrthogonalHilbertSpace, ABC

An abstract base class for Hilbert spaces of functions spaces on symmetric manifolds.

The implementation is based on the expansion of elements of the space in terms of the eigenfunctions of the Laplace-Beltrami operator.

To inherit from this base class, the user must provide methods that provide information on the eigenvalues and eigenfunctions of the Laplace-Beltrami operator, including mappings to and from the co-ordinate basis. The eigenfunction basis is necessarily orthogonal, but need not be normalised.

cluster_points(points: List[float | tuple[float, ...]], /, *, threshold: float | None = None, n_clusters: int | None = None, linkage_method: str = 'complete') → List[List[int]][source]¶

Clusters a list of points into interacting blocks based on their geodesic distance.

This is particularly useful for generating the ‘interacting_blocks’ required by localized preconditioners in Bayesian inversions.

Parameters:

points – A list of points on the symmetric manifold.
threshold – The maximum geodesic distance between points in a cluster.
n_clusters – The exact number of clusters to form (alternative to threshold).
linkage_method – The hierarchical clustering method to use (‘complete’, ‘average’, ‘single’, etc.). Defaults to ‘complete’ to ensure clusters remain geographically compact.

Returns:

A list of lists, where each sub-list contains the indices of the points belonging to a specific cluster.

property degree: int | tuple[int, ...]¶: The spectral truncation degree of the space.

abstract degree_multiplicity(degree: int) → int[source]¶: Returns the number of eigenfunctions associated with a given degree.

abstract degree_transfer_operator(target_degree: int | tuple[int, ...]) → LinearOperator[source]¶: Returns the transfer operator from this space to one with a different degree.

estimate_truncation_degree(covariance_function: Callable[[float], float], /, *, rtol: float = 1e-06, min_degree: int = 0, max_degree: int | None = None) → int[source]¶: Estimates the truncation degree required to capture the expected energy of a function drawn from a prior measure.

abstract property gaussian_curvature: float¶: Returns the Gaussian curvature (K) of the symmetric manifold.

geodesic_ball_average(center: Any, radius: float, /, *, n_points: int | None = None) → LinearForm[source]¶: Return a linear form for the average over a geodesic ball.

geodesic_ball_integral(center: Any, radius: float, /, *, n_points: int | None = None, normalize: bool = False) → LinearForm[source]¶

Return a linear form for a geodesic-ball integral.

The returned form represents $int_{B(center, radius)} u(x),dV(x)$. If normalize is true, it represents the average over the geodesic ball instead.

Subclasses with exact spectral formulas can override this method. Otherwise, a subclass-provided geodesic_ball_quadrature rule is used.

Parameters:

center – Center point of the geodesic ball.
radius – Geodesic radius in the manifold’s distance units.
n_points – Number of quadrature points requested from geodesic_ball_quadrature.
normalize – If true, divide by the quadrature measure of the ball.

geodesic_ball_quadrature(center: Any, radius: float, n_points: int) → Tuple[List[Any], ndarray][source]¶

Return quadrature points and weights for a geodesic ball.

Subclasses may implement geometry-specific quadrature rules. The weights should include the manifold volume element so that sum_i weights[i] * f(points[i]) approximates $int_{B(center, radius)} f(x),dV(x)$.

abstract geodesic_distance(p1: Any, p2: Any) → float[source]¶

Returns the shortest distance along the manifold between two points.

Parameters:

p1 – The starting point.
p2 – The end point.

abstract geodesic_quadrature(p1: Any, p2: Any, n_points: int) → Tuple[List[Any], ndarray][source]¶

Returns quadrature points and weights for a geodesic between p1 and p2.

Returns:: List of manifold coordinates. weights: Integration weights scaled by the line element.
Return type:: points

static heat_kernel(scale: float) → Callable[[float], float][source]¶: Returns a heat kernel covariance function: exp(-scale^2 * k).

heat_kernel_gaussian_measure(scale: float, /, *, expectation: Vector = None)[source]¶

Returns an invariant Gaussian measure with a heat kernel covariance equal to exp(-scale^2 * Δ).

Parameters:

scale – Scale parameter for the covariance.
expectation – The expected value for measure. Defaults to zero

identity_operator() → InvariantLinearAutomorphism[source]¶

Returns the identity operator on the space.

Overridden to return an InvariantLinearAutomorphism, ensuring that compositions with other invariant operators preserve the fast-path diagonal structure.

abstract index_to_integer(k: int | tuple[int, ...]) → int[source]¶: Maps an eigenvalue index to an integer.

property indices: Iterator[int | tuple[int, ...]]¶: Returns a list of all the eigenvalue indices

abstract integer_to_index(i: int) → int | tuple[int, ...][source]¶: Maps an integer to the eigenvalue index

invariant_automorphism(f: Callable[[float], float]) → InvariantLinearAutomorphism[source]¶

Returns an automorphism of the form f(Δ) with f a function that is well-defined on the spectrum of the Laplacian, Δ.

Parameters:: f – The function

abstract invariant_covariance_function(spectral_variances: ndarray) → Callable[[ndarray], ndarray][source]¶: Returns a vectorized function that evaluates the two-point covariance of an invariant measure as a function of geodesic distance.

invariant_gaussian_measure(f: Callable[[float], float], /, *, expectation=None) → InvariantGaussianMeasure[source]¶

Returns a Gaussian measure on the space whose covariance takes the form f(Δ) with f a function that is well-defined on the spectrum of the Laplacian, Δ.

In order for the covariance to be well-defined, the trace of the covariance operator must be finite. This condition (in the sense of convergence as the size of the approximating space being increased) is not checked.

Parameters:

f – The function
expectation – The expected value for measure. Defaults to zero

Returns:

The measure as an instance of InvariantGaussianMeasure

l2_products_operator(weighting_functions: List[Any]) → LinearOperator[source]¶

Creates an operator that computes a vector of L2 inner products.

The action on a function ‘u’ returns a vector ‘d’ where d_i = <u, w_i>_L2.

Parameters:: weighting_functions – A list of vectors representing the weighting kernels.
Returns:: A mapping from this space to EuclideanSpace(len(weighting_functions)).
Return type:: LinearOperator

property laplacian: InvariantLinearAutomorphism¶: Returns the Laplacian, Δ, as an operator on the space. Note that this is an unbounded operator and it should be used with care.

abstract laplacian_eigenvalue(k: int | tuple[int, ...]) → float[source]¶

Returns the eigenvalue of the Laplacian for a given index.

The index k can be a single integer (e.g., for a circle) or a tuple of integers (e.g., for a sphere or torus), depending on the geometry of the space.

Parameters:: k – The index of the eigenvalue to return.

abstract laplacian_eigenvector_squared_norm(k: int | tuple[int, ...]) → float[source]¶

Returns the squared norm of the eigenvalue of the Laplacian for a given index.

The index k can be a single integer (e.g., for a circle) or a tuple of integers (e.g., for a sphere or torus), depending on the geometry of the space.

Parameters:: k – The index of the eigenvalue to return.

abstract laplacian_eigenvectors_at_point(x: float | tuple[float, ...]) → ndarray[source]¶

Returns a list of the values of the eigenvectors at a given point

Parameters:: x – The evaluation point.

norm_scaled_heat_kernel_gaussian_measure(scale: float, /, *, expectation: Vector = None, std: float = 1)[source]¶

Returns an invariant Gaussian measure with a heat kernel covariance equal to exp(-scale^2 * Δ).

Parameters:

scale – Scale parameter for the covariance.
expectation – The expected value for measure. Defaults to zero
std – The desired standard deviation

norm_scaled_invariant_gaussian_measure(f: Callable[[float], float], /, *, expectation=None, std=1.0) → InvariantGaussianMeasure[source]¶

Returns a Gaussian measure on the space whose covariance takes the form f(Δ) with f a function that is well-defined on the spectrum of the Laplacian, Δ.

The measure’s covariance is scaled such that:

E[||x||^2] = std*std + ||E[x]||^2

Parameters:

f – The function
expectation – The expected value for measure. Defaults to zero
std – The desired standard deviation

Returns:

The measure as an instance of InvariantGaussianMeasure

norm_scaled_sobolev_kernel_gaussian_measure(order: float, scale: float, /, *, expectation: Vector = None, std: float = 1)[source]¶

Returns an invariant Gaussian measure with a Sobolev-type covariance proportional to (1 + scale^2 * Δ)^-order.

The measure’s covariance is scaled such that the expected value for the samples norm is equal to the given standard deviation.

Parameters:

order – Order parameter for the covariance.
scale – Scale parameter for the covariance.
expectation – The expected value for measure. Defaults to zero
std – The desired standard deviation

pairs_within_distance(points: List[float | tuple[float, ...]], max_distance: float) → Tuple[ndarray, ndarray, ndarray][source]¶: Finds all pairs of points within a given geodesic distance. Returns arrays of row indices, column indices, and their geodesic distances.

abstract project_function(f: Callable[[float | tuple[float, ...]], float | ndarray]) → Vector[source]¶

Returns an element of the space by projecting a given function.

Parameters:: f – A function that takes a point and returns a value.

abstract random_point() → Any[source]¶: Returns a single random point from the underlying symmetric space.

random_points(n: int) → List[float | tuple[float, ...]][source]¶

Returns a list of n random points.

Parameters:: n – The number of random points to generate.

abstract representative_index(degree: int) → int | tuple[int, ...][source]¶: Returns a valid eigenfunction index for the given degree.

static sobolev_kernel(order: float, scale: float) → Callable[[float], float][source]¶: Returns a Sobolev kernel covariance function: (1 + scale^2 * k)^-order.

sobolev_kernel_gaussian_measure(order: float, scale: float, /, *, expectation: Vector = None)[source]¶

Returns an invariant Gaussian measure with a Sobolev-type covariance equal to (1 + scale^2 * Δ)^-order.

Parameters:

order – Order parameter for the covariance.
scale – Scale parameter for the covariance.
expectation – The expected value for measure. Defaults to zero

property spatial_dimension: int¶: The dimension of the symetric manifold.

abstract with_degree(degree: int | tuple[int, ...]) → SymmetricHilbertSpace[source]¶: Returns a new instance of the space with a modified truncation degree.

zero_operator(codomain: HilbertSpace | None = None) → InvariantLinearAutomorphism[source]¶

Returns the zero operator.

If the codomain is the space itself (or not provided), it returns an InvariantLinearAutomorphism for optimal arithmetic operations. Otherwise, it falls back to the standard abstract zero operator.

class pygeoinf.symmetric_space.symmetric_space.SymmetricSobolevSpace(lebesgue_space: AbstractSymmetricLebesgueSpace, order: float, scale: float, /, *, safe: bool | None = True)[source]¶

Bases: MassWeightedHilbertModule, AbstractSymmetricLebesgueSpace

The Sobolev space Hˢ constructed over a symmetric Lebesgue space.

This implementation leverages the mass-weighting framework to ensure that the inner product and dual mappings correctly account for the smoothness order and scale.

Clusters a list of points into interacting blocks based on their geodesic distance.

This is particularly useful for generating the ‘interacting_blocks’ required by localized preconditioners in Bayesian inversions.

Parameters:

points – A list of points on the symmetric manifold.
threshold – The maximum geodesic distance between points in a cluster.
n_clusters – The exact number of clusters to form (alternative to threshold).
linkage_method – The hierarchical clustering method to use (‘complete’, ‘average’, ‘single’, etc.). Defaults to ‘complete’ to ensure clusters remain geographically compact.

Returns:

A list of lists, where each sub-list contains the indices of the points belonging to a specific cluster.

degree_multiplicity(degree: int) → int[source]¶: Returns the number of eigenfunctions associated with a given degree.

degree_transfer_operator(target_degree: int | tuple[int, ...]) → LinearOperator[source]¶: Returns the transfer operator from this space to one with a different degree.

dirac(point: float | tuple[float, ...]) → LinearForm[source]¶: Returns the Dirac measure at the chosen point.

dirac_representation(point) → Vector[source]¶: Returns the representation of the Dirac measure at the chosen point.

distance_localized_preconditioner(prior_measure: InvariantGaussianMeasure, points: List[float | tuple[float, ...]], max_distance: float, /, *, data_error_measure: GaussianMeasure | None = None, apply_taper: bool = False, parallel: bool = False, n_jobs: int = -1) → LinearOperator[source]¶

Builds a highly specialized, ultra-fast sparse preconditioner for point evaluation problems when the prior measure is invariant.

This exploits the property that the two-point covariance depends solely on the geodesic distance. It maps out the 1D covariance function along a deterministic path, interpolates it, optionally applies a Gaspari-Cohn taper to ensure positive-definiteness, and populates a sparse normal matrix.

If max_distance is set to 0.0, it bypasses the distance calculations and instantly builds a purely diagonal (Jacobi) preconditioner.

Parameters:

prior_measure – The invariant prior Gaussian measure.
points – The list of observation points.
max_distance – The geodesic distance beyond which covariance is assumed to be zero (enforces sparsity). Set to 0.0 for a pure diagonal.
data_error_measure – The Gaussian measure describing the data noise. Can be None if the noise is already folded into the prior_measure.
apply_taper – If true, applies Gaspari-Cohn taper.
parallel – If True, computes the error measure diagonal in parallel.
n_jobs – Number of CPU cores to use if parallel=True.

Returns:

A LinearOperator representing the inverse of the approximated sparse normal matrix.

Constructs the exact 4th-order variable-coefficient flexure operator for a floating elastic plate mapped to a Riemannian manifold.

Mathematical Formulation:¶

This operator evaluates the covariant tensor divergence of the non-Euclidean moment equations. Under the geometer’s Laplacian convention, the exact equation is:

Op(w) = Delta(D Delta w)

Delta(D_{eff}) Delta w

tr(Hess(D_{eff}) Hess(w)) + 2 K nabla(D_{eff}) cdot nabla w

K D_{eff} Delta w

rho_g w

Notation:¶

$D_{eff} = D(1 - nu)$ is the effective flexural rigidity.
$K$ is the Gaussian curvature of the manifold.

Physical & Geometric Breakdown:¶

Principal Bending: Delta(D Delta w)
The standard biharmonic flexure term.
Scalar Property Interaction: - Delta(D_{eff}) Delta w
The scalar component of the Föppl-von Kármán bilinear bracket, accounting for gradients in the effective rigidity.
Tensor Gradients & Commutator: tr(Hess(D_{eff}) Hess(w)) + 2 K nabla(D_{eff}) cdot nabla w
Accounts for the coupling of spatial rigidity gradients with surface twist. The second half of this term (the gradient-curvature coupling) arises analytically from the non-commutativity of covariant derivatives. In this implementation, both terms are elegantly generated together in a coordinate-free manner by subtracting the Bochner identity block.
Geometric Covariant Softening: - K D_{eff} Delta w
Because covariant derivatives do not commute, pure bending on a curved surface inherently penalizes the manifold’s stiffness.
Foundation Restoring Force: rho_g w
The local hydrostatic restoring force (e.g., fluid mantle buoyancy).

param flexural_rigidity:: The flexural rigidity of the plate ($D$). Can be a spatially varying field (Vector) or a constant scalar.
type flexural_rigidity:: Union[Vector, float, int]
param poisson_ratio:: The Poisson’s ratio of the material ($nu$). Can be a spatially varying field (Vector) or a constant scalar. (Used to compute the effective rigidity $D_{eff}$).
type poisson_ratio:: Union[Vector, float, int]
param buoyancy_factor:: The hydrostatic foundation restoring force coefficient ($rho_g$). Can be a spatially varying field (Vector) or a constant scalar.
type buoyancy_factor:: Union[Vector, float, int]
returns:: A formally self-adjoint linear operator representing the exact 4th-order variable-coefficient flexural equations mapped to the manifold.
rtype:: LinearOperator

property gaussian_curvature: float¶: Returns the Gaussian curvature (K) of the symmetric manifold.

geodesic_ball_integral(center: Any, radius: float, /, *, n_points: int | None = None, normalize: bool = False) → LinearForm[source]¶

Return a linear form for a geodesic-ball integral.

The returned form represents $int_{B(center, radius)} u(x),dV(x)$. If normalize is true, it represents the average over the geodesic ball instead.

Subclasses with exact spectral formulas can override this method. Otherwise, a subclass-provided geodesic_ball_quadrature rule is used.

Parameters:

center – Center point of the geodesic ball.
radius – Geodesic radius in the manifold’s distance units.
n_points – Number of quadrature points requested from geodesic_ball_quadrature.
normalize – If true, divide by the quadrature measure of the ball.

geodesic_ball_quadrature(center: Any, radius: float, n_points: int) → Tuple[List[Any], ndarray][source]¶

Return quadrature points and weights for a geodesic ball.

geodesic_distance(p1: float | tuple[float, ...], p2: float | tuple[float, ...]) → float[source]¶

Returns the shortest distance along the manifold between two points.

Parameters:

p1 – The starting point.
p2 – The end point.

geodesic_integral(p1: float | tuple[float, ...], p2: float | tuple[float, ...], /, *, n_points: int | None = None) → LinearForm[source]¶

Returns a linear functional representing the line integral of a function along a geodesic path.

This method approximates the integral $\int_{\gamma} u(s) ds$, where $\gamma$ is the shortest path (geodesic) connecting points p1 and p2. The integral is represented as a LinearForm in the dual space, constructed by summing weighted point evaluations (Dirac measures) along the path.

For Hilbert spaces with a specified scale, the method can automatically determine the required quadrature density to resolve the smooth features of the space’s sensitivity kernels.

Parameters:

p1 – The starting point of the geodesic. The type is manifold-dependent (e.g., float for Circle, tuple for Sphere).
p2 – The end point of the geodesic.
n_points (int, optional) – The number of Gauss-Legendre quadrature points. If None, it is heuristically determined as: $n = \lceil (\text{arc\_length} / \text{scale}) \times 2 \rceil$. This ensures at least two points per characteristic length-scale, providing stable sampling of the sensitivity kernel. Defaults to None.

Returns:

A linear functional whose action on a vector u computes: the approximated line integral.

Return type:

LinearForm

Raises:

NotImplementedError – If the Sobolev order $s$ is less than or equal to half the spatial dimension $n/2$.

geodesic_integral_representation(p1: float | tuple[float, ...], p2: float | tuple[float, ...], /, *, n_points: int | None = None) → Vector[source]¶

Returns the Riesz representation (sensitivity kernel) of the line integral.

This maps the LinearForm (the integral functional) back into the primal Hilbert space. Visualizing this vector reveals the “sensitivity” of the line integral to perturbations at different locations in the domain.

Parameters:

p1 – Start and end points of the geodesic.
p2 – Start and end points of the geodesic.
n_points – Number of quadrature points.

geodesic_quadrature(p1: Any, p2: Any, n_points: int) → Tuple[List[Any], ndarray][source]¶

Returns quadrature points and weights for a geodesic between p1 and p2.

Returns:: List of manifold coordinates. weights: Integration weights scaled by the line element.
Return type:: points

index_to_integer(k: int | tuple[int, ...]) → int[source]¶: Maps an eigenvalue index to an integer.

integer_to_index(i: int) → int | tuple[int, ...][source]¶: Maps an integer to the eigenvalue index

invariant_covariance_function(spectral_variances: ndarray) → Callable[[ndarray], ndarray][source]¶: Returns a vectorized function that evaluates the two-point covariance of an invariant measure as a function of geodesic distance.

Constructs the inverse of the exact variable-coefficient flexure operator.

Eigenvalue(k) = D k^2 - K D_{eff} k + rho_g

Where $D_{eff} = D(1 - nu)$ and $K$ is the Gaussian curvature.

For varying properties, it utilizes a preconditioned Conjugate Gradient solver, using the spatially averaged properties to build the exact geometric spectral inverse as an ideal preconditioner.

Parameters:

flexural_rigidity (Union[Vector, float, int]) – The flexural rigidity of the plate ($D$).
poisson_ratio (Union[Vector, float, int]) – The Poisson’s ratio of the material ($nu$). (Used to compute the effective rigidity $D_{eff}$).
buoyancy_factor (Union[Vector, float, int]) – The hydrostatic foundation restoring force coefficient ($rho_g$).
baseline_rigidity (Optional[float], optional) – A constant scalar rigidity used to construct the geometric spectral preconditioner. If None, the spatial average of flexural_rigidity is computed dynamically.
baseline_buoyancy (Optional[float], optional) – A constant scalar foundation restoring force used to construct the geometric spectral preconditioner. If None, the spatial average of buoyancy_factor is computed dynamically.
solver (Optional[IterativeLinearSolver], optional) – An instance of an iterative linear solver used to invert the varying-coefficient operator. If None, a standard Conjugate Gradient solver (CGSolver) is instantiated.

Returns:

A linear operator representing the exact geometric spectral inverse (for constant properties) or the preconditioned iterative inverse (for varying properties).

Return type:

LinearOperator

l2_products_operator(weighting_functions: List[Any]) → LinearOperator[source]¶

Creates an operator that computes a vector of L2 inner products.

The action on a function ‘u’ returns a vector ‘d’ where d_i = <u, w_i>_L2.

Parameters:: weighting_functions – A list of vectors representing the weighting kernels.
Returns:: A mapping from this space to EuclideanSpace(len(weighting_functions)).
Return type:: LinearOperator

laplacian_eigenvalue(k: int | tuple[int, ...]) → float[source]¶

Returns the eigenvalue of the Laplacian for a given index.

The index k can be a single integer (e.g., for a circle) or a tuple of integers (e.g., for a sphere or torus), depending on the geometry of the space.

Parameters:: k – The index of the eigenvalue to return.

laplacian_eigenvector_squared_norm(k: int | tuple[int, ...]) → float[source]¶

Returns the squared norm of the eigenvalue of the Laplacian for a given index.

The index k can be a single integer (e.g., for a circle) or a tuple of integers (e.g., for a sphere or torus), depending on the geometry of the space.

Parameters:: k – The index of the eigenvalue to return.

laplacian_eigenvectors_at_point(x: float | tuple[float, ...]) → List[float | ndarray][source]¶

Returns a list of the values of the eigenvectors at a given point

Parameters:: x – The evaluation point.

property order: float¶: The Sobolev order.

order_inclusion_operator(target_order: float) → LinearOperator[source]¶: Returns the inclusion operator from this space to one of a lower order.

pairs_within_distance(points: List[float | tuple[float, ...]], max_distance: float) → Tuple[ndarray, ndarray, ndarray][source]¶: Finds all pairs of points within a given geodesic distance. Returns arrays of row indices, column indices, and their geodesic distances.

path_average_operator(paths, /, *, n_points=None)[source]¶

Constructs a tomographic operator mapping a function field to its line integrals along a set of geodesic paths.

Note: Despite the name, this operator returns the line integral (the dual pairing of the function with the path functional) rather than a normalized average, unless the user manually scales the forms. This corresponds to the ‘path average’ convention often used in seismic and atmospheric tomography.

Parameters:

paths (List[Tuple[Any, Any]]) – A list of start and end point pairs defining the geodesics.
n_points (int, optional) – The number of quadrature points per path. If None, the heuristic based on the Sobolev scale is used.

Returns:

An operator mapping Space -> EuclideanSpace(len(paths)).: The adjoint of this operator performs the ‘back-projection’ mapping data residuals into the function space.

Return type:

LinearOperator

point_evaluation_operator(points: List[Any]) → LinearOperator[source]¶

Returns a linear operator that evaluates a function at a list of points.

The resulting operator maps a function (a vector in this space) to a vector in Euclidean space containing the function’s values at the specified locations. This is the primary mechanism for creating a forward operator that links a function field to a set of discrete measurements.

Parameters:: points – A list of points at which to evaluate the functions.

point_value_scaled_heat_kernel_gaussian_measure(scale: float, /, *, expectation: Vector = None, std: float = 1)[source]¶

Returns an invariant Gaussian measure with a heat-kernel covariance proportional to exp(-scale^2 * Δ).

The covariance of the operator is scaled such that the standard deviation of the point-wise values are equal to the given amplitude.

Parameters:

scale – Scale parameter for the covariance.
std – The desired standard deviation for the pointwise values.

point_value_scaled_invariant_gaussian_measure(f: Callable[[float], float], /, *, expectation: Vector = None, std: float = 1)[source]¶

Returns an invariant Gaussian measure with covariance proportional to f(Δ), where f must be such that this operator is trace-class.

The covariance of the operator is scaled such that the standard deviation of the point-wise values are equal to the given std.

Parameters:

f – A real-valued function that is well-defined on the spectrum of the Laplacian, Δ.
std – The desired standard deviation for the pointwise values.

Raises:

NotImplementedError – If the Sobolev order is less than n/2, with n the spatial dimension.

Notes

This method applies for symmetric spaces an invariant measures. As a result, the pointwise variance is the same at all points. Internally, a random point is chosen to carry out the normalisation.

point_value_scaled_sobolev_kernel_gaussian_measure(order: float, scale: float, /, *, expectation: Vector = None, std: float = 1)[source]¶

Returns an invariant Gaussian measure with a Sobolev-type covariance proportional to (1 + scale^2 * Δ)^-order.

The covariance of the operator is scaled such that the standard deviation of the point-wise values are equal to the given amplitude.

Parameters:

order – Order parameter for the covariance.
scale – Scale parameter for the covariance.
std – The desired standard deviation for the pointwise values.

random_point() → float | tuple[float, ...][source]¶: Returns a single random point from the underlying symmetric space.

random_source_receiver_paths(n_sources: int, n_receivers: int) → List[Tuple[Any, Any]][source]¶

Generates a list of source-receiver pairs by connecting every source to every receiver.

This method uses the existing random_points() logic to generate coordinates appropriate for the specific symmetric space. For a set of S sources and R receivers, this returns a list of S*R paths.

Parameters:

n_sources – The number of random source locations to generate.
n_receivers – The number of random receiver locations to generate.

Returns:

A list of tuples, where each tuple contains: a (source, receiver) pair.

Return type:

List[Tuple[Any, Any]]

representative_index(degree: int) → int | tuple[int, ...][source]¶: Returns a valid eigenfunction index for the given degree.

property safe: bool¶: True if mathematical checks are being applied.

property scale: float¶: The Sobolev length-scale.

sobolev_function(lambda_val) → float[source]¶: Returns the value of the Sobolev function associated with the space.

spatial_multiplication_operator(f: Vector) → LinearOperator[source]¶: Returns a LinearOperator corresponding to the mapping u -> f * u with f a given function that is suitably regular.

abstract with_order(order: float) → SymmetricSobolevSpace[source]¶: Returns a new instance of the exact same space but with a modified Sobolev order.

pygeoinf.symmetric_space.torus module¶

Provides concrete implementations of function spaces on the 2D Torus.

This module defines Lebesgue (L²) and Sobolev (Hˢ) spaces for doubly periodic functions defined on a rectangular domain that wraps at its boundaries. It utilizes 2D Real Fast Fourier Transforms (rfft2) for highly efficient spectral operations and coordinate basis mapping.

class pygeoinf.symmetric_space.torus.Lebesgue(kmax: int, /, *, radius_x: float = 1.0, radius_y: float = 1.0)[source]¶

Bases: AbstractSymmetricLebesgueSpace

Implementation of the Lebesgue space L² on a 2D Torus.

Functions are represented by their values on an evenly spaced 2D grid. The coordinate basis for the space is formed by 2D Fourier expansions.

degree_multiplicity(degree: int) → int[source]¶: Returns the number of eigenfunctions associated with a given degree.

degree_transfer_operator(target_degree: int) → LinearOperator[source]¶: Returns a mapping operator that safely truncates or zero-pads Fourier coefficients to transition between different spectral resolutions.

estimate_truncation_degree(covariance_function: Callable[[float], float], /, *, rtol: float = 1e-06, min_degree: int = 1, max_degree: int | None = None) → int[source]¶: Iteratively calculates the spherical truncation degree needed to capture the variance of a random field below a specified error tolerance.

property fft_factor: float¶: The scaling factor applied during forward Fourier transforms.

from_coefficient_operator(kmax: int) → LinearOperator[source]¶: Returns an operator that maps spectral coefficients to a spatial function.

from_coefficients(coeff: ndarray) → ndarray[source]¶: Maps 2D complex Fourier coefficients back to a spatial function vector.

from_components(c: ndarray) → ndarray[source]¶: Unpacks a real-valued 1D component array back into spatial data.

classmethod from_covariance(covariance_function: Callable[[float], float], /, *, radius_x: float = 1.0, radius_y: float = 1.0, rtol: float = 1e-06, min_degree: int = 1, max_degree: int | None = None, power_of_two: bool = False) → Lebesgue[source]¶

Creates an instance of the L² space with a Fourier truncation degree (kmax) automatically chosen to capture the expected energy of functions drawn from a specified prior measure.

Parameters:

covariance_function – A callable mapping a Laplacian eigenvalue to its spectral variance.
radius_x – The effective radius of the Torus in the x-direction. Defaults to 1.0.
radius_y – The effective radius of the Torus in the y-direction. Defaults to 1.0.
rtol – The relative tolerance for the energy truncation. Defaults to 1e-6.
min_degree – The absolute minimum truncation degree to return. Defaults to 1.
max_degree – An optional safety ceiling for the truncation degree.
power_of_two – If True, rounds the resulting kmax up to the nearest power of two.

Returns:

A fully instantiated L² space on the Torus with the optimal kmax.

Return type:

Lebesgue

classmethod from_heat_kernel_prior(kernel_scale: float, /, *, radius_x: float = 1.0, radius_y: float = 1.0, rtol: float = 1e-06, min_degree: int = 1, max_degree: int | None = None, power_of_two: bool = False) → Lebesgue[source]¶

Creates an instance of the L² space on the Torus, tuned to the expected energy of a Heat Kernel prior measure.

Parameters:

kernel_scale – The length-scale parameter of the heat kernel covariance.
radius_x – The effective radius of the Torus in the x-direction. Defaults to 1.0.
radius_y – The effective radius of the Torus in the y-direction. Defaults to 1.0.
rtol – The relative tolerance for the energy truncation. Defaults to 1e-6.
min_degree – The absolute minimum truncation degree to return. Defaults to 1.
max_degree – An optional safety ceiling for the truncation degree.
power_of_two – If True, rounds the resulting kmax up to the nearest power of two.

classmethod from_sobolev_kernel_prior(kernel_order: float, kernel_scale: float, /, *, radius_x: float = 1.0, radius_y: float = 1.0, rtol: float = 1e-06, min_degree: int = 1, max_degree: int | None = None, power_of_two: bool = False) → Lebesgue[source]¶

Creates an instance of the L² space on the Torus, tuned to the expected energy of a Sobolev-type prior measure.

Parameters:

kernel_order – The smoothness order of the Sobolev prior measure.
kernel_scale – The length-scale parameter of the Sobolev prior measure.
radius_x – The effective radius of the Torus in the x-direction. Defaults to 1.0.
radius_y – The effective radius of the Torus in the y-direction. Defaults to 1.0.
rtol – The relative tolerance for the energy truncation. Defaults to 1e-6.
min_degree – The absolute minimum truncation degree to return. Defaults to 1.
max_degree – An optional safety ceiling for the truncation degree.
power_of_two – If True, rounds the resulting kmax up to the nearest power of two.

property gaussian_curvature: float¶: The Gaussian curvature of the Torus (always 0.0).

geodesic_distance(p1: Tuple[float, float], p2: Tuple[float, float]) → float[source]¶: Computes the shortest path distance between two points, accounting for periodic wrapping.

geodesic_quadrature(p1: Tuple[float, float], p2: Tuple[float, float], n_points: int) → Tuple[List[Tuple[float, float]], ndarray][source]¶: Generates quadrature points and weights along a shortest-path geodesic.

index_to_integer(k: int) → int[source]¶: Maps an eigenvalue index to an integer.

inner_product(x: ndarray, y: ndarray) → float[source]¶: Computes the exact L² inner product utilizing the spectral metric.

integer_to_index(i: int) → int[source]¶: Maps an integer to the eigenvalue index

invariant_covariance_function(spectral_variances: ndarray) → Callable[[ndarray], ndarray][source]¶: Raises NotImplementedError. The anisotropic layout of the 2D Torus prevents an exact 1D distance-to-covariance mapping.

property inverse_fft_factor: float¶: The scaling factor applied during inverse Fourier transforms.

is_element(x: Any) → bool[source]¶: Checks if the given object is a valid spatial representation in this space.

property kmax: int¶: The maximum Fourier degree represented in this space.

laplacian_eigenvalue(k: int) → float[source]¶: Returns the eigenvalue of the Laplacian for the given flat index.

laplacian_eigenvector_squared_norm(k: int) → float[source]¶: Returns the squared norm of the eigenfunction for the given flat index.

laplacian_eigenvectors_at_point(pt: Tuple[float, float]) → ndarray[source]¶: Evaluates all Laplacian eigenfunctions at a given physical coordinate.

norm(x: ndarray) → float[source]¶: Computes the exact L² norm utilizing the spectral metric.

points() → Tuple[ndarray, ndarray][source]¶: Returns two flattened 1D arrays representing the X and Y coordinates of the 2D grid.

project_function(f: Callable[[Tuple[float, float]], float]) → ndarray[source]¶

Returns an element of the space by projecting a given continuous 2D function.

Parameters:: f – A callable that takes a tuple (x, y) and returns a scalar value.
Returns:: A 2D numpy array representing the function evaluated on the spatial grid.

property radius_x: float¶: The radius of the Torus in the x-direction.

property radius_y: float¶: The radius of the Torus in the y-direction.

random_point() → Tuple[float, float][source]¶: Returns a random coordinate point within the Torus domain [0, 2π] x [0, 2π].

representative_index(degree: int) → int[source]¶: Returns a valid eigenfunction index for the given degree.

spectral_projection_operator(modes: List[Tuple[int, int]]) → LinearOperator[source]¶: Returns an operator that projects a function onto a specific set of Fourier wavevector modes.

to_coefficient_operator(kmax: int) → LinearOperator[source]¶: Returns an operator that maps a function to its spectral coefficients.

to_coefficients(u: ndarray) → ndarray[source]¶: Maps a spatial function vector to its 2D complex Fourier coefficients.

to_components(x: ndarray) → ndarray[source]¶: Flattens and packs spatial data into a real-valued 1D component array.

vector_multiply(x1: ndarray, x2: ndarray) → ndarray[source]¶: Pointwise multiplication of two spatial functions.

vector_sqrt(x: ndarray) → ndarray[source]¶: Pointwise square root of a spatial function.

wavevector_indices(kx: int, ky: int) → List[int][source]¶

Retrieves the 1D component indices corresponding to a specific 2D wavevector.

Parameters:

kx – The spatial frequency in the x-direction.
ky – The spatial frequency in the y-direction.

Returns:

A list of integers representing the indices in the flattened spectral vector.

with_degree(degree: int) → Lebesgue[source]¶: Returns a new space identical to this one but with a modified truncation degree.

property zero: ndarray¶: Returns a spatial field populated entirely with zeros.

class pygeoinf.symmetric_space.torus.Sobolev(kmax: int, order: float, scale: float, /, *, radius_x: float = 1.0, radius_y: float = 1.0, safe: bool = True)[source]¶

Bases: SymmetricSobolevSpace

Implementation of the Sobolev space Hˢ on a 2D Torus.

Constructs a mass-weighted Hilbert space over the base Torus Lebesgue geometry, accounting for fractional smoothness via a spectral order parameter.

estimate_truncation_degree(covariance_function: Callable[[float], float], /, *, rtol: float = 1e-06, min_degree: int = 1, max_degree: int | None = None) → int[source]¶: Delegates the energy truncation search to the underlying Torus Lebesgue space, ensuring it loops geometrically over (kx, ky) shells.

property fft_factor: float¶: The scaling factor applied during forward Fourier transforms.

from_coefficient_operator(kmax: int) → LinearOperator[source]¶: Returns an operator that maps spectral coefficients to a spatial function.

from_coefficients(coeff: ndarray) → ndarray[source]¶: Maps 2D complex Fourier coefficients back to a spatial function vector.

from_components(c: ndarray) → ndarray[source]¶: Unpacks a 1D component array back into spatial grid data.

classmethod from_covariance(covariance_function: Callable[[float], float], order: float, scale: float, /, *, radius_x: float = 1.0, radius_y: float = 1.0, rtol: float = 1e-06, min_degree: int = 1, max_degree: int | None = None, power_of_two: bool = False, safe: bool = True) → Sobolev[source]¶

Creates an instance with a Fourier truncation degree (kmax) automatically chosen to capture the expected energy of functions drawn from a specified prior.

Parameters:

covariance_function – A callable mapping a Laplacian eigenvalue to its spectral variance.
order – The Sobolev order defining the function space.
scale – The Sobolev length-scale defining the function space.
radius_x – The radius of the Torus in the x-direction.
radius_y – The radius of the Torus in the y-direction.
rtol – The relative tolerance for the energy truncation. Defaults to 1e-6.
min_degree – The absolute minimum truncation degree to return.
max_degree – An optional safety ceiling for the truncation degree.
power_of_two – If True, rounds the resulting kmax up to the nearest power of two.
safe – If True, enables mathematical correctness checks.

classmethod from_heat_kernel_prior(kernel_scale: float, order: float, scale: float, /, *, radius_x: float = 1.0, radius_y: float = 1.0, rtol: float = 1e-06, min_degree: int = 1, max_degree: int | None = None, power_of_two: bool = False, safe: bool = True) → Sobolev[source]¶: Creates an instance of the Sobolev space tuned to the expected energy of a Heat Kernel prior measure.

classmethod from_sobolev_kernel_prior(kernel_order: float, kernel_scale: float, order: float, scale: float, /, *, radius_x: float = 1.0, radius_y: float = 1.0, rtol: float = 1e-06, min_degree: int = 1, max_degree: int | None = None, power_of_two: bool = False, safe: bool = True) → Sobolev[source]¶: Creates an instance of the Sobolev space tuned to the expected energy of a Sobolev-type prior measure.

static from_sobolev_parameters(order: float, scale: float, /, *, radius_x: float = 1.0, radius_y: float = 1.0, rtol: float = 1e-06, power_of_two: bool = False, safe: bool = True) → Sobolev[source]¶

Creates an instance with kmax automatically chosen to limit spatial truncation error.

Parameters:

order – The Sobolev order. Must be > 1.0 for valid point evaluation on a 2D domain.
scale – The Sobolev length-scale parameter.
radius_x – The radius of the Torus in the x-direction.
radius_y – The radius of the Torus in the y-direction.
rtol – Relative tolerance used in assessing truncation error.
power_of_two – If True, rounds kmax up to the next power of two.
safe – If True, enables mathematical correctness checks.

inner_product(x: ndarray, y: ndarray) → float[source]¶: Computes the exact Sobolev inner product utilizing the mass-weighted metric.

property inverse_fft_factor: float¶: The scaling factor applied during inverse Fourier transforms.

is_element(x: Any) → bool[source]¶: Checks if the given object is a valid spatial representation in this space.

property kmax: int¶: The maximum Fourier degree represented in this space.

norm(x: ndarray) → float[source]¶: Computes the exact Sobolev norm utilizing the mass-weighted metric.

points() → Tuple[ndarray, ndarray][source]¶: Returns two flattened 1D arrays representing the X and Y grid coordinates.

project_function(f: Callable[[Tuple[float, float]], float]) → ndarray[source]¶: Evaluates and projects a continuous function onto the spatial grid.

property radius_x: float¶: The radius of the Torus in the x-direction.

property radius_y: float¶: The radius of the Torus in the y-direction.

spectral_projection_operator(modes: List[Tuple[int, int]]) → LinearOperator[source]¶: Returns an operator projecting onto specific Fourier wavevector modes.

to_coefficient_operator(kmax: int) → LinearOperator[source]¶: Returns an operator that maps a function to its spectral coefficients.

to_coefficients(u: ndarray) → ndarray[source]¶: Maps a spatial function vector to its 2D complex Fourier coefficients.

to_components(x: ndarray) → ndarray[source]¶: Flattens and packs spatial data into a dense 1D component array.

wavevector_indices(kx: int, ky: int) → List[int][source]¶: Retrieves the 1D component indices corresponding to a specific 2D wavevector.

with_degree(degree: int) → Sobolev[source]¶: Returns a new Sobolev instance with a modified truncation degree.

with_order(order: float) → Sobolev[source]¶: Returns a new Sobolev instance with a modified order.

pygeoinf.symmetric_space.torus.plot(space: Lebesgue | Sobolev, u: ndarray, /, *, ax: Axes | Tuple[Axes, Any] | None = None, contour: bool = False, cmap: str = 'RdBu', symmetric: bool | float = False, contour_lines: bool = False, contour_lines_kwargs: dict | None = None, num_levels: int = 10, colorbar: bool = False, colorbar_kwargs: dict | None = None, **kwargs) → Tuple[Axes, Any][source]¶: Creates a high-quality plot of a function on the 2D Torus.

pygeoinf.symmetric_space.torus.plot_geodesic(p1: Tuple[float, float], p2: Tuple[float, float], /, *, ax: Axes | Tuple[Axes, Any] | None = None, **kwargs) → Tuple[Axes, Any][source]¶: Plots a straight-line geodesic, wrapping around the Torus boundaries if necessary.

pygeoinf.symmetric_space.torus.plot_geodesic_network(paths: List[Tuple[Tuple[float, float], Tuple[float, float]]], /, *, ax: Axes | Tuple[Axes, Any] | None = None, plot_sources: bool = True, plot_receivers: bool = True, source_kwargs: dict | None = None, receiver_kwargs: dict | None = None, **kwargs) → Tuple[Axes, List[Any]][source]¶: Plots a network of intersecting geodesic paths on the Torus.

pygeoinf.symmetric_space.torus.plot_points(points: List[Tuple[float, float]], /, *, data: List[float] | ndarray | None = None, ax: Axes | Tuple[Axes, Any] | None = None, cmap: str = 'RdBu', color: str = 'red', s: float = 20, marker: str = 'o', edgecolors: str = 'none', zorder: int = 5, symmetric: bool | float = False, colorbar: bool = False, colorbar_kwargs: dict | None = None, **kwargs) → Tuple[Axes, Any][source]¶: Plots discrete observation points wrapped onto the primary Torus domain [0, 2π] x [0, 2π].

Module contents¶

Symmetric space function spaces and utilities.

class pygeoinf.symmetric_space.AbstractSymmetricLebesgueSpace(spatial_dim: int, degree: int | tuple[int, ...], dim: int)[source]¶

Bases: SymmetricHilbertSpace, HilbertModuleMixin, ABC

A specialisation for scalar-valued L² function spaces on symmetric manifolds.

To be instantiated, such a class must provide the following additional methods:

vector_multiply vector_sqrt

Constructs the exact 4th-order variable-coefficient flexure operator for a floating elastic plate mapped to a Riemannian manifold.

Mathematical Formulation:¶

This operator evaluates the covariant tensor divergence of the non-Euclidean moment equations. Under the geometer’s Laplacian convention, the exact equation is:

Op(w) = Delta(D Delta w)

Delta(D_{eff}) Delta w

tr(Hess(D_{eff}) Hess(w)) + 2 K nabla(D_{eff}) cdot nabla w

K D_{eff} Delta w

rho_g w

Notation:¶

$D_{eff} = D(1 - nu)$ is the effective flexural rigidity.
$K$ is the Gaussian curvature of the manifold.

Physical & Geometric Breakdown:¶

Principal Bending: Delta(D Delta w)
The standard biharmonic flexure term.
Scalar Property Interaction: - Delta(D_{eff}) Delta w
The scalar component of the Föppl-von Kármán bilinear bracket, accounting for gradients in the effective rigidity.
Tensor Gradients & Commutator: tr(Hess(D_{eff}) Hess(w)) + 2 K nabla(D_{eff}) cdot nabla w
Accounts for the coupling of spatial rigidity gradients with surface twist. The second half of this term (the gradient-curvature coupling) arises analytically from the non-commutativity of covariant derivatives. In this implementation, both terms are elegantly generated together in a coordinate-free manner by subtracting the Bochner identity block.
Geometric Covariant Softening: - K D_{eff} Delta w
Because covariant derivatives do not commute, pure bending on a curved surface inherently penalizes the manifold’s stiffness.
Foundation Restoring Force: rho_g w
The local hydrostatic restoring force (e.g., fluid mantle buoyancy).

param flexural_rigidity:: The flexural rigidity of the plate ($D$). Can be a spatially varying field (Vector) or a constant scalar.
type flexural_rigidity:: Union[Vector, float, int]
param poisson_ratio:: The Poisson’s ratio of the material ($nu$). Can be a spatially varying field (Vector) or a constant scalar. (Used to compute the effective rigidity $D_{eff}$).
type poisson_ratio:: Union[Vector, float, int]
param buoyancy_factor:: The hydrostatic foundation restoring force coefficient ($rho_g$). Can be a spatially varying field (Vector) or a constant scalar.
type buoyancy_factor:: Union[Vector, float, int]
returns:: A formally self-adjoint linear operator representing the exact 4th-order variable-coefficient flexural equations mapped to the manifold.
rtype:: LinearOperator

gradient_dot_product(f: Vector, g: Vector) → Vector[source]¶: Computes the pointwise dot product of the surface gradients of two fields.

Constructs the inverse of the exact variable-coefficient flexure operator.

Eigenvalue(k) = D k^2 - K D_{eff} k + rho_g

Where $D_{eff} = D(1 - nu)$ and $K$ is the Gaussian curvature.

For varying properties, it utilizes a preconditioned Conjugate Gradient solver, using the spatially averaged properties to build the exact geometric spectral inverse as an ideal preconditioner.

Parameters:

flexural_rigidity (Union[Vector, float, int]) – The flexural rigidity of the plate ($D$).
poisson_ratio (Union[Vector, float, int]) – The Poisson’s ratio of the material ($nu$). (Used to compute the effective rigidity $D_{eff}$).
buoyancy_factor (Union[Vector, float, int]) – The hydrostatic foundation restoring force coefficient ($rho_g$).
baseline_rigidity (Optional[float], optional) – A constant scalar rigidity used to construct the geometric spectral preconditioner. If None, the spatial average of flexural_rigidity is computed dynamically.
baseline_buoyancy (Optional[float], optional) – A constant scalar foundation restoring force used to construct the geometric spectral preconditioner. If None, the spatial average of buoyancy_factor is computed dynamically.
solver (Optional[IterativeLinearSolver], optional) – An instance of an iterative linear solver used to invert the varying-coefficient operator. If None, a standard Conjugate Gradient solver (CGSolver) is instantiated.

Returns:

A linear operator representing the exact geometric spectral inverse (for constant properties) or the preconditioned iterative inverse (for varying properties).

Return type:

LinearOperator

spatial_multiplication_operator(f: Vector) → LinearOperator[source]¶: Returns a LinearOperator corresponding to the mapping u -> f * u with f a given function that is suitably regular.

class pygeoinf.symmetric_space.InvariantGaussianMeasure(domain: SymmetricHilbertSpace, spectral_variances: ndarray, /, *, expectation: Vector | None = None)[source]¶

Bases: GaussianMeasure

Transforms the measure under an affine map y = A(x) + b.

static from_function(domain: SymmetricHilbertSpace, f: Callable[[float], float], /, *, expectation=None) → InvariantGaussianMeasure[source]¶

Returns an invariant Gaussian measure on a symmetric Hilbert space whose covariance is of the form f(Δ) with f a function that is well-defined on the spectrum of the Laplacian, Δ.

Parameters:

domain – The domain of the operator
f – The function expressed in terms of the eigenvalue index
expectation – The expected value for the measure. Defaults to None which means zero.

static from_index_function(domain: SymmetricHilbertSpace, g: Callable[[int | tuple[int, ...]], float], /, *, expectation=None) → InvariantGaussianMeasure[source]¶

Returns an invariant Gaussian measure on a symmetric Hilbert space whose covariance is of the form f(Δ) with f a function that is well-defined on the spectrum of the Laplacian, Δ.

Here the function, f, is expressed implicitly as a function, g, of the eigenvalue index.

Parameters:

domain – The domain of the operator
g – The function expressed in terms of the eigenvalue index
expectation – The expected value for the measure. Defaults to None which means zero.

Computes the Hilbert-Schmidt norm of the covariance operator.

By definition, ||C||_HS = sqrt( Tr(C* C) ). Because the covariance is self-adjoint, this evaluates to sqrt( Tr(C^2) ).

Parameters:

exact – If True, computes the exact norm by sweeping the basis.
size_estimate – Initial number of stochastic probe vectors.
method – ‘variable’ to sample dynamically until ‘rtol’ is met; ‘fixed’ otherwise.
max_samples – Hard limit on stochastic samples.
rtol – Relative tolerance for the stochastic estimator.
block_size – Number of probes evaluated between convergence checks.
parallel – If True, evaluates the vectors concurrently.
n_jobs – Number of CPU cores to utilize.

Returns:

The Hilbert-Schmidt norm of the covariance operator.

Return type:

float

Notes

Optimisations within the invariant class mean that the various optional arguments are not used in this instance.

kl_divergence(other: GaussianMeasure, /, *, method: Literal['auto', 'dense', 'randomized'] = 'auto', **kwargs) → float[source]¶

Computes the exact or approximate Kullback-Leibler (KL) divergence D_KL(self || other).

If both measures are InvariantGaussianMeasures on the same domain, and method=”auto”, this calculates the divergence in ultra-fast O(N) time directly from the spectral variances.

Computes the Nuclear norm (Trace norm) of the covariance operator.

Because a covariance operator is positive semi-definite, its Nuclear norm is mathematically identical to its trace: ||C||_* = Tr(C).

Parameters:

exact – If True, computes the exact norm by extracting the full diagonal of the component matrix. This is fast for diagonal/dense matrices but O(N^2) for abstract operators. If False, uses stochastic Hutchinson trace estimation.
size_estimate – Initial number of stochastic probe vectors.
method – ‘variable’ to sample dynamically until ‘rtol’ is met; ‘fixed’ otherwise.
max_samples – Hard limit on stochastic samples.
rtol – Relative tolerance for the stochastic estimator.
block_size – Number of probes evaluated between convergence checks.
parallel – If True, evaluates the vectors concurrently.
n_jobs – Number of CPU cores to utilize.

Returns:

The Nuclear norm of the covariance operator.

Return type:

float

Notes

Optimisations within the invariant class mean that the various optional arguments are not used in this instance.

rescale_norm_variance(std: float) → InvariantGaussianMeasure[source]¶

Returns a new measure whose covariance is scaled such that

E[||x - E[x]||^2] = std^2.

The expectation of the measure is unchanged.

property spectral_variances: ndarray¶: Provides instant access to the exact eigenvalues of the covariance operator, useful for log-determinant or KL-divergence calculations.

zero_expectation() → InvariantGaussianMeasure[source]¶: Returns a new invariant measure with the same covariance, but with expectation set to zero.

class pygeoinf.symmetric_space.InvariantLinearAutomorphism(domain: SymmetricHilbertSpace, eigenvalues: ndarray)[source]¶

Bases: DiagonalSparseMatrixLinearOperator

apply_function(func: Callable[[ndarray], ndarray]) → InvariantLinearAutomorphism[source]¶

Applies a mathematical function to the operator via functional calculus.

Because this operator is invariant and diagonal in the spectral basis, the function is simply applied element-wise to the eigenvalues.

Parameters:: func – A NumPy ufunc or a callable that can be vectorized.
Returns:: A new InvariantLinearAutomorphism representing f(A).

property eigenvalues: ndarray¶: Returns the operator’s eigenvalues (the diagonal multipliers).

property exp: InvariantLinearAutomorphism¶: Returns the matrix exponential of the operator.

static from_function(domain: SymmetricHilbertSpace, f: Callable[[float], float]) → InvariantLinearAutomorphism[source]¶

Returns an invariant linear automorphism on a symmetric Hilbert space of the form f(Δ) with f a function that is well-defined on the spectrum of the Laplacian, Δ.

Parameters:

domain – The domain of the operator
f – The function

static from_index_function(domain: SymmetricHilbertSpace, g: Callable[[int | tuple[int, ...]], float]) → InvariantLinearAutomorphism[source]¶

Returns an invariant linear automorphism on a symmetric Hilbert space of the form f(Δ) with f a function that is well-defined on the spectrum of the Laplacian, Δ.

Here the function, f, is expressed implicitly as a function, g, of the eigenvalue index.

Parameters:

domain – The domain of the operator
g – The function expressed in terms of the eigenvalue index

property inverse: InvariantLinearAutomorphism¶

Returns the inverse of the invariant automorphism.

Since the operator is diagonal in the spectral basis, the inverse is formed by taking the reciprocal of the eigenvalues.

property log: InvariantLinearAutomorphism¶: Returns the matrix logarithm of the operator.

property sqrt: InvariantLinearAutomorphism¶: Returns the principal square root of the operator.

property trace: float¶: Returns the trace of the automorphism.

class pygeoinf.symmetric_space.SymmetricHilbertSpace(spatial_dim: int, degree: int | tuple[int, ...], dim: int)[source]¶

Bases: OrthogonalHilbertSpace, ABC

An abstract base class for Hilbert spaces of functions spaces on symmetric manifolds.

The implementation is based on the expansion of elements of the space in terms of the eigenfunctions of the Laplace-Beltrami operator.

Clusters a list of points into interacting blocks based on their geodesic distance.

This is particularly useful for generating the ‘interacting_blocks’ required by localized preconditioners in Bayesian inversions.

Parameters:

points – A list of points on the symmetric manifold.
threshold – The maximum geodesic distance between points in a cluster.
n_clusters – The exact number of clusters to form (alternative to threshold).
linkage_method – The hierarchical clustering method to use (‘complete’, ‘average’, ‘single’, etc.). Defaults to ‘complete’ to ensure clusters remain geographically compact.

Returns:

A list of lists, where each sub-list contains the indices of the points belonging to a specific cluster.

property degree: int | tuple[int, ...]¶: The spectral truncation degree of the space.

abstract degree_multiplicity(degree: int) → int[source]¶: Returns the number of eigenfunctions associated with a given degree.

abstract degree_transfer_operator(target_degree: int | tuple[int, ...]) → LinearOperator[source]¶: Returns the transfer operator from this space to one with a different degree.

estimate_truncation_degree(covariance_function: Callable[[float], float], /, *, rtol: float = 1e-06, min_degree: int = 0, max_degree: int | None = None) → int[source]¶: Estimates the truncation degree required to capture the expected energy of a function drawn from a prior measure.

abstract property gaussian_curvature: float¶: Returns the Gaussian curvature (K) of the symmetric manifold.

geodesic_ball_average(center: Any, radius: float, /, *, n_points: int | None = None) → LinearForm[source]¶: Return a linear form for the average over a geodesic ball.

geodesic_ball_integral(center: Any, radius: float, /, *, n_points: int | None = None, normalize: bool = False) → LinearForm[source]¶

Return a linear form for a geodesic-ball integral.

The returned form represents $int_{B(center, radius)} u(x),dV(x)$. If normalize is true, it represents the average over the geodesic ball instead.

Subclasses with exact spectral formulas can override this method. Otherwise, a subclass-provided geodesic_ball_quadrature rule is used.

Parameters:

center – Center point of the geodesic ball.
radius – Geodesic radius in the manifold’s distance units.
n_points – Number of quadrature points requested from geodesic_ball_quadrature.
normalize – If true, divide by the quadrature measure of the ball.

geodesic_ball_quadrature(center: Any, radius: float, n_points: int) → Tuple[List[Any], ndarray][source]¶

Return quadrature points and weights for a geodesic ball.

abstract geodesic_distance(p1: Any, p2: Any) → float[source]¶

Returns the shortest distance along the manifold between two points.

Parameters:

p1 – The starting point.
p2 – The end point.

abstract geodesic_quadrature(p1: Any, p2: Any, n_points: int) → Tuple[List[Any], ndarray][source]¶

Returns quadrature points and weights for a geodesic between p1 and p2.

Returns:: List of manifold coordinates. weights: Integration weights scaled by the line element.
Return type:: points

static heat_kernel(scale: float) → Callable[[float], float][source]¶: Returns a heat kernel covariance function: exp(-scale^2 * k).

heat_kernel_gaussian_measure(scale: float, /, *, expectation: Vector = None)[source]¶

Returns an invariant Gaussian measure with a heat kernel covariance equal to exp(-scale^2 * Δ).

Parameters:

scale – Scale parameter for the covariance.
expectation – The expected value for measure. Defaults to zero

identity_operator() → InvariantLinearAutomorphism[source]¶

Returns the identity operator on the space.

Overridden to return an InvariantLinearAutomorphism, ensuring that compositions with other invariant operators preserve the fast-path diagonal structure.

abstract index_to_integer(k: int | tuple[int, ...]) → int[source]¶: Maps an eigenvalue index to an integer.

property indices: Iterator[int | tuple[int, ...]]¶: Returns a list of all the eigenvalue indices

abstract integer_to_index(i: int) → int | tuple[int, ...][source]¶: Maps an integer to the eigenvalue index

invariant_automorphism(f: Callable[[float], float]) → InvariantLinearAutomorphism[source]¶

Returns an automorphism of the form f(Δ) with f a function that is well-defined on the spectrum of the Laplacian, Δ.

Parameters:: f – The function

abstract invariant_covariance_function(spectral_variances: ndarray) → Callable[[ndarray], ndarray][source]¶: Returns a vectorized function that evaluates the two-point covariance of an invariant measure as a function of geodesic distance.

invariant_gaussian_measure(f: Callable[[float], float], /, *, expectation=None) → InvariantGaussianMeasure[source]¶

Returns a Gaussian measure on the space whose covariance takes the form f(Δ) with f a function that is well-defined on the spectrum of the Laplacian, Δ.

Parameters:

f – The function
expectation – The expected value for measure. Defaults to zero

Returns:

The measure as an instance of InvariantGaussianMeasure

l2_products_operator(weighting_functions: List[Any]) → LinearOperator[source]¶

Creates an operator that computes a vector of L2 inner products.

The action on a function ‘u’ returns a vector ‘d’ where d_i = <u, w_i>_L2.

Parameters:: weighting_functions – A list of vectors representing the weighting kernels.
Returns:: A mapping from this space to EuclideanSpace(len(weighting_functions)).
Return type:: LinearOperator

property laplacian: InvariantLinearAutomorphism¶: Returns the Laplacian, Δ, as an operator on the space. Note that this is an unbounded operator and it should be used with care.

abstract laplacian_eigenvalue(k: int | tuple[int, ...]) → float[source]¶

Returns the eigenvalue of the Laplacian for a given index.

The index k can be a single integer (e.g., for a circle) or a tuple of integers (e.g., for a sphere or torus), depending on the geometry of the space.

Parameters:: k – The index of the eigenvalue to return.

abstract laplacian_eigenvector_squared_norm(k: int | tuple[int, ...]) → float[source]¶

Returns the squared norm of the eigenvalue of the Laplacian for a given index.

The index k can be a single integer (e.g., for a circle) or a tuple of integers (e.g., for a sphere or torus), depending on the geometry of the space.

Parameters:: k – The index of the eigenvalue to return.

abstract laplacian_eigenvectors_at_point(x: float | tuple[float, ...]) → ndarray[source]¶

Returns a list of the values of the eigenvectors at a given point

Parameters:: x – The evaluation point.

norm_scaled_heat_kernel_gaussian_measure(scale: float, /, *, expectation: Vector = None, std: float = 1)[source]¶

Returns an invariant Gaussian measure with a heat kernel covariance equal to exp(-scale^2 * Δ).

Parameters:

scale – Scale parameter for the covariance.
expectation – The expected value for measure. Defaults to zero
std – The desired standard deviation

norm_scaled_invariant_gaussian_measure(f: Callable[[float], float], /, *, expectation=None, std=1.0) → InvariantGaussianMeasure[source]¶

Returns a Gaussian measure on the space whose covariance takes the form f(Δ) with f a function that is well-defined on the spectrum of the Laplacian, Δ.

The measure’s covariance is scaled such that:

E[||x||^2] = std*std + ||E[x]||^2

Parameters:

f – The function
expectation – The expected value for measure. Defaults to zero
std – The desired standard deviation

Returns:

The measure as an instance of InvariantGaussianMeasure

norm_scaled_sobolev_kernel_gaussian_measure(order: float, scale: float, /, *, expectation: Vector = None, std: float = 1)[source]¶

Returns an invariant Gaussian measure with a Sobolev-type covariance proportional to (1 + scale^2 * Δ)^-order.

The measure’s covariance is scaled such that the expected value for the samples norm is equal to the given standard deviation.

Parameters:

order – Order parameter for the covariance.
scale – Scale parameter for the covariance.
expectation – The expected value for measure. Defaults to zero
std – The desired standard deviation

pairs_within_distance(points: List[float | tuple[float, ...]], max_distance: float) → Tuple[ndarray, ndarray, ndarray][source]¶: Finds all pairs of points within a given geodesic distance. Returns arrays of row indices, column indices, and their geodesic distances.

abstract project_function(f: Callable[[float | tuple[float, ...]], float | ndarray]) → Vector[source]¶

Returns an element of the space by projecting a given function.

Parameters:: f – A function that takes a point and returns a value.

abstract random_point() → Any[source]¶: Returns a single random point from the underlying symmetric space.

random_points(n: int) → List[float | tuple[float, ...]][source]¶

Returns a list of n random points.

Parameters:: n – The number of random points to generate.

abstract representative_index(degree: int) → int | tuple[int, ...][source]¶: Returns a valid eigenfunction index for the given degree.

static sobolev_kernel(order: float, scale: float) → Callable[[float], float][source]¶: Returns a Sobolev kernel covariance function: (1 + scale^2 * k)^-order.

sobolev_kernel_gaussian_measure(order: float, scale: float, /, *, expectation: Vector = None)[source]¶

Returns an invariant Gaussian measure with a Sobolev-type covariance equal to (1 + scale^2 * Δ)^-order.

Parameters:

order – Order parameter for the covariance.
scale – Scale parameter for the covariance.
expectation – The expected value for measure. Defaults to zero

property spatial_dimension: int¶: The dimension of the symetric manifold.

abstract with_degree(degree: int | tuple[int, ...]) → SymmetricHilbertSpace[source]¶: Returns a new instance of the space with a modified truncation degree.

zero_operator(codomain: HilbertSpace | None = None) → InvariantLinearAutomorphism[source]¶

Returns the zero operator.

If the codomain is the space itself (or not provided), it returns an InvariantLinearAutomorphism for optimal arithmetic operations. Otherwise, it falls back to the standard abstract zero operator.

class pygeoinf.symmetric_space.SymmetricSobolevSpace(lebesgue_space: AbstractSymmetricLebesgueSpace, order: float, scale: float, /, *, safe: bool | None = True)[source]¶

Bases: MassWeightedHilbertModule, AbstractSymmetricLebesgueSpace

The Sobolev space Hˢ constructed over a symmetric Lebesgue space.

This implementation leverages the mass-weighting framework to ensure that the inner product and dual mappings correctly account for the smoothness order and scale.

Clusters a list of points into interacting blocks based on their geodesic distance.

This is particularly useful for generating the ‘interacting_blocks’ required by localized preconditioners in Bayesian inversions.

Parameters:

points – A list of points on the symmetric manifold.
threshold – The maximum geodesic distance between points in a cluster.
n_clusters – The exact number of clusters to form (alternative to threshold).
linkage_method – The hierarchical clustering method to use (‘complete’, ‘average’, ‘single’, etc.). Defaults to ‘complete’ to ensure clusters remain geographically compact.

Returns:

A list of lists, where each sub-list contains the indices of the points belonging to a specific cluster.

degree_multiplicity(degree: int) → int[source]¶: Returns the number of eigenfunctions associated with a given degree.

degree_transfer_operator(target_degree: int | tuple[int, ...]) → LinearOperator[source]¶: Returns the transfer operator from this space to one with a different degree.

dirac(point: float | tuple[float, ...]) → LinearForm[source]¶: Returns the Dirac measure at the chosen point.

dirac_representation(point) → Vector[source]¶: Returns the representation of the Dirac measure at the chosen point.

Builds a highly specialized, ultra-fast sparse preconditioner for point evaluation problems when the prior measure is invariant.

If max_distance is set to 0.0, it bypasses the distance calculations and instantly builds a purely diagonal (Jacobi) preconditioner.

Parameters:

prior_measure – The invariant prior Gaussian measure.
points – The list of observation points.
max_distance – The geodesic distance beyond which covariance is assumed to be zero (enforces sparsity). Set to 0.0 for a pure diagonal.
data_error_measure – The Gaussian measure describing the data noise. Can be None if the noise is already folded into the prior_measure.
apply_taper – If true, applies Gaspari-Cohn taper.
parallel – If True, computes the error measure diagonal in parallel.
n_jobs – Number of CPU cores to use if parallel=True.

Returns:

A LinearOperator representing the inverse of the approximated sparse normal matrix.

Constructs the exact 4th-order variable-coefficient flexure operator for a floating elastic plate mapped to a Riemannian manifold.

Mathematical Formulation:¶

This operator evaluates the covariant tensor divergence of the non-Euclidean moment equations. Under the geometer’s Laplacian convention, the exact equation is:

Op(w) = Delta(D Delta w)

Delta(D_{eff}) Delta w

tr(Hess(D_{eff}) Hess(w)) + 2 K nabla(D_{eff}) cdot nabla w

K D_{eff} Delta w

rho_g w

Notation:¶

$D_{eff} = D(1 - nu)$ is the effective flexural rigidity.
$K$ is the Gaussian curvature of the manifold.

Physical & Geometric Breakdown:¶

Principal Bending: Delta(D Delta w)
The standard biharmonic flexure term.
Scalar Property Interaction: - Delta(D_{eff}) Delta w
The scalar component of the Föppl-von Kármán bilinear bracket, accounting for gradients in the effective rigidity.
Tensor Gradients & Commutator: tr(Hess(D_{eff}) Hess(w)) + 2 K nabla(D_{eff}) cdot nabla w
Accounts for the coupling of spatial rigidity gradients with surface twist. The second half of this term (the gradient-curvature coupling) arises analytically from the non-commutativity of covariant derivatives. In this implementation, both terms are elegantly generated together in a coordinate-free manner by subtracting the Bochner identity block.
Geometric Covariant Softening: - K D_{eff} Delta w
Because covariant derivatives do not commute, pure bending on a curved surface inherently penalizes the manifold’s stiffness.
Foundation Restoring Force: rho_g w
The local hydrostatic restoring force (e.g., fluid mantle buoyancy).

param flexural_rigidity:: The flexural rigidity of the plate ($D$). Can be a spatially varying field (Vector) or a constant scalar.
type flexural_rigidity:: Union[Vector, float, int]
param poisson_ratio:: The Poisson’s ratio of the material ($nu$). Can be a spatially varying field (Vector) or a constant scalar. (Used to compute the effective rigidity $D_{eff}$).
type poisson_ratio:: Union[Vector, float, int]
param buoyancy_factor:: The hydrostatic foundation restoring force coefficient ($rho_g$). Can be a spatially varying field (Vector) or a constant scalar.
type buoyancy_factor:: Union[Vector, float, int]
returns:: A formally self-adjoint linear operator representing the exact 4th-order variable-coefficient flexural equations mapped to the manifold.
rtype:: LinearOperator

property gaussian_curvature: float¶: Returns the Gaussian curvature (K) of the symmetric manifold.

geodesic_ball_integral(center: Any, radius: float, /, *, n_points: int | None = None, normalize: bool = False) → LinearForm[source]¶

Return a linear form for a geodesic-ball integral.

The returned form represents $int_{B(center, radius)} u(x),dV(x)$. If normalize is true, it represents the average over the geodesic ball instead.

Subclasses with exact spectral formulas can override this method. Otherwise, a subclass-provided geodesic_ball_quadrature rule is used.

Parameters:

center – Center point of the geodesic ball.
radius – Geodesic radius in the manifold’s distance units.
n_points – Number of quadrature points requested from geodesic_ball_quadrature.
normalize – If true, divide by the quadrature measure of the ball.

geodesic_ball_quadrature(center: Any, radius: float, n_points: int) → Tuple[List[Any], ndarray][source]¶

Return quadrature points and weights for a geodesic ball.

geodesic_distance(p1: float | tuple[float, ...], p2: float | tuple[float, ...]) → float[source]¶

Returns the shortest distance along the manifold between two points.

Parameters:

p1 – The starting point.
p2 – The end point.

geodesic_integral(p1: float | tuple[float, ...], p2: float | tuple[float, ...], /, *, n_points: int | None = None) → LinearForm[source]¶

Returns a linear functional representing the line integral of a function along a geodesic path.

For Hilbert spaces with a specified scale, the method can automatically determine the required quadrature density to resolve the smooth features of the space’s sensitivity kernels.

Parameters:

p1 – The starting point of the geodesic. The type is manifold-dependent (e.g., float for Circle, tuple for Sphere).
p2 – The end point of the geodesic.
n_points (int, optional) – The number of Gauss-Legendre quadrature points. If None, it is heuristically determined as: $n = \lceil (\text{arc\_length} / \text{scale}) \times 2 \rceil$. This ensures at least two points per characteristic length-scale, providing stable sampling of the sensitivity kernel. Defaults to None.

Returns:

A linear functional whose action on a vector u computes: the approximated line integral.

Return type:

LinearForm

Raises:

NotImplementedError – If the Sobolev order $s$ is less than or equal to half the spatial dimension $n/2$.

geodesic_integral_representation(p1: float | tuple[float, ...], p2: float | tuple[float, ...], /, *, n_points: int | None = None) → Vector[source]¶

Returns the Riesz representation (sensitivity kernel) of the line integral.

Parameters:

p1 – Start and end points of the geodesic.
p2 – Start and end points of the geodesic.
n_points – Number of quadrature points.

geodesic_quadrature(p1: Any, p2: Any, n_points: int) → Tuple[List[Any], ndarray][source]¶

Returns quadrature points and weights for a geodesic between p1 and p2.

Returns:: List of manifold coordinates. weights: Integration weights scaled by the line element.
Return type:: points

index_to_integer(k: int | tuple[int, ...]) → int[source]¶: Maps an eigenvalue index to an integer.

integer_to_index(i: int) → int | tuple[int, ...][source]¶: Maps an integer to the eigenvalue index

invariant_covariance_function(spectral_variances: ndarray) → Callable[[ndarray], ndarray][source]¶: Returns a vectorized function that evaluates the two-point covariance of an invariant measure as a function of geodesic distance.

Constructs the inverse of the exact variable-coefficient flexure operator.

Eigenvalue(k) = D k^2 - K D_{eff} k + rho_g

Where $D_{eff} = D(1 - nu)$ and $K$ is the Gaussian curvature.

For varying properties, it utilizes a preconditioned Conjugate Gradient solver, using the spatially averaged properties to build the exact geometric spectral inverse as an ideal preconditioner.

Parameters:

flexural_rigidity (Union[Vector, float, int]) – The flexural rigidity of the plate ($D$).
poisson_ratio (Union[Vector, float, int]) – The Poisson’s ratio of the material ($nu$). (Used to compute the effective rigidity $D_{eff}$).
buoyancy_factor (Union[Vector, float, int]) – The hydrostatic foundation restoring force coefficient ($rho_g$).
baseline_rigidity (Optional[float], optional) – A constant scalar rigidity used to construct the geometric spectral preconditioner. If None, the spatial average of flexural_rigidity is computed dynamically.
baseline_buoyancy (Optional[float], optional) – A constant scalar foundation restoring force used to construct the geometric spectral preconditioner. If None, the spatial average of buoyancy_factor is computed dynamically.
solver (Optional[IterativeLinearSolver], optional) – An instance of an iterative linear solver used to invert the varying-coefficient operator. If None, a standard Conjugate Gradient solver (CGSolver) is instantiated.

Returns:

A linear operator representing the exact geometric spectral inverse (for constant properties) or the preconditioned iterative inverse (for varying properties).

Return type:

LinearOperator

l2_products_operator(weighting_functions: List[Any]) → LinearOperator[source]¶

Creates an operator that computes a vector of L2 inner products.

The action on a function ‘u’ returns a vector ‘d’ where d_i = <u, w_i>_L2.

Parameters:: weighting_functions – A list of vectors representing the weighting kernels.
Returns:: A mapping from this space to EuclideanSpace(len(weighting_functions)).
Return type:: LinearOperator

laplacian_eigenvalue(k: int | tuple[int, ...]) → float[source]¶

Returns the eigenvalue of the Laplacian for a given index.

The index k can be a single integer (e.g., for a circle) or a tuple of integers (e.g., for a sphere or torus), depending on the geometry of the space.

Parameters:: k – The index of the eigenvalue to return.

laplacian_eigenvector_squared_norm(k: int | tuple[int, ...]) → float[source]¶

Returns the squared norm of the eigenvalue of the Laplacian for a given index.

The index k can be a single integer (e.g., for a circle) or a tuple of integers (e.g., for a sphere or torus), depending on the geometry of the space.

Parameters:: k – The index of the eigenvalue to return.

laplacian_eigenvectors_at_point(x: float | tuple[float, ...]) → List[float | ndarray][source]¶

Returns a list of the values of the eigenvectors at a given point

Parameters:: x – The evaluation point.

property order: float¶: The Sobolev order.

order_inclusion_operator(target_order: float) → LinearOperator[source]¶: Returns the inclusion operator from this space to one of a lower order.

pairs_within_distance(points: List[float | tuple[float, ...]], max_distance: float) → Tuple[ndarray, ndarray, ndarray][source]¶: Finds all pairs of points within a given geodesic distance. Returns arrays of row indices, column indices, and their geodesic distances.

path_average_operator(paths, /, *, n_points=None)[source]¶

Constructs a tomographic operator mapping a function field to its line integrals along a set of geodesic paths.

Parameters:

paths (List[Tuple[Any, Any]]) – A list of start and end point pairs defining the geodesics.
n_points (int, optional) – The number of quadrature points per path. If None, the heuristic based on the Sobolev scale is used.

Returns:

An operator mapping Space -> EuclideanSpace(len(paths)).: The adjoint of this operator performs the ‘back-projection’ mapping data residuals into the function space.

Return type:

LinearOperator

point_evaluation_operator(points: List[Any]) → LinearOperator[source]¶

Returns a linear operator that evaluates a function at a list of points.

Parameters:: points – A list of points at which to evaluate the functions.

point_value_scaled_heat_kernel_gaussian_measure(scale: float, /, *, expectation: Vector = None, std: float = 1)[source]¶

Returns an invariant Gaussian measure with a heat-kernel covariance proportional to exp(-scale^2 * Δ).

The covariance of the operator is scaled such that the standard deviation of the point-wise values are equal to the given amplitude.

Parameters:

scale – Scale parameter for the covariance.
std – The desired standard deviation for the pointwise values.

point_value_scaled_invariant_gaussian_measure(f: Callable[[float], float], /, *, expectation: Vector = None, std: float = 1)[source]¶

Returns an invariant Gaussian measure with covariance proportional to f(Δ), where f must be such that this operator is trace-class.

The covariance of the operator is scaled such that the standard deviation of the point-wise values are equal to the given std.

Parameters:

f – A real-valued function that is well-defined on the spectrum of the Laplacian, Δ.
std – The desired standard deviation for the pointwise values.

Raises:

NotImplementedError – If the Sobolev order is less than n/2, with n the spatial dimension.

Notes

This method applies for symmetric spaces an invariant measures. As a result, the pointwise variance is the same at all points. Internally, a random point is chosen to carry out the normalisation.

point_value_scaled_sobolev_kernel_gaussian_measure(order: float, scale: float, /, *, expectation: Vector = None, std: float = 1)[source]¶

Returns an invariant Gaussian measure with a Sobolev-type covariance proportional to (1 + scale^2 * Δ)^-order.

The covariance of the operator is scaled such that the standard deviation of the point-wise values are equal to the given amplitude.

Parameters:

order – Order parameter for the covariance.
scale – Scale parameter for the covariance.
std – The desired standard deviation for the pointwise values.

random_point() → float | tuple[float, ...][source]¶: Returns a single random point from the underlying symmetric space.

random_source_receiver_paths(n_sources: int, n_receivers: int) → List[Tuple[Any, Any]][source]¶

Generates a list of source-receiver pairs by connecting every source to every receiver.

This method uses the existing random_points() logic to generate coordinates appropriate for the specific symmetric space. For a set of S sources and R receivers, this returns a list of S*R paths.

Parameters:

n_sources – The number of random source locations to generate.
n_receivers – The number of random receiver locations to generate.

Returns:

A list of tuples, where each tuple contains: a (source, receiver) pair.

Return type:

List[Tuple[Any, Any]]

representative_index(degree: int) → int | tuple[int, ...][source]¶: Returns a valid eigenfunction index for the given degree.

property safe: bool¶: True if mathematical checks are being applied.

property scale: float¶: The Sobolev length-scale.

sobolev_function(lambda_val) → float[source]¶: Returns the value of the Sobolev function associated with the space.

spatial_multiplication_operator(f: Vector) → LinearOperator[source]¶: Returns a LinearOperator corresponding to the mapping u -> f * u with f a given function that is suitably regular.

abstract with_order(order: float) → SymmetricSobolevSpace[source]¶: Returns a new instance of the exact same space but with a modified Sobolev order.