hpfracc.solvers — architecture, dependencies, and maintenance

This note complements CONTRIBUTING.md, ALGORITHMS_ARCHITECTURE.md, SPECIAL_ARCHITECTURE.md, and VALIDATION_ARCHITECTURE.md. It is the detailed map for fractional ODE / PDE / SDE time-stepping, noise models, graph–temporal coupling, and the JAX/Diffrax bridge.

---

1. Design goals

1. Time-stepping first: provide reference integrators for fractional-in-time dynamics (history convolutions, GL-type weights) that are usable from examples, validation scripts, and hpfracc.ml.neural_fsde.

2. Clear separation from “solver” in other stacks: the word solver also appears in SciPy (solve_ivp), torchdiffeq, Diffrax (dfx.diffeqsolve(..., solver=...)), and generic numerical texts. This package is hpfracc.solvers — not those APIs (see §6).

3. Optional acceleration without hard cycles: FFT / array backends may consult hpfracc.ml.intelligent_backend_selector from inside lazy helpers (_get_intelligent_selector), so importing hpfracc.solvers does not require ML training code until a backend decision runs.

4. Honest modeling docs: coupled_solvers.py ships an explicit module docstring on approximation limits (spatial Euler, temporal memory model, asymmetric coupling in splitting). Treat that text as normative for interpreting coupled outputs.

---

2. Module layout (file-by-file)

Module	Role	Primary entry points
ode_solvers.py	Fractional ODE stepping, FFT-accelerated history convolution, FixedStepODESolver, solve_fractional_ode, solve_fractional_system.	Uses FractionalOrder from core.definitions; Γ via special.gamma_beta.gamma_function (lazy _get_gamma_function with SciPy fallback).
sde_solvers.py	Fractional SDE drivers, FractionalSDESolver, Euler–Maruyama / Milstein-style paths, FastHistoryConvolution, SDESolution, solve_fractional_sde, solve_fractional_sde_system.	Imports noise_models.NoiseModel; same Γ pattern as ODE; optional intelligent FFT backend.
pde_solvers.py	Spatial–temporal fractional PDE prototypes (diffusion / advection / reaction–diffusion), sparse GL-style operators, solve_fractional_pde.	SciPy sparse linear algebra; Γ from scipy.special.gamma at module level (distinct from ODE/SDE adapter path — see §5).
noise_models.py	NoiseModel hierarchy: Brownian, fractional Brownian, Lévy, coloured noise; NoiseConfig, create_noise_model, generate_noise_trajectory.	Optional NumPyro / JAX for some paths (NUMPYRO_AVAILABLE gate).
coupled_solvers.py	Graph–SDE coupling: OperatorSplittingSolver, MonolithicSolver, solve_coupled_graph_sde, CoupledSolution.	Imports ode_solvers.gamma (module-level bound function) and FastHistoryConvolution from sde_solvers. Read the module docstring before publishing results.
solvers_jax_utils.py	JAX ecosystem bridge: Diffrax VirtualBrownianTree, optional precomputed fBm via lrdbenchmark, Optax wrapper stubs. get_virtual_brownian_path is the main export used from docs / pyproject optional extra.	Not imported from hpfracc.solvers.__init__ — import hpfracc.solvers.solvers_jax_utils explicitly when wiring Diffrax. Requires diffrax, jax; fBm path needs equinox + lrdbenchmark (pip install 'hpfracc[lrdbenchmark]').
__init__.py	Re-exports ODE / PDE / SDE / noise / coupled symbols; compatibility aliases and explicit placeholders for unimplemented high-order APIs (see §4).	Keeps legacy names (FractionalODESolver = FixedStepODESolver, etc.).

---

3. Import graph (intended edges)

flowchart TB
  subgraph solvers_pkg["hpfracc.solvers"]
    ODE["ode_solvers"]
    SDE["sde_solvers"]
    PDE["pde_solvers"]
    NOISE["noise_models"]
    COUPLED["coupled_solvers"]
    JAXU["solvers_jax_utils (explicit import)"]
  end

  CORE["core.definitions.FractionalOrder"]
  SPEC["special.gamma_beta (gamma_function)"]
  MLSEL["ml.intelligent_backend_selector (lazy)"]

  ODE --> CORE
  ODE --> SPEC
  SDE --> CORE
  SDE --> NOISE
  PDE --> CORE
  COUPLED --> ODE
  COUPLED --> SDE

  ODE -.->|first FFT backend use| MLSEL
  SDE -.->|first FFT backend use| MLSEL
  PDE -.->|first array backend use| MLSEL

Downstream of solvers

• hpfracc.ml.neural_fsde imports solve_fractional_sde / FractionalSDESolver from sde_solvers.

There is no solvers → neural_fsde import at module load; the dependency is one-way ML → solvers.

---

4. Public surface and compatibility policy (__init__.py)

Supported ODE workflow: FixedStepODESolver and solve_fractional_ode with documented methods (e.g. Euler, Adams–Bashforth, Runge–Kutta, predictor–corrector where implemented).

Aliases (same implementation):

• FractionalODESolver, AdaptiveFractionalODESolver, AdaptiveFixedStepODESolver → FixedStepODESolver (adaptive behaviour is not a separate class — controlled by solver parameters where applicable).

Explicit placeholders (raise NotImplementedError on construction or call):

• AdvancedFractionalODESolver, HighOrderFractionalSolver, VariableStepPredictorCorrector, solve_advanced_fractional_ode, solve_high_order_fractional_ode.

Thin wrappers:

• solve_predictor_corrector → delegates to solve_fractional_ode.

• PredictorCorrectorSolver, AdamsBashforthMoultonSolver → FixedStepODESolver.

When extending the package, prefer adding real implementations or narrowing __all__ in a major release over growing the placeholder list.

---

5. Gamma (Γ) usage consistency

Module	Γ source	Note
ode_solvers, sde_solvers	special.gamma_beta.gamma_function with SciPy fallback via _get_gamma_function	SciPy-first, matches SPECIAL_ARCHITECTURE.md “solvers use gamma_function”.
pde_solvers	scipy.special.gamma directly	Historical / simpler path for spatial coefficients; not the JAX-aware special.gamma. If you unify behaviour later, document BC impact.

---

6. Naming collisions (read before mixing ecosystems)

Term here	Often means elsewhere	Disambiguation
solve_fractional_ode	scipy.integrate.solve_ivp	Different formulation (fractional history, not classical y'=f(t,y) with dense output API).
FractionalSDESolver	Diffrax diffeqsolve	Diffrax solves SDEs/ODEs with its own stepper objects named “solver”; hpfracc fractional SDE classes live here.
FixedStepODESolver	torchdiffeq odeint / Neural ODE dopri5	hpfracc.ml.neural_ode may use torchdiffeq when installed; that is orthogonal to this ODE module unless you explicitly bridge them.
BenchmarkSuite	hpfracc.validation.method_benchmarks / hpfracc.benchmarks	Validation benchmarks compare methods; hpfracc.benchmarks is the numerical/ML benchmark package — see BENCHMARKS_ARCHITECTURE.md.

---

7. Optional dependencies (install surface)

Feature	Typical extras / packages
NumPyro noise paths	hpfracc[probabilistic] / numpyro (see noise_models guards).
Diffrax + fBm driver	diffrax, equinox, hpfracc[lrdbenchmark] for get_virtual_brownian_path when hurst ≠ 0.5.
Intelligent FFT / array backend hints	hpfracc.ml modules; failures disable the selector and fall back to NumPy/SciPy heuristics.

---

8. Tests (where to run what)

Tree	Focus
tests/test_solvers/	API, predictor–corrector, PDE/ODE “goldmine” style coverage.
tests/solvers/	Smaller grids, accuracy, system shapes, no-placeholder regressions.
tests/test_sde_solvers/	SDE correctness, noise models, coupled dynamics.

Example focused runs from repo root:

python -m pytest tests/test_solvers/ tests/solvers/ -q
python -m pytest tests/test_sde_solvers/ -q

For JAX-only utilities, add targeted tests (or subprocess import checks) if you change solvers_jax_utils — avoid assuming GPU.

---

9. Maintenance risks

Risk	Mitigation
__init__.py eager submodule imports	import hpfracc.solvers loads ODE, PDE, SDE, noise, coupled (heavy-ish). For minimal scripts, import hpfracc.solvers.ode_solvers (or the specific submodule) directly.
Placeholder class proliferation	Document in release notes when removing symbols; tests in test_no_placeholder_solvers guard expectations.
Coupled solver semantics	Do not “fix” silently: if physics changes, update module docstring, this file §1, and tests together.
ML selector import failures	Already caught; keep fail-soft behaviour so CI without full ML stack still runs ODE/SDE smoke tests.

---

10. Related documentation

• SPECIAL_ARCHITECTURE.md — Γ / binomial / Mittag–Leffler (hpfracc.special).

• ALGORITHMS_ARCHITECTURE.md — derivative engines vs legacy integrals.

• VALIDATION_ARCHITECTURE.md — analytical/convergence checks that may call solvers but stay in validation.

• BENCHMARKS_ARCHITECTURE.md — benchmark package vs validation method_benchmarks.

• Sphinx: solvers_api.rst, neural_ode_sde_api.rst (narrative SDE / neural material; legacy docs/sde_api_reference.rst is kept in-repo but excluded from the default HTML build).

---

11. Sphinx / Read the Docs

• hpfracc.solvers.solvers_jax_utils is intentionally not run through .. automodule:: in docs/api/solvers_api.rst: optional JAX / Diffrax stacks (and JAX version quirks such as primitive_jvps) make autodoc brittle on minimal builders. Treat this file and the module docstring as the canonical export map.

• sphinx -W (warnings as errors): the default HTML build is kept strict via docs/conf.py (narrow suppress_warnings, large auxiliary Markdown trees excluded from the source set, sectional docs/api/*.rst automodule directives marked :no-index: so the flat api_reference page owns the inventory, redundant automethod stubs removed from api_reference.rst, and a few docstrings fixed so docutils does not treat |…| as substitution markup). Re-run locally as python -m sphinx -b html docs docs/_build -W -q.