eig.jl

# Eig functions
# -------------
# TODO: kwargs for sorting eigenvalues?

docs_eig_note = """
    Note that [`eig_full`](@ref) and its variants do not assume any symmetry structure on
    the input matrices, and therefore will always return complex eigenvalues and eigenvectors
    for reasons of type stability.  For the  eigenvalue decomposition of symmetric or hermitian
    operators, see [`eigh_full`](@ref).
"""

"""
    eig_full(A; kwargs...) -> D, V
    eig_full(A, alg::AbstractAlgorithm) -> D, V
    eig_full!(A, [DV]; kwargs...) -> D, V
    eig_full!(A, [DV], alg::AbstractAlgorithm) -> D, V

Compute the full eigenvalue decomposition of the square matrix `A`,
such that `A * V = V * D`, where the invertible matrix `V` contains the eigenvectors
and the diagonal matrix `D` contains the associated eigenvalues.

!!! note
    The bang method `eig_full!` optionally accepts the output structure and
    possibly destroys the input matrix `A`. Always use the return value of the function
    as it may not always be possible to use the provided `DV` as output.

!!! note
$(docs_eig_note)

See also [`eig_vals(!)`](@ref eig_vals), [`eig_trunc_no_error`](@ref eig_trunc_no_error)
and [`eig_trunc(!)`](@ref eig_trunc).
"""
@functiondef eig_full

"""
    eig_trunc(A; [trunc], kwargs...) -> D, V, ϵ
    eig_trunc(A, alg::AbstractAlgorithm) -> D, V, ϵ
    eig_trunc!(A, [DV]; [trunc], kwargs...) -> D, V, ϵ
    eig_trunc!(A, [DV], alg::AbstractAlgorithm) -> D, V, ϵ

Compute a partial or truncated eigenvalue decomposition of the matrix `A`,
such that `A * V ≈ V * D`, where the (possibly rectangular) matrix `V` contains 
a subset of eigenvectors and the diagonal matrix `D` contains the associated eigenvalues,
selected according to a truncation strategy.

The function also returns `ϵ`, the truncation error defined as the 2-norm of the 
discarded eigenvalues.

## Truncation
The truncation strategy can be controlled via the `trunc` keyword argument. This can be
either a `NamedTuple` or a [`TruncationStrategy`](@ref). If `trunc` is not provided or
nothing, all values will be kept.

### `trunc::NamedTuple`
The supported truncation keyword arguments are:

$docs_truncation_kwargs

### `trunc::TruncationStrategy`
For more control, a truncation strategy can be supplied directly.
By default, MatrixAlgebraKit supplies the following:

$docs_truncation_strategies

## Keyword Arguments
Other keyword arguments are passed to the algorithm selection procedure. If no explicit
`alg` is provided, these keywords are used to select and configure the algorithm through
[`MatrixAlgebraKit.select_algorithm`](@ref). The remaining keywords after algorithm
selection are passed to the algorithm constructor. See [`MatrixAlgebraKit.default_algorithm`](@ref)
for the default algorithm selection behavior.

When `alg` is a [`TruncatedAlgorithm`](@ref), the `trunc` keyword cannot be specified as the
truncation strategy is already embedded in the algorithm.

!!! note
    The bang method `eig_trunc!` optionally accepts the output structure and
    possibly destroys the input matrix `A`. Always use the return value of the function
    as it may not always be possible to use the provided `DV` as output.

!!! note
$docs_eig_note

See also [`eig_full(!)`](@ref eig_full), [`eig_vals(!)`](@ref eig_vals),
[`eig_trunc_no_error!`](@ref eig_trunc_no_error) and [Truncations](@ref)
for more information on truncation strategies.
"""
@functiondef eig_trunc

"""
    eig_trunc_no_error(A; [trunc], kwargs...) -> D, V
    eig_trunc_no_error(A, alg::AbstractAlgorithm) -> D, V
    eig_trunc_no_error!(A, [DV]; [trunc], kwargs...) -> D, V
    eig_trunc_no_error!(A, [DV], alg::AbstractAlgorithm) -> D, V

Compute a partial or truncated eigenvalue decomposition of the matrix `A`,
such that `A * V ≈ V * D`, where the (possibly rectangular) matrix `V` contains 
a subset of eigenvectors and the diagonal matrix `D` contains the associated eigenvalues,
selected according to a truncation strategy. The truncation error is *not* returned.

## Truncation
The truncation strategy can be controlled via the `trunc` keyword argument. This can be
either a `NamedTuple` or a [`TruncationStrategy`](@ref). If `trunc` is not provided or
nothing, all values will be kept.

### `trunc::NamedTuple`
The supported truncation keyword arguments are:

$docs_truncation_kwargs

### `trunc::TruncationStrategy`
For more control, a truncation strategy can be supplied directly.
By default, MatrixAlgebraKit supplies the following:

$docs_truncation_strategies

## Keyword Arguments
Other keyword arguments are passed to the algorithm selection procedure. If no explicit
`alg` is provided, these keywords are used to select and configure the algorithm through
[`MatrixAlgebraKit.select_algorithm`](@ref). The remaining keywords after algorithm
selection are passed to the algorithm constructor. See [`MatrixAlgebraKit.default_algorithm`](@ref)
for the default algorithm selection behavior.

When `alg` is a [`TruncatedAlgorithm`](@ref), the `trunc` keyword cannot be specified as the
truncation strategy is already embedded in the algorithm.

!!! note
    The bang method `eig_trunc!` optionally accepts the output structure and
    possibly destroys the input matrix `A`. Always use the return value of the function
    as it may not always be possible to use the provided `DV` as output.

!!! note
$docs_eig_note

See also [`eig_full(!)`](@ref eig_full), [`eig_vals(!)`](@ref eig_vals),
[`eig_trunc(!)`](@ref eig_trunc) and [Truncations](@ref) for more
information on truncation strategies.
"""
@functiondef eig_trunc_no_error

"""
    eig_vals(A; kwargs...) -> D
    eig_vals(A, alg::AbstractAlgorithm) -> D
    eig_vals!(A, [D]; kwargs...) -> D
    eig_vals!(A, [D], alg::AbstractAlgorithm) -> D

Compute the list of eigenvalues of `A`.

!!! note
    The bang method `eig_vals!` optionally accepts the output structure and
    possibly destroys the input matrix `A`. Always use the return value of the function
    as it may not always be possible to use the provided `D` as output.

!!! note
$(docs_eig_note)

See also [`eig_full(!)`](@ref eig_full) and [`eig_trunc(!)`](@ref eig_trunc).
"""
@functiondef eig_vals

# Algorithm selection
# -------------------
default_eig_algorithm(A; kwargs...) = default_eig_algorithm(typeof(A); kwargs...)
default_eig_algorithm(T::Type; kwargs...) = throw(MethodError(default_eig_algorithm, (T,)))
function default_eig_algorithm(::Type{T}; kwargs...) where {T <: YALAPACK.MaybeBlasVecOrMat}
    return LAPACK_Expert(; kwargs...)
end
function default_eig_algorithm(::Type{T}; kwargs...) where {T <: Diagonal}
    return DiagonalAlgorithm(; kwargs...)
end
function default_eig_algorithm(::Type{<:Base.ReshapedArray{T, N, A}}) where {T, N, A}
    return default_eig_algorithm(A)
end
function default_eig_algorithm(::Type{SubArray{T, N, A}}) where {T, N, A}
    return default_eig_algorithm(A)
end

for f in (:eig_full!, :eig_vals!)
    @eval function default_algorithm(::typeof($f), ::Type{A}; kwargs...) where {A}
        return default_eig_algorithm(A; kwargs...)
    end
end

for f in (:eig_trunc!, :eig_trunc_no_error!)
    @eval function select_algorithm(::typeof($f), A, alg; trunc = nothing, kwargs...)
        if alg isa TruncatedAlgorithm
            isnothing(trunc) ||
                throw(ArgumentError("`trunc` can't be specified when `alg` is a `TruncatedAlgorithm`"))
            return alg
        else
            alg_eig = select_algorithm(eig_full!, A, alg; kwargs...)
            return TruncatedAlgorithm(alg_eig, select_truncation(trunc))
        end
    end
end