Why is the return type of `eigen` a `Union`?

[albertomercurio]

Hello,

I see a small Union type when diagonalizing a matrix

julia> A = rand(10, 10);

julia> eigen(A);

julia> @code_warntype eigen(A)
MethodInstance for eigen(::Matrix{Float64})
  from eigen(A::AbstractMatrix{T}; permute, scale, sortby) where T @ LinearAlgebra ~/.julia/juliaup/julia-1.12.2+0.x64.linux.gnu/share/julia/stdlib/v1.12/LinearAlgebra/src/eigen.jl:238
Static Parameters
  T = Float64
Arguments
  #self#::Core.Const(LinearAlgebra.eigen)
  A::Matrix{Float64}
Body::Union{Eigen{ComplexF64, ComplexF64, Matrix{ComplexF64}, Vector{ComplexF64}}, Eigen{Float64, Float64, Matrix{Float64}, Vector{Float64}}}
1 ─ %1 = LinearAlgebra.:(var"#eigen#102")::Core.Const(LinearAlgebra.var"#eigen#102")
│   %2 = LinearAlgebra.eigsortby::Core.Const(LinearAlgebra.eigsortby)
│   %3 = (%1)(true, true, %2, #self#, A)::Union{Eigen{ComplexF64, ComplexF64, Matrix{ComplexF64}, Vector{ComplexF64}}, Eigen{Float64, Float64, Matrix{Float64}, Vector{Float64}}}
└──      return %3

I think this happens because it is checking if the matrix is Hermitian at runtime here

LinearAlgebra.jl/src/eigen.jl

Lines 247 to 254 in 189d49d

	function _eigen(A::AbstractMatrix{T}; permute=true, scale=true, sortby=eigsortby) where {T}
	isdiag(A) && return eigen(Diagonal{eigtype(T)}(diag(A)); sortby)
	if ishermitian(A)
	eigen!(eigencopy_oftype(Hermitian(A), eigtype(T)); sortby)
	else
	eigen!(eigencopy_oftype(A, eigtype(T)); permute, scale, sortby)
	end
	end

And this returns the Union type between real or complex numbers.

I was wondering if one should remove that check, and implement specific methods for Hermitian or Symmetric matrices. If so, I can do a PR.

[dkarrasch]

That method exists. It's the one that gets called in line 250 above. So, in case surrounding code knows by type that the argument is Hermitian, then there is no type instability. The type instability occurs only for generic matrices.

[albertomercurio]

Ok, should one just return complex numbers for the general case?

[dkarrasch]

I'm sure this has been written on purpose, and that there's been discussion around that. If you know your matrix will never be Hermitian (and not rely on optimized methods in the backend), you could create a copy yourself and directly call eigen!.

[araujoms]

I think it should indeed just return complex numbers for the general case, but this has been discussed to death already and I don't think it is likely to change.