Validated 2D particle solvers for soil–fluid interaction analysis.
This is the scientific core behind
ParticlesInteraction
(the interactive sandbox), rebuilt validation-first: SI units end-to-end,
no rendering, and a benchmark suite (validate.py) that checks every solver
against analytical solutions and published experimental data. Nothing here is
tuned "to look right" — if a behaviour matters, there is a benchmark for it.
particlesci.fluid.PBFluid— 2D Position-Based Fluids (Macklin & Müller 2013) in SI units (m, s, kg/m of depth; y-up). Jacobi constraint solver with explicit under-relaxation, one-sided (compression-only) density constraint, CFL-limited internal sub-stepping at ≤1/240 s, self-calibrated CFM epsilon from the rest lattice, XSPH smoothing. Static boundary particles (e.g. a grain bed) enter the density constraint.particlesci.granular.Granular— 2D grains (pymunk) in SI units: metre-scale collision slop, Coulomb friction, rolling resistance as exponential angular damping (mu_roll, the standard circular-DEM surrogate).
python3 -m pip install -r requirements.txt
python3 validate.py # full suite (a few minutes, CPU)
python3 validate.py --quick # coarser resolution
python3 validate.py --only hydrostatic| # | Benchmark | Reference | What is checked |
|---|---|---|---|
| 1 | hydrostatic |
analytical | interior density = ρ₀ (incompressibility), volume conservation, liquid at rest stays at rest, linear pressure profile with depth |
| 2 | dambreak |
Martin & Moyce (1952) experiments | dimensionless front position Z = x/a at T = t·√(2g/a) = 2 |
| 3 | repose |
Coulomb friction | poured-pile angle of repose ≈ arctan(μ) |
| 4 | infiltration |
Ergun equation | superficial velocity of water through a static coarse-grain bed |
Current results: 6/7 checks pass — recorded in
docs/VALIDATION.md — including the failures and what
they mean. Exit code = number of failed checks, so the suite is CI-ready.
Documented in detail in docs/VALIDATION.md; the headline lessons from
building the suite:
- PBD stiffness is timestep-dependent. At a 1/60 s solver step the hydrostatic interior density error is ~10 %; at 1/240 s it is ~1 %. The solver therefore sub-steps internally at ≤ 1/240 s.
- Jacobi constraint projection needs under-relaxation (ω ≈ 0.5): both ends of every pair are corrected simultaneously, and full-strength corrections diverge. (The predecessor project achieved the same effect accidentally through an oversized ε.)
- The density constraint must be one-sided (compression only): surface and wall particles have truncated kernel support, and a two-sided constraint turns that deficiency into spurious attraction.
- ✅ SI solvers + validation harness (this)
- Calibrated coupling: boundary pseudo-mass and drag derived from grain geometry (Ergun closure), not tuned
- GPU port (Taichi): DFSPH + proper DEM contacts, 10⁵–10⁶ particles, 3D
- Science modules, each with its own benchmark: erosion/corrosion, sediment transport (Shields), unsaturated flow (Green–Ampt)
- Experiment configs, parameter sweeps, reproducible reports; interactive front-end reattached from ParticlesInteraction
particlesci/
├── particlesci/
│ ├── fluid.py # PBF solver (SI)
│ └── granular.py # granular solver (SI)
├── validate.py # benchmark suite (exit code = #fails)
├── scripts/
│ └── calibrate.py # parameter-sweep experiments behind the defaults
└── docs/
└── VALIDATION.md # current benchmark results & analysis