comments

change some input tests
remove redundant comment
include ComplexF16 in tests
fix unchanged test names and docs
improve allocations

Author: sanderdemeyer

src/implementations/exponential.jl

@@ -8,12 +8,20 @@ copy_input(::typeof(exponential), A::Diagonal) = copy(A)
 
 function check_input(::typeof(exponential!), A::AbstractMatrix, expA::AbstractMatrix, alg::AbstractAlgorithm)
     m, n = size(A)
-    m == n || throw(DimensionMismatch("square input matrix expected"))
-    @assert expA isa AbstractMatrix
+    m == n || throw(DimensionMismatch("square input matrix expected. Got ($m,$n)"))
     @check_size(expA, (m, m))
     return @check_scalar(expA, A)
 end
 
+function check_input(::typeof(exponential!), A::AbstractMatrix, expA::AbstractMatrix, ::DiagonalAlgorithm)
+    m, n = size(A)
+    @assert m == n && isdiag(A)
+    @assert expA isa Diagonal
+    @check_size(expA, (m, m))
+    @check_scalar(expA, A)
+    return nothing
+end
+
 # Outputs
 # -------
 function initialize_output(::typeof(exponential!), A::AbstractMatrix, ::AbstractAlgorithm)
@@ -22,6 +30,10 @@ function initialize_output(::typeof(exponential!), A::AbstractMatrix, ::Abstract
     return expA
 end
 
+function initialize_output(::typeof(exponential!), A::Diagonal, ::DiagonalAlgorithm)
+    return similar(A)
+end
+
 # Implementation
 # --------------
 function exponential!(A::AbstractMatrix, expA::AbstractMatrix, alg::ExponentialViaLA)
@@ -31,20 +43,24 @@ end
 
 function exponential!(A::AbstractMatrix, expA::AbstractMatrix, alg::ExponentialViaEigh)
     D, V = eigh_full(A, alg.eigh_alg)
-    copyto!(expA, V * Diagonal(exp.(diagview(D))) * inv(V))
+    iV = inv(V)
+    map!(exp, diagview(D), diagview(D))
+    mul!(expA, rmul!(V, D), iV)
     return expA
 end
 
 function exponential!(A::AbstractMatrix, expA::AbstractMatrix, alg::ExponentialViaEig)
     D, V = eig_full(A, alg.eig_alg)
-    copyto!(expA, V * Diagonal(exp.(diagview(D))) * inv(V))
+    iV = inv(V)
+    map!(exp, diagview(D), diagview(D))
+    mul!(expA, rmul!(V, D), iV)
     return expA
 end
 
 # Diagonal logic
 # --------------
 function exponential!(A::Diagonal, expA, alg::DiagonalAlgorithm)
     check_input(exponential!, A, expA, alg)
-    copyto!(expA, Diagonal(LinearAlgebra.exp.(diagview(A))))
+    map!(exp, diagview(expA), diagview(A))
     return expA
 end

src/interface/decompositions.jl

@@ -325,16 +325,15 @@ const ROCSOLVER_SVDAlgorithm = Union{ROCSOLVER_QRIteration, ROCSOLVER_Jacobi}
 """
     ExponentialViaLA()
 
-Algorithm type to denote finding the LQ decomposition of `A` by computing the QR decomposition of `Aᵀ`.
-The `qr_alg` specifies which QR-decomposition implementation to use.
+Algorithm type to denote finding the exponential of `A` via the implementation of `LinearAlgebra`.
 """
 @algdef ExponentialViaLA
 
 """
     ExponentialViaEigh()
 
-Algorithm type to denote finding the LQ decomposition of `A` by computing the QR decomposition of `Aᵀ`.
-The `qr_alg` specifies which QR-decomposition implementation to use.
+Algorithm type to denote finding the exponential `A` by computing the hermitian eigendecomposition of `A`.
+The `eigh_alg` specifies which hermitian eigendecomposition implementation to use.
 """
 struct ExponentialViaEigh{A <: AbstractAlgorithm} <: AbstractAlgorithm
     eigh_alg::A
@@ -348,8 +347,8 @@ end
 """
     ExponentialViaEig()
 
-Algorithm type to denote finding the LQ decomposition of `A` by computing the QR decomposition of `Aᵀ`.
-The `qr_alg` specifies which QR-decomposition implementation to use.
+Algorithm type to denote finding the exponential `A` by computing the eigendecomposition of `A`.
+The `eig_alg` specifies which eigendecomposition implementation to use.
 """
 struct ExponentialViaEig{A <: AbstractAlgorithm} <: AbstractAlgorithm
     eig_alg::A

src/interface/exponential.jl

@@ -1,7 +1,6 @@
 # Exponential functions
 # --------------
 @functiondef exponential
-# @algdef exponential!
 
 # Algorithm selection
 # -------------------

test/exponential.jl

@@ -6,9 +6,9 @@ using MatrixAlgebraKit: diagview
 using LinearAlgebra
 
 BLASFloats = (Float32, Float64, ComplexF32, ComplexF64)
-GenericFloats = (Float16, BigFloat, Complex{BigFloat})
+GenericFloats = (Float16, ComplexF16, BigFloat, Complex{BigFloat})
 
-@testset "exp! for T = $T" for T in BLASFloats
+@testset "exponential! for T = $T" for T in BLASFloats
     rng = StableRNG(123)
     m = 2
 
@@ -30,7 +30,7 @@ GenericFloats = (Float16, BigFloat, Complex{BigFloat})
     end
 end
 
-@testset "svd for Diagonal{$T}" for T in (BLASFloats..., GenericFloats...)
+@testset "exponential! for Diagonal{$T}" for T in (BLASFloats..., GenericFloats...)
     rng = StableRNG(123)
     atol = sqrt(eps(real(T)))
     m = 54