Unwrap AutoSpecializeCallable for Enzyme AD and fix Runic formatting

Enzyme cannot differentiate through FunctionWrappers' llvmcall, causing
EnzymeMutabilityException in all IIP Vector{Float64} tests with AutoEnzyme.
Unwrap the function in construct_jacobian_cache when the AD backend is
Enzyme-based (including AutoSparse(AutoEnzyme(...))), so DI sees the raw
user function. Also apply Runic formatting to SCCNonlinearSolve files.

Co-Authored-By: Chris Rackauckas <accounts@chrisrackauckas.com>
Co-Authored-By: Claude Opus 4.6 <noreply@anthropic.com>

Author: ChrisRackauckas

lib/NonlinearSolveBase/src/autospecialize.jl

@@ -81,6 +81,17 @@ Otherwise return `f` unchanged.
 get_raw_f(f) = f
 get_raw_f(f::AutoSpecializeCallable) = f.orig
 
+"""
+    _uses_enzyme_ad(ad) -> Bool
+
+Return `true` if `ad` is an Enzyme-based AD backend (possibly wrapped in `AutoSparse`).
+Enzyme cannot differentiate through FunctionWrappers' `llvmcall`, so
+`AutoSpecializeCallable` must be unwrapped before passing to DI with Enzyme.
+"""
+_uses_enzyme_ad(::ADTypes.AutoEnzyme) = true
+_uses_enzyme_ad(ad::AutoSparse) = _uses_enzyme_ad(ADTypes.dense_ad(ad))
+_uses_enzyme_ad(_) = false
+
 # Default dispatch assumes no ForwardDiff loaded.
 # The ForwardDiff extension overrides these with dual-aware versions.
 

lib/NonlinearSolveBase/src/jacobian.jl

@@ -53,6 +53,11 @@ function construct_jacobian_cache(
         end
         autodiff = standardize_forwarddiff_tag(autodiff, prob)
         autodiff = construct_concrete_adtype(f, autodiff)
+        # Enzyme cannot differentiate through FunctionWrappers' llvmcall.
+        # Unwrap AutoSpecializeCallable so DI sees the raw user function.
+        if is_fw_wrapped(f.f) && _uses_enzyme_ad(autodiff)
+            f = @set f.f = get_raw_f(f.f)
+        end
         di_extras = if SciMLBase.isinplace(f)
             DI.prepare_jacobian(f, fu_cache, autodiff, u, Constant(p), strict = Val(false))
         else

lib/SCCNonlinearSolve/src/SCCNonlinearSolve.jl

@@ -85,9 +85,11 @@ function solve_single_scc(alg, prob, explicitfun, sols; kwargs...)
         )
     else
         sol = SciMLBase.solve(prob, alg.nlalg; kwargs...)
-        SciMLBase.strip_solution(SciMLBase.build_solution(
-            prob, nothing, sol.u, sol.resid, retcode = sol.retcode
-        ))
+        SciMLBase.strip_solution(
+            SciMLBase.build_solution(
+                prob, nothing, sol.u, sol.resid, retcode = sol.retcode
+            )
+        )
     end
 
     return _sol
@@ -101,20 +103,24 @@ function iteratively_build_sols(alg, probs::AbstractVector, explicitfuns::Abstra
     uType = typeof(probvec(prob1))
     T = eltype(uType)
     rType = uType  # resid has same type as u for nonlinear problems
-    ST = SciMLBase.NonlinearSolution{T, 1, uType, rType,
-        NamedTuple{(:p,), Tuple{Nothing}}, Nothing, Nothing, Nothing, Nothing, Nothing}
+    ST = SciMLBase.NonlinearSolution{
+        T, 1, uType, rType,
+        NamedTuple{(:p,), Tuple{Nothing}}, Nothing, Nothing, Nothing, Nothing, Nothing,
+    }
     sols = Vector{ST}(undef, length(probs))
     for i in eachindex(probs)
-        sols[i] = solve_single_scc(alg, probs[i], explicitfuns[i], view(sols, 1:i-1); kwargs...)
+        sols[i] = solve_single_scc(alg, probs[i], explicitfuns[i], view(sols, 1:(i - 1)); kwargs...)
     end
     return sols
 end
 
 @generated function iteratively_build_sols(alg, probs::Tuple, explicitfuns::Tuple, ::Val{N}; kwargs...) where {N}
     return quote
         Base.Cartesian.@nexprs $N i -> begin
-            prob_i = solve_single_scc(alg, probs[i], explicitfuns[i],
-                Base.Cartesian.@ntuple((i - 1), j -> prob_j); kwargs...)
+            prob_i = solve_single_scc(
+                alg, probs[i], explicitfuns[i],
+                Base.Cartesian.@ntuple((i - 1), j -> prob_j); kwargs...
+            )
         end
         return Base.Cartesian.@ntuple $N i -> prob_i
     end

lib/SCCNonlinearSolve/test/core_tests.jl

@@ -248,8 +248,10 @@ end
     # Wrap explicitfuns with FunctionWrapper for type unification.
     # The stripped solution type is deterministic — compute it from u0 type.
     uType = Vector{Float64}
-    SSol = SciMLBase.NonlinearSolution{Float64, 1, uType, uType,
-        NamedTuple{(:p,), Tuple{Nothing}}, Nothing, Nothing, Nothing, Nothing, Nothing}
+    SSol = SciMLBase.NonlinearSolution{
+        Float64, 1, uType, uType,
+        NamedTuple{(:p,), Tuple{Nothing}}, Nothing, Nothing, Nothing, Nothing, Nothing,
+    }
     SolsView = SubArray{SSol, 1, Vector{SSol}, Tuple{UnitRange{Int64}}, true}
     EFW = FunctionWrapper{Nothing, Tuple{Vector{Float64}, SolsView}}
     ef1_wrapped = EFW(explicitfun1_raw)