Solvers API Reference 

Mathematical basis:

conv(C, Y) = ifft(fft(C) * fft(Y))

Performance:

Direct summation: O(N²)
FFT-based: O(N log N)
Backend selection: < 0.001 ms overhead

hpfracc.solvers.ode_solvers._fast_history_sum(coeffs, f_hist, reverse=True, verbose=False)[source]

Compute weighted sum of history using FFT-based convolution.

This is the key optimization for fractional ODE solvers, replacing:: sum_{j=0}^{n} coeffs[n-j] * f_hist[j]

with an O(N log N) FFT-based approach instead of O(N²) direct summation.

Parameters:

coeffs (ndarray) – Coefficient array of length N
f_hist (ndarray) – History array of shape (N,) or (N, m) where m is state dimension
reverse (bool) – If True, apply coefficients in reverse order (typical for convolution)
verbose (bool)

Returns:

Weighted sum result (scalar or array of length m)

Return type:

class hpfracc.solvers.ode_solvers.FixedStepODESolver(derivative_type='caputo', method='predictor_corrector', adaptive=True, tol=1e-06, max_iter=1000, *, fractional_order=None, rtol=None, atol=None, order=1, min_h=None, max_h=None, min_step=None, max_step=None)[source]

Bases: object

Base class for fixed-step fractional ODE solvers.

Provides common functionality for solving fractional ordinary differential equations of the form:

D^α y(t) = f(t, y(t))

where D^α is a fractional derivative operator.

Parameters:

derivative_type (str)
method (str)
adaptive (bool)
tol (float)
max_iter (int)
fractional_order (float | FractionalOrder | None)
rtol (float | None)
atol (float | None)
order (int)
min_h (float | None)
max_h (float | None)
min_step (float | None)
max_step (float | None)

__init__(derivative_type='caputo', method='predictor_corrector', adaptive=True, tol=1e-06, max_iter=1000, *, fractional_order=None, rtol=None, atol=None, order=1, min_h=None, max_h=None, min_step=None, max_step=None)[source]

Initialize fractional ODE solver.

Parameters:

derivative_type (str) – Type of fractional derivative (“caputo”, “riemann_liouville”, “grunwald_letnikov”)
method (str) – Numerical method (“predictor_corrector”, “adams_bashforth”, “runge_kutta”)
adaptive (bool) – Use adaptive step size control
tol (float) – Tolerance for convergence
max_iter (int) – Maximum number of iterations
fractional_order (float | FractionalOrder | None)
rtol (float | None)
atol (float | None)
order (int)
min_h (float | None)
max_h (float | None)
min_step (float | None)
max_step (float | None)

_validate_alpha(alpha)[source]

Validate the fractional order.

Parameters:: alpha (float | FractionalOrder)

solve(f, t_span, y0, alpha, h=None, **kwargs)[source]

Solve fractional ODE.

Parameters:

f (Callable) – Right-hand side function f(t, y)
t_span (Tuple[float, float]) – Time interval (t0, tf)
y0 (float | ndarray) – Initial condition(s)
alpha (float | FractionalOrder) – Fractional order
h (float | None) – Step size (None for adaptive)
**kwargs – Additional solver parameters

Returns:

Tuple of (t_values, y_values)

Return type:

_solve_predictor_corrector(f, t0, tf, y0, alpha, h, **kwargs)[source]

Solve using Adams-Bashforth-Moulton predictor-corrector for Caputo fractional ODEs.

Based on the Volterra integral formulation: y(t) = y_0 + (1/Γ(α)) ∫_0^t (t-τ)^(α-1) f(τ, y(τ)) dτ

Reference: Diethelm, K. (2010). The Analysis of Fractional Differential Equations.

_solve_adams_bashforth(f, t0, tf, y0, alpha, h, **kwargs)[source]

Solve using Adams-Bashforth method.

Parameters:

f (Callable) – Right-hand side function
t0 (float) – Initial time
tf (float) – Final time
y0 (float | ndarray) – Initial condition
alpha (float | FractionalOrder) – Fractional order
h (float) – Step size
**kwargs – Additional parameters

Returns:

Tuple of (t_values, y_values)

Return type:

_solve_runge_kutta(f, t0, tf, y0, alpha, h, **kwargs)[source]

Solve using fractional Runge-Kutta method.

Parameters:

f (Callable) – Right-hand side function
t0 (float) – Initial time
tf (float) – Final time
y0 (float | ndarray) – Initial condition
alpha (float | FractionalOrder) – Fractional order
h (float) – Step size
**kwargs – Additional parameters

Returns:

Tuple of (t_values, y_values)

Return type:

_solve_euler(f, t0, tf, y0, alpha, h, **kwargs)[source]

Solve using fractional Euler method.

Parameters:

f (Callable) – Right-hand side function
t0 (float) – Initial time
tf (float) – Final time
y0 (float | ndarray) – Initial condition
alpha (float | FractionalOrder) – Fractional order
h (float) – Step size
**kwargs – Additional parameters

Returns:

Tuple of (t_values, y_values)

Return type:

_euler_step(f, t_values, y_values, n, alpha, h)[source]

Minimal fractional Euler update for 0<α≤1.

y_n = y_{n-1} + h^α / Γ(α+1) * f(t_{n-1}, y_{n-1})

Parameters:

f (Callable)
t_values (ndarray)
y_values (ndarray)
n (int)
alpha (float | FractionalOrder)
h (float)

Return type:

_compute_fractional_coefficients(alpha, N)[source]

Compute fractional derivative coefficients.

Parameters:

alpha (float | FractionalOrder) – Fractional order
N (int) – Number of time steps

Returns:

Array of coefficients

Return type:

hpfracc.solvers.ode_solvers.solve_fractional_ode(f, t_span, y0, alpha, derivative_type='caputo', method='predictor_corrector', adaptive=False, h=None, **kwargs)[source]

Solve fractional ODE.

Parameters:

f (Callable)
t_span (Tuple[float, float])
y0 (float | ndarray)
alpha (float | FractionalOrder)
derivative_type (str)
method (str)
adaptive (bool)
h (float | None)

Return type:

hpfracc.solvers.ode_solvers.solve_fractional_system(f, t_span, y0, alpha, derivative_type='caputo', method='predictor_corrector', **kwargs)[source]

Solve system of fractional ODEs.

Parameters:

f (Callable) – Right-hand side function f(t, y)
t_span (Tuple[float, float]) – Time interval (t0, tf)
y0 (ndarray) – Initial conditions
alpha (float | ndarray) – Fractional orders (scalar or array)
derivative_type (str) – Type of fractional derivative
method (str) – Numerical method
**kwargs – Additional solver parameters

Returns:

Tuple of (t_values, y_values)

Return type:

PDE Solvers

Fractional Partial Differential Equation Solvers

This module provides comprehensive solvers for fractional PDEs including finite difference methods, spectral methods, and adaptive mesh refinement.

Performance Note (v2.1.0): - Intelligent backend selection for sparse matrix operations - Optimal array operations based on problem size

hpfracc.solvers.pde_solvers._get_intelligent_selector()[source]: Get intelligent backend selector instance for PDE solvers.

hpfracc.solvers.pde_solvers._select_array_backend(data_size, operation_type='element_wise')[source]

Select optimal backend for array operations in PDE solving.

Parameters:

data_size (int) – Number of elements to process
operation_type (str) – Type of operation (e.g., “element_wise”, “matrix_multiply”)

Returns:

“numpy”, “numba”, or “jax”

Return type:

Backend name

class hpfracc.solvers.pde_solvers.FractionalPDESolver(pde_type='diffusion', method='finite_difference', spatial_order=2, temporal_order=1, adaptive=False, *, fractional_order=None, boundary_conditions=None)[source]

Bases: object

Base class for fractional PDE solvers.

Provides common functionality for solving fractional partial differential equations of various types.

Parameters:

pde_type (str)
method (str)
spatial_order (int)
temporal_order (int)
adaptive (bool)
fractional_order (float | FractionalOrder | None)
boundary_conditions (str | None)

__init__(pde_type='diffusion', method='finite_difference', spatial_order=2, temporal_order=1, adaptive=False, *, fractional_order=None, boundary_conditions=None)[source]

Initialize fractional PDE solver.

Parameters:

pde_type (str) – Type of PDE (“diffusion”, “advection”, “reaction_diffusion”)
method (str) – Numerical method (“finite_difference”, “spectral”, “finite_element”)
spatial_order (int) – Order of spatial discretization
temporal_order (int) – Order of temporal discretization
adaptive (bool) – Use adaptive mesh refinement
fractional_order (float | FractionalOrder | None)
boundary_conditions (str | None)

_validate_orders(alpha, beta)[source]

Validate fractional orders for PDE solver.

Parameters:

alpha (float)
beta (float)

solve(t_span, x_span, initial_condition, boundary_conditions, alpha, beta, **kwargs)[source]

Solve the fractional PDE.

Parameters:

t_span (Tuple[float, float])
x_span (Tuple[float, float])
initial_condition (Callable)
boundary_conditions (Dict)
alpha (float)
beta (float)

Return type:

Dict[str, ndarray]

class hpfracc.solvers.pde_solvers.FractionalDiffusionSolver(method='finite_difference', spatial_order=2, temporal_order=1, derivative_type='caputo')[source]

Bases: FractionalPDESolver

Solver for fractional diffusion equations.

Solves equations of the form: ∂^α u/∂t^α = D ∂^β u/∂x^β + f(x, t, u)

where α and β are fractional orders.

Parameters:

method (str)
spatial_order (int)
temporal_order (int)
derivative_type (str)

__init__(method='finite_difference', spatial_order=2, temporal_order=1, derivative_type='caputo')[source]

Initialize fractional diffusion solver.

Parameters:

method (str) – Numerical method
spatial_order (int) – Order of spatial discretization
temporal_order (int) – Order of temporal discretization
derivative_type (str) – Type of fractional derivative

solve(x_span, t_span, initial_condition, boundary_conditions, alpha, beta, diffusion_coeff=1.0, source_term=None, nx=100, nt=100, **kwargs)[source]

Solve fractional diffusion equation.

Parameters:

x_span (Tuple[float, float]) – Spatial interval (x0, xf)
t_span (Tuple[float, float]) – Time interval (t0, tf)
initial_condition (Callable) – Initial condition u(x, 0)
boundary_conditions (Tuple[Callable, Callable]) – Boundary conditions (left_bc, right_bc)
alpha (float | FractionalOrder) – Temporal fractional order
beta (float | FractionalOrder) – Spatial fractional order
diffusion_coeff (float) – Diffusion coefficient
source_term (Callable | None) – Source term f(x, t, u)
nx (int) – Number of spatial grid points
nt (int) – Number of temporal grid points
**kwargs – Additional solver parameters

Returns:

Tuple of (t_values, x_values, solution)

Return type:

_solve_spatial_step(solution, n, x_values, t_values, alpha, beta, source_term, nx, nt, pde_params)[source]

Solve spatial problem at current time step.

Parameters:

solution (ndarray) – Solution array
n (int) – Current time step
x_values (ndarray) – Spatial grid
t_values (ndarray) – Time grid
alpha (float | FractionalOrder) – Temporal fractional order
beta (float | FractionalOrder) – Spatial fractional order
source_term (Callable | None) – Source term
nx (int) – Number of spatial grid points
nt (int) – Number of temporal grid points
pde_params (dict) – Dictionary of parameters (diffusion_coeff, dx, dt)

Returns:

Solution at interior points

Return type:

_finite_difference_step(solution, n, x_values, t_values, alpha, beta, source_term, nx, nt, pde_params)[source]

Finite difference step.

Parameters:

solution (ndarray) – Solution array
n (int) – Current time step
x_values (ndarray) – Spatial grid
t_values (ndarray) – Time grid
alpha (float | FractionalOrder) – Temporal fractional order
beta (float | FractionalOrder) – Spatial fractional order
source_term (Callable | None) – Source term
nx (int) – Number of spatial grid points
nt (int) – Number of temporal grid points
pde_params (dict) – Dictionary of parameters (diffusion_coeff, dx, dt)

Returns:

Solution at interior points

Return type:

_spectral_step(solution, n, x_values, t_values, alpha, beta, source_term, nx, nt, pde_params)[source]: Correct, efficient, diagonal spectral solve for periodic BCs.

_compute_temporal_derivative(solution, n, alpha, dt, **kwargs)[source]

Computes the RHS vector of the system, based on previous time steps.

Parameters:

solution (ndarray)
n (int)
alpha (float | FractionalOrder)
dt (float)

Return type:

_compute_spatial_derivative(u, beta, dx)[source]

Compute spatial fractional derivative.

Parameters:

u (ndarray) – Solution at current time
beta (float | FractionalOrder) – Spatial fractional order
dx (float) – Spatial step size

Returns:

Spatial derivative

Return type:

_grunwald_letnikov_coeffs(order, n_points)[source]

Compute Grünwald-Letnikov coefficients.

Parameters:

order (float)
n_points (int)

Return type:

_validate_orders(alpha, beta)[source]

Validate fractional orders for PDE solver.

Parameters:

alpha (float)
beta (float)

_build_spatial_matrix(nx, beta, diffusion_coeff, dx, dt, alpha, side='lhs')[source]

Build spatial discretization matrix.

Parameters:

nx (int) – Number of spatial points
beta (float | FractionalOrder) – Spatial fractional order
diffusion_coeff (float) – Diffusion coefficient
dx (float) – Spatial step size
dt (float) – Temporal step size
alpha (float | FractionalOrder) – Temporal fractional order
side (str) – ‘lhs’ for left-hand side of the equation, ‘rhs’ for right-hand side

Returns:

Sparse matrix for spatial discretization

Return type:

spmatrix

_get_spatial_operator(nx, beta, dx)[source]

Computes the matrix for the Grünwald-Letnikov spatial derivative.

Parameters:

nx (int)
beta (float | FractionalOrder)
dx (float)

Return type:

class hpfracc.solvers.pde_solvers.FractionalAdvectionSolver(method='finite_difference', spatial_order=2, temporal_order=1, derivative_type='caputo')[source]

Bases: FractionalPDESolver

Solver for fractional advection-diffusion equations.

Solves equations of the form:: D_t^α u = v * D_x^β u

Currently, only integer orders (alpha=1, beta=1) are supported.

Parameters:

method (str)
spatial_order (int)
temporal_order (int)
derivative_type (str)

__init__(method='finite_difference', spatial_order=2, temporal_order=1, derivative_type='caputo')[source]

Initialize fractional PDE solver.

Parameters:

pde_type – Type of PDE (“diffusion”, “advection”, “reaction_diffusion”)
method (str) – Numerical method (“finite_difference”, “spectral”, “finite_element”)
spatial_order (int) – Order of spatial discretization
temporal_order (int) – Order of temporal discretization
adaptive – Use adaptive mesh refinement
derivative_type (str)

solve(t_span, x_span, initial_condition, boundary_conditions, alpha, beta, **kwargs)[source]

Solve the fractional advection equation.

Parameters:

t_span (Tuple[float, float])
x_span (Tuple[float, float])
initial_condition (Callable)
boundary_conditions (Dict)
alpha (float)
beta (float)

Return type:

Dict[str, ndarray]

_build_spatial_matrix(nx, beta, diffusion_coeff, dx, dt, alpha, side='lhs')[source]

Build spatial discretization matrix.

Parameters:

nx (int) – Number of spatial points
beta (float | FractionalOrder) – Spatial fractional order
diffusion_coeff (float) – Diffusion coefficient
dx (float) – Spatial step size
dt (float) – Temporal step size
alpha (float | FractionalOrder) – Temporal fractional order
side (str) – ‘lhs’ for left-hand side of the equation, ‘rhs’ for right-hand side

Returns:

Sparse matrix for spatial discretization

Return type:

spmatrix

_get_spatial_operator(nx, beta, dx)[source]

Computes the matrix for the Grünwald-Letnikov spatial derivative.

Parameters:

nx (int)
beta (float | FractionalOrder)
dx (float)

Return type:

class hpfracc.solvers.pde_solvers.FractionalReactionDiffusionSolver(method='finite_difference', spatial_order=2, temporal_order=1, derivative_type='caputo')[source]

Bases: FractionalPDESolver

Solver for fractional reaction-diffusion equations.

Solves equations of the form: ∂^α u/∂t^α = D ∂^β u/∂x^β + R(u) + f(x, t, u)

where α and β are fractional orders and R(u) is the reaction term.

Parameters:

method (str)
spatial_order (int)
temporal_order (int)
derivative_type (str)

__init__(method='finite_difference', spatial_order=2, temporal_order=1, derivative_type='caputo')[source]

Initialize fractional reaction-diffusion solver.

Parameters:

method (str) – Numerical method
spatial_order (int) – Order of spatial discretization
temporal_order (int) – Order of temporal discretization
derivative_type (str) – Type of fractional derivative

solve(x_span, t_span, initial_condition, boundary_conditions, alpha, beta, diffusion_coeff=1.0, reaction_term=None, source_term=None, nx=100, nt=100, **kwargs)[source]

Solve fractional reaction-diffusion equation.

Parameters:

x_span (Tuple[float, float]) – Spatial interval (x0, xf)
t_span (Tuple[float, float]) – Time interval (t0, tf)
initial_condition (Callable) – Initial condition u(x, 0)
boundary_conditions (Tuple[Callable, Callable]) – Boundary conditions (left_bc, right_bc)
alpha (float | FractionalOrder) – Temporal fractional order
beta (float | FractionalOrder) – Spatial fractional order
diffusion_coeff (float) – Diffusion coefficient
reaction_term (Callable | None) – Reaction term R(u)
source_term (Callable | None) – Source term f(x, t, u)
nx (int) – Number of spatial grid points
nt (int) – Number of temporal grid points
**kwargs – Additional solver parameters

Returns:

Tuple of (t_values, x_values, solution)

Return type:

hpfracc.solvers.pde_solvers.solve_fractional_diffusion(x_span, t_span, initial_condition, boundary_conditions, alpha, beta, diffusion_coeff=1.0, source_term=None, method='finite_difference', nx=100, nt=100, **kwargs)[source]

Solve fractional diffusion equation.

Parameters:

x_span (Tuple[float, float]) – Spatial interval (x0, xf)
t_span (Tuple[float, float]) – Time interval (t0, tf)
initial_condition (Callable) – Initial condition u(x, 0)
boundary_conditions (Tuple[Callable, Callable]) – Boundary conditions (left_bc, right_bc)
alpha (float | FractionalOrder) – Temporal fractional order
beta (float | FractionalOrder) – Spatial fractional order
diffusion_coeff (float) – Diffusion coefficient
source_term (Callable | None) – Source term f(x, t, u)
method (str) – Numerical method
nx (int) – Number of spatial grid points
nt (int) – Number of temporal grid points
**kwargs – Additional solver parameters

Returns:

Tuple of (t_values, x_values, solution)

Return type:

hpfracc.solvers.pde_solvers.solve_fractional_advection(x_span, t_span, initial_condition, boundary_conditions, alpha, beta, velocity=1.0, source_term=None, method='finite_difference', nx=100, nt=100, **kwargs)[source]

Solve fractional advection equation.

Parameters:

x_span (Tuple[float, float]) – Spatial interval (x0, xf)
t_span (Tuple[float, float]) – Time interval (t0, tf)
initial_condition (Callable) – Initial condition u(x, 0)
boundary_conditions (Tuple[Callable, Callable]) – Boundary conditions (left_bc, right_bc)
alpha (float | FractionalOrder) – Temporal fractional order
beta (float | FractionalOrder) – Spatial fractional order
velocity (float) – Advection velocity
source_term (Callable | None) – Source term f(x, t, u)
method (str) – Numerical method
nx (int) – Number of spatial grid points
nt (int) – Number of temporal grid points
**kwargs – Additional solver parameters

Returns:

Tuple of (t_values, x_values, solution)

Return type:

hpfracc.solvers.pde_solvers.solve_fractional_reaction_diffusion(x_span, t_span, initial_condition, boundary_conditions, alpha, beta, diffusion_coeff=1.0, reaction_term=None, source_term=None, method='finite_difference', nx=100, nt=100, **kwargs)[source]

Solve fractional reaction-diffusion equation.

Parameters:

x_span (Tuple[float, float]) – Spatial interval (x0, xf)
t_span (Tuple[float, float]) – Time interval (t0, tf)
initial_condition (Callable) – Initial condition u(x, 0)
boundary_conditions (Tuple[Callable, Callable]) – Boundary conditions (left_bc, right_bc)
alpha (float | FractionalOrder) – Temporal fractional order
beta (float | FractionalOrder) – Spatial fractional order
diffusion_coeff (float) – Diffusion coefficient
reaction_term (Callable | None) – Reaction term R(u)
source_term (Callable | None) – Source term f(x, t, u)
method (str) – Numerical method
nx (int) – Number of spatial grid points
nt (int) – Number of temporal grid points
**kwargs – Additional solver parameters

Returns:

Tuple of (t_values, x_values, solution)

Return type:

hpfracc.solvers.pde_solvers.solve_fractional_pde(x_span, t_span, initial_condition, boundary_conditions, alpha, beta, equation_type='diffusion', **kwargs)[source]

Generic solver for fractional PDEs.

Parameters:

x_span (Tuple[float, float]) – Spatial interval (x0, xf)
t_span (Tuple[float, float]) – Time interval (t0, tf)
initial_condition (Callable) – Initial condition u(x, 0)
boundary_conditions (Tuple[Callable, Callable]) – Boundary conditions (left_bc, right_bc)
alpha (float | FractionalOrder) – Temporal fractional order
beta (float | FractionalOrder) – Spatial fractional order
equation_type (str) – Type of PDE (“diffusion”, “advection”, “reaction_diffusion”)
**kwargs – Additional solver parameters

Returns:

Tuple of (t_values, x_values, solution)

Return type:

SDE Solvers

Fractional Stochastic Differential Equation Solvers

This module provides comprehensive solvers for fractional SDEs including various numerical methods, adaptive step size control, and error estimation.

Performance Note: - Uses FFT-based convolution for O(N log N) history summation instead of O(N²) - Intelligent backend selection for optimal performance - Multi-backend support (PyTorch, JAX, NumPy/Numba)

hpfracc.solvers.sde_solvers._get_gamma_function()[source]: Get gamma function through adapter system.

hpfracc.solvers.sde_solvers._get_intelligent_selector()[source]: Get intelligent backend selector instance.

hpfracc.solvers.sde_solvers._fft_convolution(coeffs, values, axis=0)[source]

Fast convolution using FFT for O(N log N) performance.

Parameters:

coeffs (ndarray) – Coefficient array (1D)
values (ndarray) – Value array (can be 1D or 2D)
axis (int) – Axis along which to perform convolution

Returns:

Convolution result (same shape as values)

Return type:

class hpfracc.solvers.sde_solvers.SDESolution(t, y, fractional_order, method, drift_func, diffusion_func, metadata=None)[source]

Bases: object

Solution object for SDE solvers.

Parameters:

t (ndarray)
y (ndarray)
fractional_order (float | FractionalOrder)
method (str)
drift_func (Callable)
diffusion_func (Callable)
metadata (Dict[str, Any])

t: ndarray

y: ndarray

fractional_order: float | FractionalOrder

method: str

drift_func: Callable

diffusion_func: Callable

metadata: Dict[str, Any] = None

__init__(t, y, fractional_order, method, drift_func, diffusion_func, metadata=None)

Parameters:

t (ndarray)
y (ndarray)
fractional_order (float | FractionalOrder)
method (str)
drift_func (Callable)
diffusion_func (Callable)
metadata (Dict[str, Any] | None)

Return type:

None

class hpfracc.solvers.sde_solvers.FractionalSDESolver(fractional_order, definition='caputo', adaptive=False, rtol=1e-05, atol=1e-08)[source]

Bases: ABC

Base class for fractional SDE solvers.

A fractional SDE takes the form:

D^α X(t) = f(t, X(t)) dt + g(t, X(t)) dW(t)

where:

α is the fractional order
f is the drift function
g is the diffusion function
W(t) is a Wiener process

Parameters:

fractional_order (float | FractionalOrder)
definition (str)
adaptive (bool)
rtol (float)
atol (float)

__init__(fractional_order, definition='caputo', adaptive=False, rtol=1e-05, atol=1e-08)[source]

Initialize fractional SDE solver.

Parameters:

fractional_order (float | FractionalOrder) – Fractional order (0 < α < 2)
definition (str) – Type of fractional derivative (“caputo” or “riemann_liouville”)
adaptive (bool) – Use adaptive step size
rtol (float) – Relative tolerance for adaptive stepping
atol (float) – Absolute tolerance for adaptive stepping

abstract solve(drift, diffusion, x0, t_span, **kwargs)[source]

Solve fractional SDE.

Parameters:

drift (Callable[[float, ndarray], ndarray]) – Drift function f(t, x) -> R^d
diffusion (Callable[[float, ndarray], ndarray]) – Diffusion function g(t, x) -> R^(d x m) where m is dimension of noise
x0 (ndarray) – Initial condition
t_span (Tuple[float, float]) – Time interval (t0, tf)
**kwargs – Additional solver parameters

Returns:

SDESolution object containing trajectory

Return type:

class hpfracc.solvers.sde_solvers.FastHistoryConvolution(alpha, num_steps, dim)[source]

Bases: object

Helper class for efficient history convolution in fractional SDE solvers.

Currently uses optimized dot product (O(N^2)), but structured to support FFT-based block convolution (O(N log N)) in future updates.

Parameters:

alpha (float)
num_steps (int)
dim (int)

__init__(alpha, num_steps, dim)[source]

Parameters:

alpha (float)
num_steps (int)
dim (int)

update(value)[source]

Add new value to history.

Parameters:: value (ndarray)

convolve()[source]

Compute convolution of weights with history.

Return type:: ndarray

class hpfracc.solvers.sde_solvers.FractionalEulerMaruyama(*args, **kwargs)[source]

Bases: FractionalSDESolver

Fractional Euler-Maruyama method for solving fractional SDEs.

This is a first-order strong convergence method.

__init__(*args, **kwargs)[source]

Initialize fractional SDE solver.

Parameters:

fractional_order – Fractional order (0 < α < 2)
definition – Type of fractional derivative (“caputo” or “riemann_liouville”)
adaptive – Use adaptive step size
rtol – Relative tolerance for adaptive stepping
atol – Absolute tolerance for adaptive stepping

solve(drift, diffusion, x0, t_span, num_steps=100, seed=None, **kwargs)[source]

Solve fractional SDE using Euler-Maruyama method.

Parameters:

drift (Callable) – Drift function f(t, x)
diffusion (Callable) – Diffusion function g(t, x)
x0 (ndarray) – Initial condition
t_span (Tuple[float, float]) – Time interval (t0, tf)
num_steps (int) – Number of time steps
seed (int | None) – Random seed for Wiener process

Returns:

SDESolution object

Return type:

class hpfracc.solvers.sde_solvers.FractionalMilstein(*args, **kwargs)[source]

Bases: FractionalSDESolver

Fractional Milstein method for solving fractional SDEs.

This is a second-order strong convergence method with improved accuracy.

__init__(*args, **kwargs)[source]

Initialize fractional SDE solver.

Parameters:

fractional_order – Fractional order (0 < α < 2)
definition – Type of fractional derivative (“caputo” or “riemann_liouville”)
adaptive – Use adaptive step size
rtol – Relative tolerance for adaptive stepping
atol – Absolute tolerance for adaptive stepping

solve(drift, diffusion, x0, t_span, num_steps=100, seed=None, **kwargs)[source]

Solve fractional SDE using Milstein method.

Parameters:

drift (Callable) – Drift function f(t, x)
diffusion (Callable) – Diffusion function g(t, x)
x0 (ndarray) – Initial condition
t_span (Tuple[float, float]) – Time interval (t0, tf)
num_steps (int) – Number of time steps
seed (int | None) – Random seed for Wiener process

Returns:

SDESolution object

Return type:

hpfracc.solvers.sde_solvers.solve_fractional_sde(drift, diffusion, x0, t_span, fractional_order=0.5, method='euler_maruyama', num_steps=100, seed=None, **kwargs)[source]

Solve a fractional SDE.

Parameters:

drift (Callable) – Drift function f(t, x) -> R^d
diffusion (Callable) – Diffusion function g(t, x) -> R^(d x m)
x0 (ndarray) – Initial condition
t_span (Tuple[float, float]) – Time interval (t0, tf)
fractional_order (float | FractionalOrder) – Fractional order (0 < α < 2)
method (str) – Solver method (“euler_maruyama”, “milstein”, “predictor_corrector”)
num_steps (int) – Number of time steps
seed (int | None) – Random seed for Wiener process
**kwargs – Additional solver parameters

Returns:

SDESolution object

Return type:

Example

>>> def drift(t, x):
...     return -0.5 * x
>>> def diffusion(t, x):
...     return 0.2 * np.eye(1)
>>> x0 = np.array([1.0])
>>> sol = solve_fractional_sde(drift, diffusion, x0, (0, 1), 0.5, num_steps=100)
>>> print(sol.y[-1])

hpfracc.solvers.sde_solvers.solve_fractional_sde_system(drift, diffusion, x0, t_span, fractional_order, method='euler_maruyama', noise_type='additive', num_steps=100, seed=None, **kwargs)[source]

Solve a system of coupled fractional SDEs.

Parameters:

drift (Callable) – Drift function f(t, x) -> R^d
diffusion (Callable) – Diffusion function g(t, x) -> R^(d x m)
x0 (ndarray) – Initial condition
t_span (Tuple[float, float]) – Time interval (t0, tf)
fractional_order (float | FractionalOrder | List[float]) – Fractional order(s) for system
method (str) – Solver method
noise_type (str) – Type of noise (“additive” or “multiplicative”)
num_steps (int) – Number of time steps
seed (int | None) – Random seed
**kwargs – Additional parameters

Returns:

SDESolution object

Return type:

Noise Models

Stochastic Noise Models for Fractional SDEs

This module provides various noise models including Brownian motion, fractional Brownian motion, Lévy processes, and coloured noise.

class hpfracc.solvers.noise_models.NoiseModel[source]

Bases: ABC

Base class for stochastic noise models.

abstract generate_increment(t, dt, size=(), seed=None)[source]

Generate noise increment.

Parameters:

t (float) – Current time
dt (float) – Time step
size (Tuple[int, ...]) – Output shape
seed (int | None) – Random seed

Returns:

Noise increment array

Return type:

class hpfracc.solvers.noise_models.BrownianMotion(scale=1.0)[source]

Bases: NoiseModel

Standard Brownian motion (Wiener process).

Generates independent Gaussian increments with variance dt.

Parameters:: scale (float)

__init__(scale=1.0)[source]

Initialize Brownian motion.

Parameters:: scale (float) – Scaling factor for noise

generate_increment(t, dt, size=(), seed=None)[source]

Generate Wiener increment.

Parameters:

t (float)
dt (float)
size (Tuple[int, ...])
seed (int | None)

Return type:

variance(dt)[source]

Get variance of increment.

Parameters:: dt (float)
Return type:: float

class hpfracc.solvers.noise_models.FractionalBrownianMotion(hurst=0.5, scale=1.0)[source]

Bases: NoiseModel

Fractional Brownian motion (fBm).

A Gaussian process with long-range dependence characterized by the Hurst exponent H (0 < H < 1).

Parameters:

hurst (float)
scale (float)

__init__(hurst=0.5, scale=1.0)[source]

Initialize fractional Brownian motion.

Parameters:

hurst (float) – Hurst exponent (H=0.5 gives standard Brownian motion)
scale (float) – Scaling factor for noise

generate_increment(t, dt, size=(), seed=None)[source]

Generate fractional Brownian motion increment.

Note: This is a simplified implementation. For full fBm, need to account for covariance structure.

Parameters:

t (float)
dt (float)
size (Tuple[int, ...])
seed (int | None)

Return type:

property is_standard_bm: bool: Check if this is standard Brownian motion.

class hpfracc.solvers.noise_models.LevyNoise(alpha=1.5, beta=0.0, scale=1.0, location=0.0)[source]

Bases: NoiseModel

Lévy noise for jump diffusions.

Uses stable distributions to model heavy-tailed noise.

Parameters:

alpha (float)
beta (float)
scale (float)
location (float)

__init__(alpha=1.5, beta=0.0, scale=1.0, location=0.0)[source]

Initialize Lévy noise.

Parameters:

alpha (float) – Stability parameter (0 < α ≤ 2)
beta (float) – Skewness parameter (-1 ≤ β ≤ 1)
scale (float) – Scale parameter
location (float) – Location parameter

generate_increment(t, dt, size=(), seed=None)[source]

Generate Lévy noise increment.

Uses stable distribution sampling (simplified implementation).

Parameters:

t (float)
dt (float)
size (Tuple[int, ...])
seed (int | None)

Return type:

class hpfracc.solvers.noise_models.ColouredNoise(correlation_time=1.0, amplitude=1.0, seed=None)[source]

Bases: NoiseModel

Coloured noise (Ornstein-Uhlenbeck process).

Gaussian noise with exponential autocorrelation.

Parameters:

correlation_time (float)
amplitude (float)
seed (int | None)

__init__(correlation_time=1.0, amplitude=1.0, seed=None)[source]

Initialize coloured noise.

Parameters:

correlation_time (float) – Correlation time constant
amplitude (float) – Noise amplitude
seed (int | None) – Random seed for state initialization

generate_increment(t, dt, size=(), seed=None)[source]

Generate coloured noise increment.

Parameters:

t (float)
dt (float)
size (Tuple[int, ...])
seed (int | None)

Return type:

reset()[source]: Reset the noise process state.

class hpfracc.solvers.noise_models.NoiseConfig(noise_type='brownian', hurst=0.5, scale=1.0, alpha=1.5, beta=0.0, correlation_time=1.0, amplitude=1.0)[source]

Bases: object

Configuration for noise models.

Parameters:

noise_type (str)
hurst (float)
scale (float)
alpha (float)
beta (float)
correlation_time (float)
amplitude (float)

noise_type: str = 'brownian'

hurst: float = 0.5

scale: float = 1.0

alpha: float = 1.5

beta: float = 0.0

correlation_time: float = 1.0

amplitude: float = 1.0

__init__(noise_type='brownian', hurst=0.5, scale=1.0, alpha=1.5, beta=0.0, correlation_time=1.0, amplitude=1.0)

Parameters:

noise_type (str)
hurst (float)
scale (float)
alpha (float)
beta (float)
correlation_time (float)
amplitude (float)

Return type:

None

hpfracc.solvers.noise_models.create_noise_model(config)[source]

Create a noise model from configuration.

Parameters:: config (NoiseConfig) – Noise configuration
Returns:: NoiseModel instance
Return type:: NoiseModel

hpfracc.solvers.noise_models.generate_noise_trajectory(noise_model, t_span, num_steps, size=(), seed=None)[source]

Generate a complete noise trajectory.

Parameters:

noise_model (NoiseModel) – Noise model to use
t_span (Tuple[float, float]) – Time interval (t0, tf)
num_steps (int) – Number of steps
size (Tuple[int, ...]) – Shape of noise increments
seed (int | None) – Random seed

Returns:

Tuple of (time array, noise increments array)

Return type: