Native sparse balancing#

The native backend compiles the graph directly into a labelled, float64 sparse quadratic program and sends it to Gurobi through gurobipy. It avoids constructing Pyomo variables and expression trees. Native Gurobi is the default backend; Pyomo remains an explicit alternative and independent equation reference.

Install the package and provide a valid Gurobi license:

pip install -e packages/bonsai-io

Select the backend#

Use the normal registered workflow tasks and set options on Config:

from bonsai.io import Config, Workflow

config = Config()
# Defaults: native Gurobi, barrier, eight threads and a work limit.
workflow = Workflow(config=config)
workflow.execute(
    start="convert_sut_to_graph", stop="balance_sut", force_run=True
)

This example uses already assembled SUT accounts. Rebuild those accounts first when upstream supply/use data change. Rebuild the graph when the technical constraint CSV changes: the serialized graph carries the constraints consumed by balancing. Set balance_backend="pyomo" to use the existing Pyomo path. Neither backend is silently substituted for the other after a failure.

The remaining native settings are balance_residual_tolerance (default 2e-5) and balance_diagnostic_time_limit (default 60 seconds). Solver parameters in balance_solver_options override the numerical and work-limit defaults. The native path requires Gurobi; the legacy Pyomo path retains its Ipopt option.

The independent acceptance tolerance is separate from Gurobi’s FeasibilityTol, which remains 1e-6. The acceptance check divides each row violation by max(1, sum(abs(A_row) * abs(factors))) for the normalised matrix. Consequently, 2e-5 is not a universal percentage error in physical or monetary quantities. Original-unit residuals remain diagnostic, and the chosen acceptance tolerance is recorded in balance_diagnostics.json. Set balance_residual_tolerance=1e-6 explicitly to retain the stricter acceptance policy.

For a graph already available in memory:

balanced = graph.run_balance(
    backend="native", solver="gurobi", solver_options={"Threads": 1}
)

Responsibilities and numerical representation#

SUTGraph -> GraphBalanceCompiler -> SparseBalanceProblem
                                      |
                               GurobiSparseSolver
                                      |
                                BalanceSolution
                                      |
                            checked update of SUTGraph

GraphBalanceCompiler reads graph data and emits the current balance equations. It does not mutate the graph or write files. SparseBalanceProblem contains:

  • One column per supply/use flow layer, labelled by account, endpoints, layer and unit. Its variable is a nonnegative multiplicative adjustment factor.

  • CSR constraint rows with lower and upper bounds, labelled by rule family, node IDs and units. Human-readable names and locations are stored separately from compact solver names. Technical ratios retain both flow and reference units.

  • Initial quantities, variable bounds and diagonal objective weights. The objective remains sum(sqrt(abs(initial)) * (factor - 1)**2).

  • Row scales and aggregated skip counts. Existing TW, waste-service product, near-zero and extreme-imbalance exclusions are retained.

Each row and its bounds are always divided by its largest absolute coefficient. Zero-coefficient rows use scale one. There is no row-scaling option. The row scale is retained to report violations in the original row units. Scaling changes the interpretation of an absolute solver feasibility tolerance; inspect the independent residual report rather than relying on solver status alone. The original-unit maximum spans different units and is diagnostic, not a global mass or monetary residual.

The compiler preserves the legacy transfer rule’s reversed input-key test. Normal directed inputs therefore remain skipped and are counted under transfer_coefficient_constraint:legacy_reverse_key_absent. Fixing this existing indexing defect would enable additional constraints and is a separate modelling change, rather than an implicit consequence of selecting a new backend.

GurobiSparseSolver owns the solver environment and model for one solve. It uses bulk sparse constraints and a diagonal quadratic objective, without dense matrices, through Gurobi’s matrix API. It is not a persistent worker across separate workflow executions. No automatic warm-start speedup is claimed for Gurobi barrier.

A candidate with OPTIMAL or SUBOPTIMAL status must pass independent row and bound checks before application. Suboptimal acceptance additionally requires deterministic solver settings, is logged as a warning, and records accepted_suboptimal=True and optimality_proven=False. Replicability does not waive feasibility checks. Other termination states and missing candidates remain errors. All source flow bindings and quantities are checked before any graph mutation. Signed initial inventories retain their signs because factors remain nonnegative. Native export keeps validated small factors; it does not round factors below 1e-8 to zero after validation. Factors are saved on the existing Property.warm_start field for provenance; this does not activate warm starting. Small solver-tolerance excursions outside a variable bound are projected onto that bound only after the raw solution passes validation. The projected point is checked again, and the largest correction is recorded. This prevents small negative factors from reversing the sign of an inventory or supply flow.

Repeatability#

The compiler sorts node IDs, flow endpoints and layer keys before constructing the matrix. Its result does not depend on dictionary insertion order or Python hash randomisation. Rows are always normalised before solver ingestion. The production benchmark observed loss of small coefficients under this scaling; retain independent residual checks and inspect numerical diagnostics.

Native defaults are Method=2 (barrier), Threads=8, Seed=0 and WorkLimit=8000. There is no default wall-clock cutoff. Work units are not seconds; this changes the stopping budget from the earlier 8,000-second limit. Finite TimeLimit, Threads=0 and methods outside 0, 1, 2, 4 disqualify suboptimal acceptance. A caller can still supply these options, but such a suboptimal result is rejected.

Gurobi documents deterministic results for the same ordered model, parameters, hardware and solver version. Work limits provide deterministic stopping, whereas time limits can stop at different points under different machine load. See Gurobi’s determinism guidance. This is not a promise of identical results across versions or hardware.

Diagnostics retain the model fingerprint, SHA-256 of the checked factor vector, Gurobi version, options, work consumed and effective repeatability settings. Matching configuration is not an empirical repeat-run test. The regression suite checks exact matrices and cold-solve factors after graph reordering; full-data compilation checks and historical production optimisation results are described in the comparison report.

Diagnostics and failures#

The balance_sut task owns these files in the configured Balanced run directory:

  • gurobi_native.log: native solver log, unless overridden in solver options.

  • balance_diagnostics.json: dimensions, options, skip counts, phase timings, termination, objective and independent residual results.

  • infeasible_constraints.json: solver-computed IIS rows and variable bounds mapped back to graph labels. An empty list means no IIS was produced for this run; consult the diagnostic status rather than interpreting it as feasibility.

Numerical errors and ambiguous infeasibility may trigger a bounded diagnostic solve of the same feasible region with a zero objective and dual simplex. This is an infeasibility diagnostic, not a fallback balancing solution. Its result never authorizes balanced graph export. Time limits, unacceptable statuses and failed residual checks raise BalanceSolveError before balanced artifacts are written. The exception carries diagnostics and any IIS for programmatic callers; graph/compiler/solver helpers do not write those reports themselves.

Validation and rollout#

The first full pre-balance production run on 29 October 2026 did not produce an exportable native solution with explicit maximum-coefficient row scaling. Gurobi returned SUBOPTIMAL (status 13) after 139 barrier iterations and about 40 minutes of optimisation. The adapter rejected the candidate and wrote no balanced graph. This is not a proof of infeasibility: status 13 means that optimality tolerances were not satisfied; see the Gurobi status definitions. The full Pyomo run also returned status 13, but its legacy acceptance policy exported the candidate. The completed native control without extra row scaling also returned status 13 and was rejected. Independent audits found 428 native and 394 Pyomo rows above either recorded residual threshold. The full equations matched within floating-point precision, but neither candidate passed all checks. Native used 38.5% less peak committed memory and reached its stopping point in 28.8% less core time; these are not two certified solutions. Those measurements precede the user’s subsequent choice of native as the default and acceptance of repeatable suboptimal results. The earlier candidate is not automatically exportable: it failed the recorded residual checks. No full-production repeat solve with the new canonical ordering has been performed.

The regression suite compares every active row, variable bound and objective coefficient with the independent legacy Pyomo model. It also checks native solves, signed quantities, reference units, skipped rules, stale graph detection and labelled infeasibility reports. Gurobi integration tests require a license.

python -m pytest tests/test_sparse_balance.py

The earlier synthetic construction benchmark is evidence for a pilot, not a production speed guarantee. Keep the same input graph, objective and tolerances when comparing backends. Compare compilation/build/optimisation/export time and peak memory together with objective, residuals, changed flows, constraint coverage and downstream accounting. Do not overwrite an accepted production snapshot merely to benchmark the new implementation.

Implementation verification on 29 October 2026#

The implemented compiler also passed a read-only compilation of the complete current graph: 4,148,481 variables, 2,962,106 rows and 14,079,749 nonzeros. Compilation took 61.4 seconds; peak working memory, including the loaded graph, was 6.25 GiB. This check did not solve or overwrite the production snapshot.

Isolated solves used Python 3.12.10, NumPy 2.0.0, SciPy 1.15.3, Pyomo 6.9.2 and Gurobi 13.0.2. Both adapters used one thread, barrier, NumericFocus=3, ScaleFlag=2, feasibility/optimality tolerances of 1e-7 and BarConvTol=1e-10. The Pyomo comparison used its persistent Gurobi interface. These are individual pilot measurements, not repeated production benchmarks.

Input

Variables

Preparation: Pyomo / native

Optimisation: Pyomo / native

Peak working memory: Pyomo / native

Replicated recycling/market toy graph

30,000

4.72 / 1.10 s

0.69 / 0.76 s

399 / 295 MiB

Luxembourg activity slice

43,677

4.23 / 0.93 s

1.57 / 4.18 s

429 / 321 MiB

Preparation includes compilation/model construction and native solver transfer. Memory includes the input graph and independent checking. Native preparation was faster, while optimisation took longer on this real slice. Explicit row scaling and different row ordering can change solver iterations and timing even when the equations are mathematically equivalent.

The Luxembourg slice comes from the already-balanced snapshot and retains complete inputs/outputs of LU activities, including their existing technical bounds. It excludes other countries’ producers, so some boundary product balances are skipped. This is an integration check rather than a forecast for a new pre-balance solve. Every active row, variable bound and objective coefficient matched the legacy model. Both solutions passed independent checks, with scaled row violations below 1.2e-15; their objectives differed by 2.9e-12 in relative terms. The largest mass-flow difference was 0.0039 tonnes on an approximately 265,000-tonne flow.

Zero initial quantities have zero objective weight and an unidentified adjustment factor. Of 408 factors differing by more than 1e-5, 406 belonged to zero quantities, whose exported values remain zero. Compare balanced quantities, objective and residuals together rather than requiring identical factors in these flat directions. The regression suite includes this semidefinite case.