Fall back to ExpMethodGeneric for immutable matrices in exponential!#237
Merged
ChrisRackauckas merged 1 commit intoJul 3, 2026
Conversation
exponential!(A) with no method previously always selected the in-place ExpMethodHigham2005, which throws a setindex! error on immutable inputs such as StaticArrays' SMatrix. ExpMethodGeneric already has a fully non-mutating code path for immutable matrices, so route them there. This also makes ForwardDiff.jacobian(exponential!, ::SMatrix) work, including for Float32 where StaticArrays' own exp overflows its Pade coefficients through the 3x3 cofactor solve and returns NaN; the normalized Pade coefficients of ExpMethodGeneric avoid that overflow. See https://discourse.julialang.org/t/how-to-use-forwarddiff-jl-with-inv-and-matrix-exponential/137880 Co-Authored-By: Chris Rackauckas <accounts@chrisrackauckas.com>
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
2 participants
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.
Note
This PR should be ignored until reviewed by @ChrisRackauckas.
What
exponential!(A)with no method argument now falls back toExpMethodGeneric()whenAis immutable (perArrayInterface.ismutable, the same predicateExpMethodGenericuses internally), instead of always selecting the in-placeExpMethodHigham2005.Why
exponential!(A)on a StaticArraysSMatrixcurrently throws:because the default
ExpMethodHigham2005path mutates. MeanwhileExpMethodGenericalready contains a fully non-mutating code path for immutable inputs (theismutable(x)guard insrc/exp_generic.jl), so the capability exists — it just isn't reachable from the no-method entry point.This came out of https://discourse.julialang.org/t/how-to-use-forwarddiff-jl-with-inv-and-matrix-exponential/137880, where the recommendation to use the (now-removed name)
exp_genericon static arrays failed with the error above.Two additional data points verified locally while preparing this:
ForwardDiff.jacobian(exponential!, ::SMatrix)works through the new fallback and matches finite differences (~1e-7) for 3×3/4×4, Float64 and Float32.SMatrix{3,3,Float32}, StaticArrays' ownexpgives NaN Jacobians in ~1/3 of random samples (its Base-style Padé coefficients, up to ~3.2e16, overflow Float32 through the 3×3 cofactor solve — as diagnosed by mikmoore on Discourse).ExpMethodGenericgave 0 NaNs in 1000 samples because its Padé coefficients are normalized (≤ 1). So the fallback is also the numerically safer path for Float32 static arrays, not just the one that doesn't error.Changes
src/exp.jl: default method selection routes immutable matrices toExpMethodGeneric(); docstring note added.test/basictests.jl: new testset coveringexponential!(::SMatrix)values vsexp(Matrix(A))andForwardDiff.jacobianvs theMatrixreference, for 3×3/4×4 × Float64/Float32, including NaN checks.Behavior note: scalars also now work through the no-method entry point (
ismutable(::Number) == false→ExpMethodGeneric), where previouslyexponential!(1.0)threw aMethodError. Purely additive.Testing
GROUP=Core Pkg.test()run locally on Julia 1.11.9:Core/basictests.jl | 349 passed, 1 broken (pre-existing kiops @test_broken), 0 failures. Runic applied to the changed files.🤖 Generated with Claude Code
https://claude.ai/code/session_01WQWjLDJ228E8w9631ZeTdM