Architecture

liulu enforces a strict separation between physics (astrophysical content) and numerics (pure numerical routines), wired together by a single driver. Respecting this split is the main way the codebase stays extensible — violating it is the main way to break it.

liulu/
├── model.py            # StreamingModel: the driver (wires physics + numerics)
├── physics/            # astrophysical content only
│   ├── base.py             # ABCs: RealSpaceCorrelation, PairwiseVelocityMoments, VelocityPDF
│   ├── real_space.py       # xi(r) providers: TabulatedXi, LinearTheoryXi, HalofitXi
│   ├── pairwise_velocity.py# TabulatedMoments (r-t pairwise velocity moments)
│   ├── velocity_pdf.py     # GaussianPDF, SkewTPDF, NIGPDF, GeneralizedHyperbolicPDF, TabulatedGHPDF
│   ├── los_projection.py   # project_to_los: r-t moments -> LOS moments (geometry)
│   └── cosmology.py        # Cosmology(FlatLambdaCDM) helper
├── numerics/           # no cosmology knowledge; operates on callables & arrays
│   ├── interpolation.py    # log_linear / log_spline / loglog interpolators
│   └── hankel.py           # pk_to_xi, xi_to_pk
└── io/                 # readers / writers for tables

The two design rules

Physics layer — no raw numerics

Classes in liulu/physics/ expose __call__ / evaluate and contain astrophysical content only. No np.trapz/scipy.integrate, no FFTs, no grid construction. New physics models are new classes inheriting the ABCs in physics/base.py.

Numerics layer — no cosmology

Functions in liulu/numerics/ are pure numerical routines operating on callables and arrays. They take signatures like f, a, b, tol — never a Cosmology or a moments object, and they have no knowledge of what \(\xi\) or \(\mathcal{P}\) mean.

The driver StreamingModel is the only place the two layers meet, via constructor injection: you build the physics components, hand them to the model, and the model applies the numerical integration strategy.

Data flow

 user supplies                physics layer                 driver
┌───────────────┐   ┌──────────────────────────────┐   ┌──────────────┐
│ xi_real       │──▶│                              │   │              │
│ (RealSpace…)  │   │                              │   │              │
│               │   │ PairwiseVelocityMoments      │   │ StreamingModel│
│ moments       │──▶│   --project_to_los-->        │──▶│  integrates   │──▶ xi^s(s_perp,s_par)
│ (Pairwise…)   │   │   LOS moments                │   │  over y       │    multipoles xi_0,2,4
│               │   │   --VelocityPDF.__call__-->  │   │              │
│               │   │   P(v_los | r_perp, r_par)   │   │              │
└───────────────┘   └──────────────────────────────┘   └──────────────┘

1. You provide a real-space xi_real and a PairwiseVelocityMoments.
2. The physics layer projects radial-tangential moments to the line of sight via project_to_los.
3. A VelocityPDF evaluates \(\mathcal{P}(v_\text{los}\mid r_\perp,r_\parallel)\).
4. StreamingModel integrates \([1+\xi(r)]\,\mathcal{P}\) over \(y\) and projects to multipoles.

Extending the package

• New \(\xi(r)\) source → subclass RealSpaceCorrelation.
• New velocity model → subclass VelocityPDF (and/or PairwiseVelocityMoments).
• New numerical strategy (interpolation, quadrature) → add to numerics/ with a callable/array signature, then inject it.

project_to_los is not an ABC — it is geometry (a binomial expansion), not a model choice. Extend it by adding higher-order \(c_{ij}\) methods to the moments class, never by hand-unrolling per-order formulae.