Skip to content

Commit

Permalink
Copy parent in copy_similar for Symmetric/Hermitian (#54473)
Browse files Browse the repository at this point in the history
This should be equivalent to copying the `Symmetric`/`Hermitian` matrix,
as only one triangular half is used. However, this reduces latency when
the parent is a `Matrix` (and possibly others too, although I've not
checked this).

```julia
julia> using LinearAlgebra

julia> A = rand(2,2); H = Hermitian(A);

julia> @time eigen(H);
  0.342207 seconds (440.61 k allocations: 22.926 MiB, 99.95% compilation time) # nightly
  0.260913 seconds (326.62 k allocations: 16.926 MiB, 99.93% compilation time) # This PR
```
Copying the parent is also marginally faster, although this doesn't
really matter in `eigen`:
```julia
julia> A = rand(1000,1000); H = Hermitian(A);

julia> @Btime LinearAlgebra.copy_similar($H, Float64);
  1.839 ms (3 allocations: 7.63 MiB) # nightly
  1.760 ms (3 allocations: 7.63 MiB) # This PR
```

---------

Co-authored-by: Daniel Karrasch <daniel.karrasch@posteo.de>
  • Loading branch information
jishnub and dkarrasch authored May 20, 2024
1 parent b5e91af commit 047e699
Showing 1 changed file with 3 additions and 2 deletions.
5 changes: 3 additions & 2 deletions stdlib/LinearAlgebra/src/symmetriceigen.jl
Original file line number Diff line number Diff line change
@@ -1,8 +1,9 @@
# This file is a part of Julia. License is MIT: https://julialang.org/license

# preserve HermOrSym wrapper
eigencopy_oftype(A::Hermitian, S) = Hermitian(copy_similar(A, S), sym_uplo(A.uplo))
eigencopy_oftype(A::Symmetric, S) = Symmetric(copy_similar(A, S), sym_uplo(A.uplo))
# Call `copytrito!` instead of `copy_similar` to only copy the matching triangular half
eigencopy_oftype(A::Hermitian, S) = Hermitian(copytrito!(similar(parent(A), S, size(A)), A.data, A.uplo), sym_uplo(A.uplo))
eigencopy_oftype(A::Symmetric, S) = Symmetric(copytrito!(similar(parent(A), S, size(A)), A.data, A.uplo), sym_uplo(A.uplo))

# Eigensolvers for symmetric and Hermitian matrices
eigen!(A::RealHermSymComplexHerm{<:BlasReal,<:StridedMatrix}; sortby::Union{Function,Nothing}=nothing) =
Expand Down

0 comments on commit 047e699

Please sign in to comment.