Release Notes

0.9.0

Built-in Radau IIA(5) DAE solver, SuiteSparse KLU as the default sparse linear-solver backend, and a silent-correctness fix for LoopEqn Jacobians that multiply a Var by a TimeSeriesParam.

New

  • Radau — 3-stage 5th-order Radau IIA fully-implicit Runge-Kutta DAE integrator. Ported from SciML’s RadauIIA5 (OrdinaryDiffEqFIRK) and adapted to Solverz’s nDAE API. Stiffly accurate (y_{n+1} = y_n + Z_3), with simplified Newton iteration in the W-coordinate system that decouples each step into one real and one complex sparse LU factorisation. Adaptive step size uses the Hairer / SciML predictive (Gustafsson) controller. Exposed as Solverz.Radau.

    from Solverz import Radau, Opt
    sol = Radau(ndae, [0.0, 20.0], y0, Opt(rtol=1e-6, atol=1e-9))
    

    Suitable for stiff or oscillatory DAEs that need higher order than Rodas. Dense output and terminal/non-terminal event detection both run on the Radau IIA(5) collocation polynomial, so event times are resolved to machine precision regardless of how large the controller lets h grow on smooth trajectories.

Changed

  • SuiteSparse KLU is now the default sparse-LU backend. Newton / Rosenbrock / BDF solves factor sparse Jacobians through a pure-Python ctypes binding to SuiteSparse KLU, with scipy SuperLU as the automatic fallback when libklu is missing or the matrix is complex or dense. KLU is roughly a factor of two faster than SuperLU on power-system / IEGS network Jacobians and needs no build step. Override per solve with Opt(linsolver='superlu'), for a script with set_linsolver('superlu'), for a block with with linsolver('superlu'):, or globally with SOLVERZ_LINSOLVER=superlu. KLU requires SuiteSparse installed (brew install suite-sparse or conda install -c conda-forge suitesparse). To reproduce older SuperLU-only results exactly, force SuperLU, since KLU can change adaptive step sequences at the tolerance level.

Fixed

  • LoopEqn Jacobian silently baked TimeSeriesParam contributions. The matrix-derivative constant classifier in the LoopEqn Jacobian pipeline (is_constant_matrix_deri) did not inspect whether a Para free symbol was a TimeSeriesParam. Any LoopEqn body multiplying a Var by a TimeSeriesParam factor (e.g. a short-circuit G_shunt[i] * ux[i] term) produced a Jacobian whose F_ correctly evaluated get_v_t(t) but whose J_ froze at the build-time TimeSeriesParam value and never updated, corrupting every short-circuit / leak / valve-trip / step-change simulation built that way. The classifier is now TimeSeriesParam-aware, and loop_jac_to_solverz_expr falls through to a per-cell dynamic kernel that re-evaluates the TimeSeriesParam every J call.

  • Saturation / In now preserve input shape under @njit. The Numba-compiled forms could return a scalar for a length-1 vector input; they now return an array matching the input shape.

0.8.6

Mat_Mul fallback diagnostic warnings (#125). When a Mat_Mul placeholder or a mutable Jacobian block is forced onto the slower scipy.sparse fallback path, the module printer now emits a UserWarning that tells the user what expression broke the fast path, why, and how to rewrite it.

New

  • Layer 1 Mat_Mul placeholder fallback warnings. The module printer’s print_F now emits one UserWarning per fallback placeholder, with specialised messages for negation, scalar multiplication, sum-of-matrices, multi-argument Mat_Mul folding, sparse dim=2 Param references in the operand, and demotion cascades (root-cause traceback).

  • Layer 2 mutable Jacobian block fallback warnings. When a Jacobian block term doesn’t match the Diag / row-scale / col-scale / biscale fast-path shapes, the module printer emits a UserWarning naming the equation, variable, and offending term.

  • matrix_calculus.md diagnostic warnings section listing every warning category, what triggers it, and the canonical fix.

Changed

  • Warnings only fire in module mode (not inline), because inline mode has no fast/fallback split. The existing _warn_dense_matmul_params warning (from FormJac time) still fires in both modes.

Fixed

  • Mat_Mul printer associativity. When the matrix argument was a non-atomic expression such as A + B or 2 * A, Mat_Mul’s numpy / sympystr / octave printers emitted A + B@x instead of (A + B)@x, so Mat_Mul(A+B, x) lambdified to a 2-D ndarray and failed Array(..., dim=1) validation in create_instance(). Each printer now wraps any non-Symbol/Function operand in parentheses.

  • Layer 2 classifier no longer conflates structurally distinct fallback shapes. The previous classifier short-circuited on any Diag-bearing term and unconditionally emitted "multiple Diag nodes", misdiagnosing single-Diag and element-wise Mul cases. The new classifier dispatches on the post-sign-extraction core and produces distinct messages for biscale, single-Diag, no-Diag Mat_Mul, element-wise Mul, bare Para, and other shapes.

  • Layer 1 Add+non-Para suggestion previously contained literal Python source (' + '.join(f'Mat_Mul({arg}, <operand>)' for arg in ...)) instead of a usable distributed-form rewrite.

How to silence

import warnings
warnings.simplefilter('ignore', UserWarning)

0.8.5

LoopEqn prototype close-out: new Set primitive for subset iteration (#129), Phase J2 translator Patterns 1 / 2 / 4 landings (#133), safety guards (#130, #132), and a full LoopEqn documentation rewrite.

New

  • Pyomo-style Set primitive (#129). Solverz.Set (internal class IndexSet) replaces the three-object subset-iteration pattern — auxiliary Param(..., dtype=int) + bounded Idx + every m.Var[m.subset_param[i]] indirection — with a one-line declaration:

    m.PVPQ = Set('PVPQ', pv_pq_arr)
    m.Bus  = Set('Bus', nb)
    i_p = m.PVPQ.idx('i_p')
    j   = m.Bus.idx('j')
    body_P = m.Vm[i_p] * Sum(m.Vm[j] * m.Gbus[i_p, j] * ..., j) + ...
    m.P_eqn = LoopEqn('P_eqn', outer_index=i_p, body=body_P, model=m)
    

    The body rewriter inserts the indirect gather (m.Vm[PVPQ[i_p]]) only when the target array is strictly larger than the set — subset-aligned storage whose length equals the set is indexed bare. Same indirect-outer path the sparsity analyzer already handles, no engine changes required. Partial- unknown Vars (Pyomo’s Var(m.PQ) pattern) are tracked separately in #136.

  • Model.add name-collision detection (#130). Raises ValueError when merging a sub-model whose attribute would overwrite an existing Param / Var / Eqn with a non-equal value. Previously the clash was silent — a real IES integration bug (gas- vs heat-network pipe_from_node parameters silently overwriting each other, producing |F| 1.5e6 residuals) motivated the guard. Identical values (shared object or value-equal Params like a common Cp) still merge without error.

  • LoopEqn + inline / partial-Jacobian path guard (#132, C2 guard from #128). made_numerical and Equations.FormPartialJac now raise NotImplementedError when the equation system contains any LoopEqn instance. LoopEqn kernels emit Python for-loops that are 3–7× slower than the lambdify-vectorised legacy Eqn path without Numba JIT; the guard redirects users to Solverz.module_printer(..., jit=True), which is LoopEqn’s design target.

Performance

  • Phase J2 translator Patterns 1 / 2 / 4 (#133 completion). LoopEqn Jacobian blocks of the following shapes now land on the constant-matrix fast path with zero per-Newton-iteration J cost:

    • Pattern 1 DiagSelectTerm — Diag(u) @ SelectMat @ Para.

    • Pattern 2 bilinear mixing — Diag(u) @ Para @ Diag(v) (new biscale term, flat / nested / -1-wrapped forms).

    • Pattern 4 Sum-KD — identity and non-identity column maps, bare and indirect outer. All three cases keep is_constant_matrix_deri = True so Value0 is baked into _data_ at module-build time.

  • Pattern 3 identity-map indirect KD → DiagTerm (#133). A top-level KD(diff, map[outer]) * Sum(...) with identity map (and n_outer == n_diff) folds into the existing direct-KD Diag(Mat_Mul(Para, iVar)) branch.

  • Mutable-matrix analyzer extensions_LoopJacSelectMat accepted as matrix operand alongside Para, including Mat_Mul(SelectMat, Para) composites; variadic Mat_Mul chains accepted in _classify_matmul; _extract_sign_and_core strips any numeric scalar coefficient.

  • Measured impact (Big IES, 8 996-DOF, Rodas3 Mode II) — wall 13.76 s vs 52.25 s for legacy (3.80× faster), F-eval 709.8 µs vs 4 348.5 µs, J-eval 4 634 µs vs 15 832 µs. Step counts identical (757 accepted) confirming trajectory match.

Documentation

  • docs/src/loopeqn.md — full rewrite as a beginner walkthrough. 10 sections from motivation to troubleshooting; two Cookbook worked examples (polar power-flow, integrated energy system) instead of synthetic snippets; Set-based subset iteration; the module_printer(jit=True) requirement; LoopOde; performance numbers.

  • docs/src/loopeqn_translator.md (new) — developer-facing appendix that defines every translator jargon term (KD, Mat_Mul, Phase J1/J2/J3, CSR walker, SelectMat, biscale) at first use and houses the Phase J2 coverage table, the supported-body-shape table, and the not-supported list.

  • docs/src/conf.py — enable myst_enable_extensions = ['dollarmath', 'amsmath', 'colon_fence'] so $$...$$ and :::{math} blocks render typeset in .md pages.

  • API reference gains IndexSet autoclass entry.

Internal

  • SymbolExtractor leaf-symbol fix — Para / iVar / iAliasVar used as indices now register in SymInIndex by name, so lambdified bodies that gather via a Param (e.g. iVar('Ts')[Para('fht_nl')]) no longer raise NameError at runtime.

  • Inline printer no longer registers LoopEqn.NUM_EQN or LoopEqnDiff._kernel_func_name callables — that code path is unreachable behind the guard and has been pruned.

  • Value0 shape is verified against (n_outer, n_diff) before Pattern 4 emits its translated Param expression, preventing a latent IndexError when a sparse-pattern row falls outside the LoopEqn’s equation-block range.

0.8.4

Never released. Content merged into 0.8.5.

0.8.3

Tooling-only release. No Python source-code, test, or runtime behaviour changes. The PyPI wheel for 0.8.3 is bit-identical to 0.8.2 modulo the new in-tree skill files (which are not bundled into the wheel).

Tooling

  • New Claude Code skill at .claude/skills/solverz-modeling/ that teaches Claude how to use Solverz for symbolic modeling and numerical simulation. Bundles:

    • SKILL.md — 4-step workflow (equation type → build → compile → solve), Var / Param / Eqn / Ode idioms, inline vs module_printer decision, every built-in solver with when-to-use notes, Mat_Mul fast vs fallback path with the full rewrite table, common pitfalls table, and a quick-reference card.

    • references/ecosystem.md — chapter map of the Solverz Cookbook, every reusable block in SolMuseum (gt, pv, st, eb, eps_network, heat_network, gas_network, pde.heat, pde.gas), and every helper in SolUtil (PowerFlow, DhsFlow, GasFlow, DhsFaultFlow).

    • references/examples/ — six canonical end-to-end runnable examples covering AE / DAE / FDAE / mutable Jacobian / events / AliasVar / Mat_Mul / model.add() composition:

      • bouncing-ball.md — minimal DAE with event handling

      • power-flow.md — canonical AE with Mat_Mul (case30, rectangular coordinates)

      • heat-flow.md — AE with mutable-matrix Jacobian

      • m3b9-dynamics.md — DAE with TimeSeriesParam fault scenario

      • gas-characteristics.md — FDAE with AliasVar (method of characteristics)

      • integrated-energy-system.md — multi-domain DAE composition using all three SolUtil flow solvers + the major SolMuseum DAE / AE blocks via model.add(...)

    • README.md — install command (symlink) and sync rule for contributors.

    Install on a contributor’s machine:

    ln -sfn "$(pwd)/.claude/skills/solverz-modeling" \
            ~/.claude/skills/solverz-modeling
    

    After the symlink is in place, git pull on a Solverz checkout updates the skill content automatically — no re-install step. Inside a Solverz checkout itself the skill auto-loads even without the global symlink, because Claude Code reads <cwd>/.claude/skills/.

The skill is not bundled into the PyPI wheel — it lives only in the source tree. Pip-installed users who want the skill should clone the repo separately.

The contribution guide in .claude/skills/solverz-modeling/README.md documents the sync rule for Solverz contributors: when a PR changes Solverz’s public API, the same PR should update the relevant skill files. Reviewers will check both.

0.8.2

Documentation-only release. No source-code, test, or behaviour changes. The PyPI wheel for 0.8.2 is bit-identical to 0.8.1 modulo the in-tree documentation.

The release rewrites the Matrix-Vector Calculus chapter to address eight reviewer comments about misleading language, stale code references, missing terminology, and content duplication with the Solverz Cookbook power-flow chapter.

Documentation

  • New Glossary section at the top of the matrix calculus chapter, defining SpMV, SpMM, CSC / CSR, @njit, scatter-add, fancy indexing, fast path, fallback path, lambdify, hot F / hot J, cold compile, LICM, inline mode, and module printer mode. Every body usage of these terms now hyperlinks to the glossary entry.

  • Supported Operations table clarified. The third column was renamed from “Derivative” to “Jacobian block (∂/∂x for vector x)” and a paragraph was added above the table explaining that Solverz first computes the elementwise vector derivative (e.g. cos(x) for sin(x)) and only inserts the diag(...) matrix wrapper at Jacobian assembly time. The previous wording risked making vector-equation users think sin(x) produces a matrix directly.

  • Stale Mat_Mul + scipy.sparse claim removed. The first note block under ## Supported Operations previously said “Mat_Mul uses scipy.sparse directly”, which was true for 0.8.0 but false for 0.8.1 (where the fast path moves the matvec into SolCF.csc_matvec inside inner_F). The note now describes the fast / fallback split correctly and points users at the Layer 1 discussion.

  • Newton-step language replaced with solver-neutral phrasing. matrix_calculus.md is the API reference for Solverz’s matrix calculus engine, which is shared across AE / FDAE / DAE / ODE. Wherever the previous text said “every Newton step assembles Jacobian block data”, it now says “every J_(y, p) call” or “every solver step”, with one explicit paragraph noting that the cost model applies to all J_ consumers, not just algebraic-equation Newton solvers. The remaining “Newton” references are intentional and concrete (e.g. naming Newton-Raphson as one specific solver among several).

  • Layer 1 narrative refreshed for the second-pass review fixes. References to the obsolete _is_csc_matvec_fast_path helper were updated to the current _classify_matmul_placeholders (with inner shape predicate _shape_is_fast). A new paragraph describes the dependency-aware demotion introduced in the second-pass review fix: a fast-path candidate consumed by a fallback placeholder’s matrix or operand expression is demoted to fallback so the wrapper never emits a reference to a not-yet-materialised placeholder.

  • Fallback path subsections rewritten. The earlier “Why not fancy indexing?” subsection title and prose suggested that fancy indexing was an option Solverz “still supports”. The new title is “Two paths: scatter-add (fast) and fancy indexing (fallback)”, and the body makes explicit that both paths always co-exist in every generated module — the runtime picks per Jacobian block based on what the symbolic classifier recognised at code-gen time. There is no on/off switch.

  • Benchmark environment block added to the Performance subsection: hardware (Apple M4), OS (macOS 26.4), Python (3.11.13), library versions (numpy 2.3.3, scipy 1.16.0, numba 0.65.0, sympy 1.13.3, Solverz 0.8.1+), and methodology (10 warm-up + 5000–20000 timed iterations, median of three repeats). Without this metadata the absolute numbers in the doc were unreproducible.

  • “When to use Mat_Mul” section slimmed down, deferring the full case30 decision matrix and the known hot-F regression case to the Solverz Cookbook power-flow chapter. The Solverz-dev doc is the API reference and now keeps just the 3-bullet API guidance plus the “matrix shapes that fall out of the fast path” reference list. The Cookbook is the right place for the case-driven narrative; this avoids contradictions when one of the two docs drifts.

The companion Solverz-Cookbook v0.8.2 release does the same three things in pf.md (terminology cleanup, benchmark environment block, refined “Newton step” wording on the one line where it overgeneralised). The Cookbook’s heat flow chapter is unchanged.

0.8.1

Hotfix release addressing correctness findings in the 0.8.0 Mat_Mul mutable-matrix Jacobian code generation. Users of Mat_Mul should upgrade — 0.8.0 has a latent miscompilation when two independent Diag(...) terms land on the same (i,i) positions (see Bug Fixes). The release also brings a large runtime-performance improvement to the Mat_Mul hot-F path (see Performance below).

Performance

  • Mat_Mul hot F: 4.4× speedup on small networks. Before 0.8.1 every Mat_Mul(A, v) precompute executed as a scipy.sparse SpMV in the Python F_ wrapper (_sz_mm_N = A @ v). On small systems this was dispatch-bound: each SpMV crossed Python → scipy → C → back and cost ≈ 1.5 µs per call of pure overhead, regardless of how little arithmetic the matrix actually had. Eight SpMVs per power-flow F_() call added up to ≈ 12 µs of scipy dispatch on case30, compared to well under 1 µs of actual matvec work.

    The 0.8.1 code generator now recognises plain sparse dim=2 Para matrix operands and emits the matvec inside inner_F using the existing SolCF.csc_matvec Numba helper (Solverz/num_api/custom_function.py). The CSC decomposition (<name>_data / <name>_indices / <name>_indptr / <name>_shape0) that Solverz/model/basic.py already emits for every sparse dim=2 Param is the direct input for the helper — no new setting entries, no new helper function.

    Case30 power-flow benchmark (per F_() call, Apple M4):

    Formulation

    0.8.1 baseline

    0.8.1 fast path

    Speedup

    Mat_Mul (rectangular)

    14.1 µs

    3.23 µs

    4.4×

    For-loop (polar) reference

    1.40 µs

    1.11 µs

    The Mat_Mul / polar hot-F ratio drops from 10.1× to 2.9×. J call, cold compile, module render, and every other phase are unchanged. The remaining 2.9× gap is the structural cost of 8 SolCF.csc_matvec calls + 3 sub-function calls + dispatcher in Mat_Mul vs 53 inlined scalar kernels in polar; closing it further would require SpMV fusion or switching to CSR format, neither of which is in this release. See When to use (and not use) Mat_Mul in the matrix calculus guide for a full decision matrix.

  • Fallback pathMat_Mul placeholders whose matrix operand is not a plain sparse Para (negated matrices Mat_Mul(-A, x), nested matrix expressions, dense dim=2 params) keep the old scipy SpMV path in the wrapper. They are functionally correct but do not benefit from the fast path; users who hit performance regressions should check whether they can rewrite -Mat_Mul(A,x) instead of Mat_Mul(-A,x), or declare matrices with sparse=True.

  • print_F dead-load cleanup. With the fast path in place, the F_ wrapper no longer loads A = p_["A"] for sparse dim=2 matrices that are used only as fast-path Mat_Mul operands — previously that line was dead but still emitted. The filter inspects each placeholder’s matrix_arg and only keeps the matrix load if at least one fallback Mat_Mul needs it.

Bug Fixes

  • Multiple Diag terms now accumulate correctly. Equations whose Jacobians contain two or more independent Diag terms sharing output positions — e.g. x*(A@x) + x*(B@x) producing diag(A@x) + diag(B@x) + diag(x)@A + diag(x)@B — previously had one of the two diagonal contributions silently overwritten in the module-printer path. The scatter-add kernel now uses += for every term, including diag. Inline mode was already correct. Regression test: test_multi_diag_accumulation.

  • Model no longer crashes on dense dim=2 parameters. model.create_instance() unconditionally tried to decompose every dim=2 parameter into sparse CSC flat arrays (.data, .indices, .indptr), which for a dense ndarray fed a memoryview into Array and raised TypeError: Unsupported array type <class 'memoryview'>. Decomposition is now restricted to sparse and dim == 2; dense matrices pass through untouched.

  • Selective @njit gating respects sparse parameter content. The inner_F / inner_J helpers are now generated without @njit when the generated parameter list contains a sparse dim=2 param, a triggerable param, or a TimeSeriesParam — objects Numba cannot lower. Pure element-wise models are unaffected.

  • FormJac and JacBlock agree on the mutable-matrix predicate. The “is this block a mutable-matrix block?” decision is now made in a single place using the same criterion on both sides — matrix-valued derivative that is not a plain Para / -Para. Previously FormJac additionally required the expression to contain both Mat_Mul and Diag, while JacBlock.is_mutable_matrix did not. The divergence would have let a derivative like Diag(x) skip the SpDiag perturbation step and produce a shrunken CooRow / CooCol at a flat start — the downstream scatter-add kernel would then write to fewer output positions than the runtime expected. No known model hit this in practice but the fix closes the corner case.

API Changes

  • Time-varying sparse dim=2 Params are now rejected at construction. Declaring a Param(..., dim=2, sparse=True, triggerable=True, ...) or a TimeSeriesParam(..., dim=2, sparse=True) raises NotImplementedError at the point of construction, regardless of whether the parameter is ever referenced in an equation.

    Every Solverz code path that consumes a sparse dim=2 Param caches its CSC decomposition (<name>_data, <name>_indices, <name>_indptr, <name>_shape0) at model-build time: the legacy MatVecMul pipeline, the new 0.8.1 Mat_Mul SolCF.csc_matvec fast path, and the mutable-matrix Jacobian scatter-add kernel all read the frozen flats. A runtime trigger_fun firing or a get_v_t(t) update simply gets lost, and the Newton iteration either diverges or silently converges to the wrong solution.

    Earlier 0.8.1 drafts narrowed this check to “sparse dim=2 time-varying Para used as the matrix operand of a Mat_Mul”, catching the combination only at FormJac time. That was a loophole: legacy MatVecMul usage, or any other code path that consumes the CSC flats, slipped through. The 0.8.1 policy closes the loophole by rejecting the shape at the exact line where it is declared — the error message now points at the user’s source, not a deep internal.

    The check lives in ParamBase.__init__ and TimeSeriesParam.__init__, with a backstop in Equations._check_no_timevar_sparse_matrices (runs on every FormJac call) for the edge case where a tainted Param was built via __new__ + attribute assignment, bypassing the __init__ guard.

    Allowed alternatives: triggerable / time-series vectors or scalars sitting next to a Mat_Mul, element-wise formulations where per-row coefficients are 1-D time-varying parameters, and dense dim=2 parameters (sparse=False, which takes the MutableMatJacDataModule fallback path that re-evaluates the full block expression on every J_() call and tolerates runtime updates). See the Restrictions section of the Matrix Calculus guide for full workarounds.

  • Reserved symbol prefixes _sz_mm_ and _sz_mb_. Any user symbol (Var, Param, iVar, Para, …) whose name starts with either prefix is rejected at construction time with a ValueError. These prefixes are used by the code generator for Mat_Mul precompute helpers and mutable-matrix Jacobian block helpers. The check is bypassed internally via internal_use=True.

  • Dense dim=2 params in Mat_Mul emit a one-shot UserWarning. A parameter declared with dim=2, sparse=False and used inside a Mat_Mul works correctly via the fallback path but forfeits the scatter-add fast path. FormJac now warns once per offending parameter to flag the performance cost. See the new Restrictions and reserved names section of the matrix calculus guide for migration guidance.

Documentation

  • Extended Matrix-Vector Calculus:

    • New Restrictions and reserved names section documenting the three API boundaries introduced in 0.8.1.

    • Code examples updated to show the _sz_mm_ / _sz_mb_ helper names actually emitted by the code printer.

    • Explicit note about the += accumulation rule for diagonal scatter-add terms.

    • New explicit list of the cases that push a mutable-matrix block onto the fallback path — useful when a model is slower than expected and you want to know whether a block is hitting the fast or slow path.

    • The immutability warning now spells out that the restriction only strictly applies to matrices used in the vectorised Mat_Mul fast path (the fallback path re-evaluates the full expression each call and would reflect mutations), but still recommends treating all sparse matrix params as immutable.

Internal cleanup

  • Removed the dead include_sparse_in_list parameter from print_F / print_inner_F / print_J / print_inner_J and from _has_sparse_in_param_list. With the Mat_Mul precompute architecture every caller hard-coded False and the parameter was vestigial.

  • Dropped the unused _var_base_name / _var_access_expr helpers from mutable_mat_analyzer; they were leftovers from an earlier draft of print_inner_J.

Full Changelog

0.8.0…0.8.1

0.8.0

Highlights

Complete Matrix-Vector Calculus — Solverz now fully supports symbolic differentiation of mixed matrix-vector equations. Write equations like e*(G@e - B@f) + f*(B@e + G@f) - P and get analytical Jacobians automatically.

Unified Mat_Mul InterfaceMat_Mul(A, x) replaces the legacy MatVecMul as the standard matrix-vector product. It uses scipy.sparse directly (faster than the old csc_matvec) and supports full matrix calculus.

New Features

  • Matrix calculus operators: exp, sin, cos, ln, power (**), transpose, and Diag now work inside matrix-vector expressions with automatic differentiation.

    from Solverz import Mat_Mul, Var, Param, Model, Eqn
    from Solverz.sym_algebra.functions import exp
    
    m = Model()
    m.A = Param('A', [[1, 0], [0, 2]], dim=2, sparse=True)
    m.x = Var('x', [1, 1])
    m.f = Eqn('f', exp(Mat_Mul(m.A, m.x)))  # Jacobian: diag(exp(A@x)) @ A
    
  • Mutable matrix Jacobian: Variable-dependent matrix derivatives (e.g., diag(e)@G + diag(f)@B from power flow equations) are now evaluated dynamically at each Newton step. The sparsity pattern is determined at initialization and remains fixed; only the data values are updated.

  • Selective Numba @njit: In module mode, Numba compilation is applied selectively — equations using Mat_Mul run with scipy.sparse (fast C-level sparse operations), while pure element-wise equations retain @njit acceleration.

  • atan2 symbolic function: Added atan2(y, x) for computing the two-argument arctangent in symbolic equations.

  • Plugin-based module discovery: Third-party numerical modules (e.g., SolMuseum) are now discovered via entry_points(group='solverz.num_api') instead of hard-coded imports. Packages register via pyproject.toml [project.entry-points."solverz.num_api"]. Closes #118.

  • Improved solution stats: Solvers now record more detailed statistics and profiling information in the solution object.

Bug Fixes

  • Stabilized solution slicing and incidence matrix helpers.

Deprecations

  • MatVecMul is deprecated — use Mat_Mul(A, x) instead. MatVecMul will emit a DeprecationWarning when used. It will be removed in a future release.

    # Before (deprecated):
    from Solverz import MatVecMul
    m.f = Eqn('f', MatVecMul(m.A, m.x) - m.b)
    
    # After (recommended):
    from Solverz import Mat_Mul
    m.f = Eqn('f', Mat_Mul(m.A, m.x) - m.b)
    

Documentation

  • New: Matrix-Vector Calculus — functionality, mathematical background, and application examples (power flow, heat network, nonlinear equations).

  • New: Extending Matrix Calculus — developer guide for adding new operations to the matrix calculus module.

  • Updated: Getting Started — matrix equation examples now use Mat_Mul.

Full Changelog

0.7.2…0.8.0