Skip to content

feat: pass the build solver to the interpolator before setup; add a native claudius script - #362

Merged
lachlangrose merged 56 commits into
Loop3D:masterfrom
lachlangrose:feat/native-claudius-script
Oct 28, 2026
Merged

lachlangrose merged 56 commits into
Loop3D:masterfrom
lachlangrose:feat/native-claudius-script

Conversation

@lachlangrose

Copy link
Copy Markdown
Member

Depends on #354. This branch is made from the head of #354, because it uses the native solver. Until #354 merges, this PR also shows the commits of #354. The change of this PR is only the last commit, 81f640d8.

Summary

  • Fix: GeologicalFeatureBuilder.build gave the solver build argument only to solve_system. Setup ran with the default solver. So with solver="native", the FDI setup still built the sparse rows, the memory saving was lost, and the solver logged solver="native" saves memory only when it is set before setup. The builder now sets interpolator.solver from the build argument before setup.
  • Script: scripts/native_claudius.py builds the claudius stratigraphy with the FDI interpolator and solver="native".
    • --solver both (the default) also builds the model with cg and compares the two fields.
    • --plot shows the 3D view. --screenshot PATH saves it to a PNG.
    • Other options: --nelements, --tol, --maxiter, --verbose.
  • Tests: tests/unit/modelling/test_builder_solver.py
    • A solver build argument reaches interpolator.solver.
    • With solver="native", the setup builds no sparse rows and logs no warning. This test skips when loop_interpolation_core is not installed.

Both tests fail without the fix and pass with it.

Results of the script

Claudius, 58 × 72 × 49 = 204,624 nodes:

native cg
Build time, tol 1e-8 7.4 – 10 s 41 s
Memory core peak 8.2 MB (40 bytes per node) 1,147,630 × 204,624 sparse matrix
Field difference, tol 1e-8 3.3e-6 of the range
Field difference, tol 1e-11 3.2e-9 of the range

The difference becomes smaller with the tolerance, so the two solvers find the same field.

Test plan

  • tests/unit/modelling/test_builder_solver.py: 2 passed (with the native core installed)
  • tests/unit: 700 passed, 7 failed. All 7 failures are in tests/unit/interpolator/test_2d_p1_p2_support.py, and the same 7 fail without this change.
  • scripts/native_claudius.py --solver native --verbose: no "saves memory" warning
  • ruff check on the new files

🤖 Generated with Claude Code

…tive extra

Wave 1 of spec 001: release-please entries, uv workspace exclude, Cargo and
maturin scaffolds, 20x20x20 reference fixture, native extra, and a failing
test for the native solver choice.
…native solver tests

Wave 2 of spec 001: copy cg, grid, stencil, interpolator and edge-case tests
from loop-rs, accept the native solver name, and add Python tests that skip
until the core is installed.
…y and D8 stubs

Baseline for the 100^3 reference problem R1 is 84.6 bytes per dof.
…nput and memory tracker tests

After table cuts, R1 uses 79.1 bytes per dof (baseline 84.6).
T013 also removes the Clone derive from SolveSession (approved by the user)
because Problem now holds the interpolator.
…lazy grid tables

Constraint record sizes follow the real 144-byte stored row. The scaling test
now keeps the row density constant (approved by the user). The table size test
builds the tables before it measures them.
R1 uses 63.3 bytes per dof through the session (target 75). Also adds the
session_r1 measurement test, makes every session_cg test take the allocator
lock, and implements Field::values_f64 for double precision.
… tests and the fixture parity test

The parity test fails with a difference of 24 percent of the range. T026
stays open until the mapping error is found.
…ure, border and gradient basis

Curvature rows are interior on all three axes, border rows are centred
second derivatives on layers 1 and n-2, and the gradient strike and dip use the
horizontal strike of the Python solver. The fixture parity test now passes.
Four ported unit tests that depended on the old one-sided border are updated.
…y pre-solve

ADMM stops only when the largest violation is at most 1e-4 of the field range
(SC-004). A pre-solve without bounds ends the solve when no bound is active.
Infeasible problems return the 20 worst rows after a stall test.
R3 uses 1.10 times the peak of R1 (limit 1.2). ADMM rebuilds the data
right-hand side each outer iteration and the constraint rows are reserved
before they are added. Adds the inequality scale tests.
…neric over Scalar

The native solver ignores column scaling (plan answer), so the column
scaling row is removed from the unsupported settings test.
The native tests use the same tolerance and iteration limit as the reference
solver. The inequality test uses a well-posed problem and turns column
scaling off for the Python ADMM reference, which does not scale the bounds.
…ive solver guide

Single precision uses 0.56 times the peak of double precision on the 50^3
test. The final residual is recomputed in double precision.
…e memory estimate

The estimate follows the measured peak within 12 percent from 10^3 to 200^3
(SC-009). R1 uses 63.2 MB (SC-001) and 250^3 in single precision uses 456 MB
(SC-008). The scaling test range is 5.5 to 8.5 because the shell scratch grows
with the surface (approved by the user).
…el feature

The wheel is about 3 times faster than the matrix solver on a 100^3 grid
(single thread: about 10 percent slower). The browser build stays single thread.
…ild script

The module is 92.5 KB after gzip (limit 2 MB). The fixture test and the
responsiveness tests pass in Chromium.
…b README

TypeScript is a dev dependency of the browser test package, so the type check
command of the task works.
… the solver comparison

The browser tests pass in Chromium.
The browser result equals the Rust result bit for bit. Memory 65.8 MB.
Set the version to 0.1.0 in the release manifest, Cargo.toml and
Cargo.lock. Add Cargo.toml to the extra files of the package, so that
release-please bumps it with pyproject.toml.

Review finding 9 (R010).
…ssion

When finalize returns an error, store it in the session. Later calls to
step, best_so_far and finish return the stored error. Before, a second
step could report a zero-length field as converged.

The f32 buffer is not reserved before the f64 residual check. That
would raise the f32 peak by 4 bytes per dof and break the 60 % peak
limit. The stored error is the guard for that failure point.

Review finding 7 (R006).
…rted input

- Reject interpolator classes that the core cannot represent. Only
  FiniteDifferenceInterpolator, or a subclass that sets
  supports_native = True in its own body, is accepted (finding 1).
- Reject sparse rows that setup did not make. Setup records the names of
  its rows; any other name, and add_inequality_feature, raise
  NotImplementedError (finding 1).
- A native setup sets _native_setup. A later solve with a matrix solver
  raises ValueError before the field is changed (finding 2).
- Clear solver_history and reset the timing on the native path
  (finding 10).
- Keep the old field when the native solve does not converge. The
  partial field is in latest_native_report["field"] (finding 8).
- Cast the field to float64 and check its length (finding 6).

Review findings 1, 2, 6, 8, 10 (R001-R005).
…ative adapter

- Drop data rows like the matrix path: a non-finite position, value or
  vector, a weight of 0 or NaN, and a gradient or tangent of length 0.
  A negative weight is sent as its absolute value, because least
  squares squares the weight (finding 3).
- Interface family: drop NaN ids and non-finite positions, then groups
  with fewer than two points. A negative ipw is used as abs(ipw)
  (finding 4).
- Inequality pairs: the new setup keyword inequality_pair_max_rows
  (default 100 000 for each pair) caps the rows. Above the cap, the
  adapter samples the pairs at random with a fixed seed
  (inequality_pair_seed, default 0) and logs a warning. Rows are built
  in chunks of 100 000 with index arithmetic, never the full cross
  product. Rows with a non-finite position are dropped (finding 5).

Review findings 3, 4, 5 (R007-R009).
…o solve_system

solve_system reads ridge_factor, add_ridge_regularisation, warm_start and
x0 from solver_kwargs. cg, lsmr and admm accept a start value. The start
value is read before _pre_solve clears the previous field. With column
scaling the start value is divided by the scaling diagonal.
InterpolatorBuilder.solve has ridge_factor and warm_start parameters.

The native path still receives the raw kwargs, so check_supported
rejects these keys until the core supports them.
…oblem

Problem.set_ridge_factor sets the ridge term of the normal operator. 0
removes it. Problem.set_initial_guess sets the start value of the solve.
CG starts with r = rhs - N x0, and the ADMM pre-solve starts from x0.
The residual stays relative to |rhs|. The memory estimate adds 8 bytes
for each dof when a start value is set.

The Python and wasm bindings expose both methods.
…tive solver

solve_system gives x0 (from x0 or warm_start) and the effective ridge
factor to _native.solve, which sends them to the core. The native
solver accepts ridge_factor, add_ridge_regularisation, warm_start and
x0. The user guide lists the new options.
Use as_chunks for fixed-size chunks (chunks_exact_to_as_chunks), and put
getter_with_clone only on the String and Vec fields of StepResult
(clone_on_copy). The JavaScript API does not change.
The core is not a workspace member, so the extra loop-interpolation[native]
had no source and uv could not resolve the workspace. A path source lets
uv read its static metadata. The core is built only when the extra is used.
…ields

macOS arm64 differs by about 1e-6 from the golden field because of the
summation order. The solver tolerance does not resolve 1e-6.
A fork pull request gets a read-only GITHUB_TOKEN, so the push of the
lint/style-fixes branch failed with 403. Push and master runs are not
changed.
…ance

The parallel sums of the core do not add in a fixed order, so two runs
can differ in the last bits (SC-003). On macOS the difference was 1.4e-10
and the bit-equal check failed. The sampled rows are still compared bit
for bit.
…s stencils

The row of the normal operator at a node depends only on the class of
each coordinate (first, second, third, interior, and the mirrored three),
so 343 kernels cover the curvature and border families. Each kernel is
built from the family stencils and weights. apply_normal_operator_into
does one sweep: edge nodes with their class kernel, interior runs as
contiguous multiply-add passes with the 25-tap kernel, then the data
rows. This replaces the fused-interior, shell-correction and border
passes and their scratch buffers. Grids with an axis of fewer than 7
nodes, and directional regularisation, keep the old path. The unfused
path stays the reference in the tests.

One apply on 100^3: 46.4 ms -> 7.7 ms on one thread, 40.8 ms -> 8.1 ms
with rayon. The peak is 40 bytes per dof (was about 53), and the memory
estimate has a class-stencil model (ratio 1.09 to 1.14 of the peak).
The f32 residual check needs no rhs vector on this path.

Rename the start value API to the spec names: Problem::set_start,
CgState::restart_from, Python set_start, wasm setStart. Add tests: the
class stencil equals the reference to 2e-15 on 9 grid shapes, a start
at the solution needs 0 iterations, and nested 10^3 -> 20^3 equals a
cold 20^3 solve within SC-003.
…MM solves

ADMM stops at its tolerance, so a different summation order can change
the field by about that tolerance. On macOS two runs differed by 2.5e-8,
above the 1e-9 limit of the previous commit. SC-003 allows 1e-6 of the
field range.
core-wheels.yml builds abi3 wheels for Linux x86_64 and aarch64, macOS
arm64 and x86_64 and Windows x64, and an sdist, with maturin-action. Each
wheel that the runner can run is installed and solves a small problem.
It runs on pull requests that change the crate, and pypi.yml calls it.

pypi.yml uploads the core with PyPI trusted publishing. No other package
waits for it. The core pyproject gets authors, licence, classifiers and
project URLs.
…ative claudius script

GeologicalFeatureBuilder.build gave the 'solver' build argument only to
solve_system. Setup ran with the default solver, so with solver='native'
the FDI setup still built the sparse rows, the memory saving was lost,
and the solver logged a warning. The builder now sets interpolator.solver
from the build argument before setup.

scripts/native_claudius.py builds the claudius stratigraphy with the
native solver, can compare it with cg, and can show or save the 3D view.
@lachlangrose
lachlangrose merged commit 5024e4f into Loop3D:master Oct 28, 2026
2 checks passed
@lachlangrose
lachlangrose deleted the feat/native-claudius-script branch October 28, 2026 09:47
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

1 participant