Jacobian Verification

When implementing analytical Jacobians (OdeJacobian(), OutputPortJacobian(), InputPortJacobian()), it is easy to introduce sign errors or incorrect entries. FhSim provides a built-in verification mode that compares the analytical/hybrid Jacobian against a purely numerical Jacobian computed via forward finite differences.

<SIMULATION>
    <Integrator Method="BDF">
        <Diagnostics>
            <Settings>
                <JacobianVerify/>
            </Settings>
        </Diagnostics>
    </Integrator>
</SIMULATION>

Jacobian verification is controlled by <Diagnostics><Triggers>:

Configuration	Description
no `<JacobianVerify/>`	No verification (default)
`<JacobianVerify/>`	Verification enabled
`<Triggers Once="true"/>`	Verify once, then disable
`<Triggers Once="false"/>`	Verify whenever trigger conditions are met

The verification uses a relative error threshold of 1e-4 (typical for finite differences with O(ε^0.875) perturbation). When mismatches are detected, the output includes:

Summary: max absolute error, max relative error, total mismatch count
Per-entry details for up to 100 mismatches: row/column indices, analytical value, numerical value, relative error, and state names (e.g., Body1::pos_x, Body2::vel_y)

Example output:

Jacobian verification: max_abs_err=0.0234, max_rel_err=0.891, mismatches=3 out of 144 entries (tol=0.0001)

J[2,5]: analytical=-1.0, numerical=0.123, rel_err=0.891 [Body1::vel_z, Body2::pos_x]

Sparsity Pattern Visualization

To understand the Jacobian structure and verify port connections are producing the expected sparsity pattern, FhSim can dump the sparsity pattern at model initialisation. Both Engine and Sundials paths support the same Jacobian diagnostics and linear-solver options via <Integrator Method="..."> and related sub-elements.

<SIMULATION>
    <Integrator Method="BDF">
        <Diagnostics>
            <Settings>
                <JacobianSparsityDump Mode="full"/>
            </Settings>
        </Diagnostics>
    </Integrator>
</SIMULATION>

JacobianSparsityDump Mode values

Mode	Description
`off`	No dump (default)
`block`	SimObject-level block summary written to `sparsity_block.txt`
`full`	ASCII grid to log (n ≤ 80 only) + `jacobian_sparsity.mtx` + `jacobian_sparsity.png`

Image-size tuning for sparsity PNG output is currently internal (no dedicated XML attribute in the unified diagnostics schema).

Output files ("full" mode)

**jacobian_sparsity.mtx** — Matrix Market coordinate format, pattern general type (no values):

%%MatrixMarket matrix coordinate pattern general
% Jacobian sparsity pattern: 50 x 50, nnz = 148
50 50 148
1 1
1 2
2 1
2 2
2 3
...

This file can be loaded by MATLAB (spconvert), SciPy (scipy.io.mmread), or any Matrix Market reader for further analysis.

Note: The .mtx file is only written when all on-diagonal blocks and port links provide analytical sparsity information. If any SimObject returns HasJacobians()=false or any port link lacks HasPortJacobians(), FhSim cannot guarantee the pattern is tight (it would over-estimate nnz with dense blocks). In this case the .mtx output is suppressed and only the PNG image is written (using the dense-fallback pattern described below).

**jacobian_sparsity.png** — Grayscale image of the sparsity pattern:

n = 50 (1:1 mapping)        n = 5000 (downsampled 3x3 → 1)
┌──────────────────┐         ┌──────────────────┐
│█░░░░░░░░░░░░░░░░░│         │▓▓░░░░░░░░░░░░░░░░│
│██░░░░░░░░░░░░░░░░│         │▓▓▓░░░░░░░░░░░░░░░│
│░██░░░░░░░░░░░░░░░│   →     │░▓▓▓░░░░░░░░░░░░░░│
│░░██░░░░░░░░░░░░░░│         │░░▓▓▓░░░░░░░░░░░░░│
│...               │         │...               │
└──────────────────┘         └──────────────────┘
(black = non-zero)            (darker = denser block)

Dense-fallback image: When the analytical sparsity pattern is incomplete (some SimObjects or port links lack Jacobians), FhSim generates a conservative fallback pattern before rendering the PNG:

Non-analytical on-diagonal blocks → filled as solid (fully dense) black squares.
Non-analytical port links → filled as solid black rectangles at the appropriate row/column offsets.
Analytical entries are shown at their true sparsity.

This fallback pattern is a superset of the true sparsity — every structural non-zero is visible, but some true zeros may appear black. The fallback image is always generated (even when .mtx is suppressed), so you can still visualize the overall coupling structure of a mixed analytical/numerical model. The log will contain a message identifying which SimObjects triggered the dense-fallback.

Block mode output

"block" mode writes sparsity_block.txt with a per-SimObject breakdown, including the half-bandwidth of each on-diagonal block and the column span of each port link:

=== Jacobian Block Sparsity (block mode) ===
  Total states: 24
  SimObjects with states:
    Body1: offset=0, size=6, half_bw=5
    Body2: offset=6, size=6, half_bw=5
    Cable: offset=12, size=12, half_bw=2
  Port Jacobian links (off-diagonal blocks):
    Body1 -> Body2 [6:12, 0:6] analytical=yes  cols=6
    Cable -> Body2 [6:12, 12:24] analytical=no  cols=12

half_bw is the maximum |i-j| over all non-zero entries within the on-diagonal block. For a block-tridiagonal single-SimObject this equals bs-1 (block size minus one). Port link cols is the number of columns in the off-diagonal block.

Note: For SimObjects that do not provide analytical Jacobians (HasJacobians()=false), their on-diagonal block is reported as fully dense (half_bw = size-1). Similarly, port links where either endpoint lacks HasPortJacobians()=true are reported as analytical=no with full column span.

Troubleshooting

Issue	Cause	Solution
Simulation uses numerical Jacobian despite implementing `OdeJacobian()`	Some SimObject returns `false` from `HasJacobians()`	Ensure all SimObjects with states override `HasJacobians()` to return `true`
Newton iteration fails to converge	Jacobian is incorrect	Verify Jacobian using finite difference comparison (see above)
Wrong results with implicit method	Sign error in Jacobian	Check that J[i][j] = dF[i]/dX[j], not dX[i]/dF[j]
Slow simulation with sparse Jacobian	System too small for sparse overhead	Return `-1` from `GetJacobianSparsity()` to use dense format
"Analytical Jacobian required but not available" error	`<Jacobian Policy="analytical">` but some SimObject lacks `OdeJacobian()`	Change to `Policy="auto"` or `"hybrid"`, or implement `OdeJacobian()` for all SimObjects
"Self-connection detected" error	A SimObject's output port is connected to its own input port with implicit integration active	Remove the self-connection or switch to explicit integration
"singular Jacobian produced NaN/Inf" error	The Jacobian matrix is singular at the current state	Check for redundant equations, equilibrium points, or numerical overflow in Jacobian computation
Off-diagonal Jacobian blocks are zero despite port connections	SimObjects do not implement `HasPortJacobians()`	Implement `OutputPortJacobian()` and `InputPortJacobian()`, or rely on hybrid numerical assembly
Sparse solver slower than dense	System too small for sparse overhead (n < 30)	Use `<LinearSolver Type="dense">` or verify sparsity pattern is correct
Iterative solver not converging	Preconditioner ineffective or tolerance too tight	Increase `<LinearSolver MaxIter="...">`, relax `<LinearSolver Tol="...">`, or switch to `<LinearSolver Type="sparse">`

Table of Contents