COMP: Prefer the tiered itk::bridge::Math spelling for ITK's Eigen numerics - #2019
Open
hjmjohnson wants to merge 6 commits into
Open
COMP: Prefer the tiered itk::bridge::Math spelling for ITK's Eigen numerics#2019hjmjohnson wants to merge 6 commits into
hjmjohnson wants to merge 6 commits into
Conversation
hjmjohnson
force-pushed
the
enh-tier-unsupported-math-6620
branch
from
August 20, 2026 21:48
1419483 to
853c0d4
Compare
ITK moved the Eigen-backed convenience numerics into itk::bridge and renamed the headers to itkBridge*.h. Probe the tiered header first so a post-move ITK stays on the Eigen path, and keep the pre-move itkMathSVD.h arm so older ITK checkouts still resolve. The itk::bridge tier is explicitly a migration aid rather than supported ITK API, so the ANTS_ITK_MATH_NS indirection is what limits the exposure: the spelling is named once here, not at each call site.
Route the symmetric 3x3 diffusion-tensor inverse in GetMetricTensorCost
(both overloads) through a gated matHelper::SymmetricInverse helper. When a
new-enough ITK provides itk::Math::SolveSymmetric (capability macro
ITK_MATH_HAS_SOLVE_SYMMETRIC), DT^-1 is built column-by-column from Eigen
LDLT solves of DT X = I; otherwise the legacy vnl_matrix_inverse is used.
Guarded by __has_include plus the macro, so it builds unchanged against the
currently pinned ITK (which takes the legacy path).
Math justification:
The diffusion tensor DT is symmetric positive-(semi)definite, so DT^-1 is
symmetric and the metric cost v^T (DT^-1)^p v is well posed. LDLT
(Bunch-Kaufman symmetric pivoting) is the natural factorization for a
symmetric system; solving DT X = I column by column yields DT^-1 exactly and
handles near-singular tensors more gracefully than a general inverse. (For
p = 1 the full inverse is unnecessary, since v^T DT^-1 v = v . solve(DT, v),
and for even p the form reduces to ||(DT^-1)^(p/2) v||^2 via repeated solves;
the drop-in keeps the existing matrix-power squaring and so still materializes
the inverse, for a minimal, behavior-preserving diff.)
Validation:
Cross-checked against vnl_matrix_inverse on random SPD 3x3 tensors over
matrixpower in {1,2,3} (18 cases): worst relative difference 2.6e-15 for the
inverse-via-solve form used here, and 1.7e-15 for the quadratic-form shortcut
-- machine precision. The gated header compiles against the pinned ITK via the
legacy path (tensor consumers rebuilt clean).
Part of the vnl -> itk::Math de-vendoring effort; this symmetric-inverse hot
path was identified in the joint-label-fusion solver analysis.
… (gated) Route the symmetric left-Cauchy-Green tensor inverse in VectorFieldGradientImageFunction::EvaluateEulerianStrainTensor (point and index overloads) through a gated vfg_detail::SymmetricInverse helper. When a new-enough ITK provides itk::Math::SolveSymmetric (capability macro ITK_MATH_HAS_SOLVE_SYMMETRIC), the inverse is built column-by-column from Eigen LDLT solves of B X = I; otherwise the legacy vnl_matrix_inverse is used. Guarded by __has_include plus the macro, so it builds unchanged against the currently pinned ITK. Companion to the diffusion-tensor swap in this branch. Math justification: The Eulerian (Almansi) strain is E = 1/2 (I - (F F^T)^-1), where F is the deformation gradient and B = F F^T is the left Cauchy-Green tensor. B is symmetric positive-definite by construction, so its inverse is symmetric and LDLT (Bunch-Kaufman symmetric pivoting) is the appropriate factorization; solving B X = I column by column yields B^-1 exactly. The full inverse matrix is required here (E subtracts it from the identity element-wise), so this is the form-the-inverse-via-solve case rather than a single solve. Validation: Cross-checked the full strain tensor E against the vnl_matrix_inverse path on random near-identity deformation gradients at the default float precision (8 cases): worst elementwise |dE| = 2.4e-7, i.e. float machine precision -- and E = 1/2(I - B^-1) subtracts the inverse from I with no averaging, so the error is undamped. The gated header compiles against the pinned ITK via the legacy path (ImageMathHelper consumers rebuilt clean). Part of the vnl -> itk::Math de-vendoring effort; this symmetric-inverse hot path (per-point strain evaluation) was identified in the forest audit.
…mmetric Replace the open-coded column-by-column inverse loops in the gated SymmetricInverse helpers (TensorFunctions.h, itkVectorFieldGradientImageFunction.hxx) with a single itk::Math::InverseSymmetric(A) call. The prior loop re-factorized A once per column (O(n^4)); InverseSymmetric factorizes once (O(n^3)). Behavior and the #ifdef ITK_MATH_HAS_SOLVE_SYMMETRIC gate are unchanged; the legacy vnl_matrix_inverse fallback is retained and still compiles against the pinned ITK. Requires an ITK carrying itk::Math::InverseSymmetric on the new path.
…ted)
Route two symmetric-matrix inverses in antsSCCANObject through a gated
sccan_detail::SymmetricInverse helper (itk::Math::InverseSymmetric, Eigen LDLT,
when available; vnl_matrix_inverse otherwise):
- VNLPseudoInverse: regularization*I + b*b^T (SPD Gram), previously inverted
via vnl_svd(mat).inverse().
- The PCG preconditioner chollow*diaginv*chollow^T (symmetric), previously
vnl_matrix_inverse.
Both operands are symmetric, so LDLT is the appropriate factorization. Gated by
__has_include plus ITK_MATH_HAS_SOLVE_SYMMETRIC; the sccan tool rebuilds clean
against the pinned ITK via the legacy path. Sites from the forest audit.
ITK moved the Eigen-backed convenience numerics into itk::bridge and renamed the headers to itkBridge*.h. Probe the tiered header first so a post-move ITK stays on the Eigen path, and keep the pre-move arm so older ITK checkouts still resolve; the vnl fallback is unchanged. The itk::bridge tier is explicitly a migration aid rather than supported ITK API. These call sites are guarded by ITK_BRIDGE_MATH_HAS_SOLVE_SYMMETRIC and fall back to vnl, so ANTs does not hard-depend on that tier.
hjmjohnson
force-pushed
the
enh-tier-unsupported-math-6620
branch
from
August 21, 2026 22:52
853c0d4 to
0202eb6
Compare
hjmjohnson
marked this pull request as ready for review
August 21, 2026 23:39
hjmjohnson
added a commit
to hjmjohnson/itk_forest_build_testbed
that referenced
this pull request
Aug 22, 2026
ITK variants re-minted off current upstream tips, patch sets re-derived rather
than copied:
itk-main slicer-v6.0.0-2026-08-22-fee22c1e224 (main + namespace patch)
itk-release-5.4 slicer-v5.4.7-2026-08-21-7fefad582bf (5.4.7 + namespace,
itkFileTools)
main needs one patch and release-5.4 needs two; the itkFileTools fix has landed
upstream in main but not in 5.4. Verified per docs/slicer-itk-policy.md rather
than carrying the previous variant's count forward.
Both ANTs sub-pins -- the ANTs that BRAINSTools and SlicerANTs build internally
-- move from a bare SHA to hjmjohnson/ANTs@integration/itk-forest-2026-08-22,
which is identical to ANTsX/ANTs#2019 head. Tracking the branch keeps the inner
ANTs in step with the PR as it evolves; a bare SHA silently drifts, which is how
an earlier sweep ended up building a 2024-era ANTs.
#2019 is the third piece of the itk::bridge::Math alignment: ITK #6768 provides
the namespace, BRAINSTools f34f9f75 already consumes it, and without #2019 ANTs
is the only consumer left on the old spelling.
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
Add this suggestion to a batch that can be applied as a single commit.This suggestion is invalid because no changes were made to the code.Suggestions cannot be applied while the pull request is closed.Suggestions cannot be applied while viewing a subset of changes.Only one suggestion per line can be applied in a batch.Add this suggestion to a batch that can be applied as a single commit.Applying suggestions on deleted lines is not supported.You must change the existing code in this line in order to create a valid suggestion.Outdated suggestions cannot be applied.This suggestion has been applied or marked resolved.Suggestions cannot be applied from pending reviews.Suggestions cannot be applied on multi-line comments.Suggestions cannot be applied while the pull request is queued to merge.Suggestion cannot be applied right now. Please check back later.
Keeps ANTs on ITK's Eigen-backed numerics now that InsightSoftwareConsortium/ITK#6768
has merged, moving the convenience numerics to
itk::bridge::Math(InsightSoftwareConsortium/ITK#6620). Chained
__has_includeguards keep both thepost-move and pre-move ITK on the Eigen path, so this builds against either vintage.
Affects
Utilities/antsNearestRotation.h(nearest-rotation / Wahba solve),Tensor/TensorFunctions.h,Utilities/antsSCCANObject.h, andUtilities/itkVectorFieldGradientImageFunction.hxx.Why this is needed even though nothing breaks
The pre-existing guards were bare
__has_includeprobes of the pre-move header:After the ITK rename that probe simply fails. There is no compile error — the
capability macro stays undefined and the code silently falls back to
vnl_matrix_inverse/vnl_svd, with nothing in any log to say so.Chaining the probe keeps either ITK vintage on the Eigen path:
Because both arms are retained, this is safe against pre-move and post-move ITK alike.
Test plan — built and measured against merged ITK
Validated against ITK
mainat313b5234d52(i.e. after #6768 merged), installedfrom a dedicated build with the remote modules ANTs requires
(
GenericLabelInterpolator,AdaptiveDenoising,ITKReview,ITKMINC,ITKIOTransformMINC,ITKTBB).1. Guard-arm probe. A syntax-only compile of each guard against the installed
ITK, with one
#errorper arm, reports which branch the preprocessor actually takes:itkBridgeMathLDLT.hchainTensorFunctions.h,antsSCCANObject.h,itkVectorFieldGradientImageFunction.hxxARM=BRIDGE(eigen)itkBridgeMathSVD.hchainantsNearestRotation.hARM=BRIDGE(eigen)Neither falls through to the
vnlarm. This matters because a successful build isconsistent with all three arms — the
vnlfallback compiles fine, it just quietlyabandons Eigen. The probe is what distinguishes them.
2. Full build.
ninjaover the whole project: 417/417 targets, 0 errors,0
FAILEDlines.Not run: the ANTs test suite. This change alters which numerical backend is selected,
not the algorithms, and the Eigen/vnl equivalence was measured previously
(<=2.6e-15 on random SPD 3x3 tensors).