Skip to content

Fix \ for Dual #236

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
merged 18 commits into from
Apr 11, 2025
34 changes: 34 additions & 0 deletions src/overloads/arrays.jl
Original file line number Diff line number Diff line change
Expand Up @@ -122,6 +122,19 @@ function LinearAlgebra.:\(
return Ainv * B
end

function LinearAlgebra.:\(
A::AbstractMatrix, B::AbstractVecOrMat{T}
) where {T<:AbstractTracer}
return inv(A) * B
end

function LinearAlgebra.:\(
A::AbstractMatrix{T1}, B::AbstractVecOrMat{T2}
) where {T1<:AbstractTracer,T2<:AbstractTracer}
Ainv = LinearAlgebra.pinv(A)
return Ainv * B
end

## Exponential
function Base.exp(A::AbstractMatrix{T}) where {T<:AbstractTracer}
LinearAlgebra.checksquare(A)
Expand Down Expand Up @@ -187,6 +200,27 @@ function LinearAlgebra.logabsdet(A::AbstractMatrix{D}) where {D<:Dual}
t1, t2 = LinearAlgebra.logabsdet(tracers)
return (D(p1, t1), D(p2, t2))
end
function LinearAlgebra.:\(A::AbstractMatrix{<:Dual}, B::AbstractVector)
primals, tracers = split_dual_array(A)
p = primals \ B
t = tracers \ B
return Dual.(p, t)
end
function LinearAlgebra.:\(A::AbstractMatrix, B::AbstractVector{D}) where {D<:Dual}
primals, tracers = split_dual_array(B)
p = A \ primals
t = A \ tracers
return Dual.(p, t)
end
function LinearAlgebra.:\(
A::AbstractMatrix{D1}, B::AbstractVector{D2}
) where {D1<:Dual,D2<:Dual}
A_primals, A_tracers = split_dual_array(A)
B_primals, B_tracers = split_dual_array(B)
p = A_primals \ B_primals
t = A_tracers \ B_tracers
return Dual.(p, t)
end

#==============#
# SparseArrays #
Expand Down
3 changes: 3 additions & 0 deletions test/test_arrays.jl
Original file line number Diff line number Diff line change
Expand Up @@ -399,6 +399,9 @@ end

x = rand(2)
b = A \ x

A = rand(2, 2)
x = [t1; t2]
@test all(t -> sameidx(t, s_out), b)
end

Expand Down
22 changes: 22 additions & 0 deletions test/test_constructors.jl
Original file line number Diff line number Diff line change
Expand Up @@ -225,3 +225,25 @@ end
end
end
end

# https://github.com/adrhill/SparseConnectivityTracer.jl/issues/235
@testset "Left division on Dual with BigFloat" begin
@testset "$T" for T in union(GRADIENT_TRACERS, HESSIAN_TRACERS)
P = IndexSetGradientPattern{Int,BitSet}
T = GradientTracer{P}

p = P(BitSet(2))
t_full = T(p)
A = rand(BigFloat, 2, 2)
B = rand(BigFloat, 2)

dualized_A = Dual.(A, t_full)
dualized_B = Dual.(B, t_full)
dualized_x = dualized_A \ dualized_B

@test_nowarn dualized_x = dualized_A \ dualized_B
@test_nowarn A \ dualized_B
@test_nowarn dualized_A \ B
@test eltype(dualized_x) == eltype(dualized_A)
end
end
Loading