Core of Exact Diagonalization
ExactDiagonalization.EDCore.edtimer
— Constantconst edtimer = TimerOutput()
The default shared timer for all exact diagonalization methods.
ExactDiagonalization.EDCore.ED
— TypeED(
lattice::AbstractLattice,
hilbert::Hilbert,
terms::Union{Term, Tuple{Term, Vararg{Term}}},
quantumnumbers::Union{AbelianNumber, Tuple{AbelianNumber, Vararg{AbelianNumber}}},
dtype::Type{<:Number}=commontype(terms),
boundary::Boundary=plain;
neighbors::Union{Nothing, Int, Neighbors}=nothing,
timer::TimerOutput=edtimer,
)
Construct the exact diagonalization method for a quantum lattice system.
ExactDiagonalization.EDCore.ED
— TypeED(
lattice::AbstractLattice, hilbert::Hilbert, terms::Union{Term, Tuple{Term, Vararg{Term}}}, targetspace::TargetSpace=TargetSpace(hilbert), dtype::Type{<:Number}=commontype(terms), boundary::Boundary=plain;
neighbors::Union{Nothing, Int, Neighbors}=nothing, timer::TimerOutput=edtimer
)
Construct the exact diagonalization method for a quantum lattice system.
ExactDiagonalization.EDCore.ED
— TypeED{K<:EDKind, L<:AbstractLattice, G<:OperatorGenerator, M<:Image} <: Frontend
Exact diagonalization method of a quantum lattice system.
ExactDiagonalization.EDCore.EDEigen
— TypeEDEigen{V<:Number, T<:Number, S<:Sector} <: Factorization{T}
Eigen decomposition in exact diagonalization method.
Compared to the usual eigen decomposition Eigen
, EDEigen
contains a :sectors
attribute to store the sectors of Hilbert space in which the eigen values and eigen vectors are computed. Furthermore, given that in different sectors the dimensions of the sub-Hilbert spaces can also be different, the :vectors
attribute of EDEigen
is a vector of vector instead of a matrix.
ExactDiagonalization.EDCore.EDKind
— TypeEDKind{K}
The kind of the exact diagonalization method applied to a quantum lattice system.
ExactDiagonalization.EDCore.EDMatrix
— TypeEDMatrix{S<:Sector, M<:SparseMatrixCSC} <: OperatorPack{M, Tuple{S, S}}
Matrix representation of quantum operators between a ket and a bra Hilbert space.
ExactDiagonalization.EDCore.EDMatrix
— MethodEDMatrix(sector::Sector, m::SparseMatrixCSC)
EDMatrix(braket::NTuple{2, Sector}, m::SparseMatrixCSC)
Construct a matrix representation when
- the ket and bra spaces share the same bases;
2-3) the ket and bra spaces may be different.
ExactDiagonalization.EDCore.EDMatrixRepresentation
— TypeEDMatrixRepresentation{D<:Number, S<:Sector, T} <: MatrixRepresentation
Exact matrix representation of a quantum lattice system on a target Hilbert space.
ExactDiagonalization.EDCore.EDMatrixRepresentation
— MethodEDMatrixRepresentation{D}(target::TargetSpace, table) where {D<:Number}
Construct a exact matrix representation.
ExactDiagonalization.EDCore.Sector
— Typeabstract type Sector <: OperatorUnit
A sector of the Hilbert space which forms the bases of an irreducible representation of the Hamiltonian of a quantum lattice system.
ExactDiagonalization.EDCore.SectorFilter
— TypeSectorFilter{S} <: LinearTransformation
Filter the target bra and ket Hilbert spaces.
ExactDiagonalization.EDCore.TargetSpace
— TypeTargetSpace{S<:Sector} <: VectorSpace{S}
The target Hilbert space in which the exact diagonalization method is performed, which could be the direct sum of several sectors.
ExactDiagonalization.EDCore.TargetSpace
— MethodTargetSpace(hilbert::Hilbert)
TargetSpace(hilbert::Hilbert, quantumnumbers::Union{AbelianNumber, Tuple{AbelianNumber, Vararg{AbelianNumber}}})
Construct a target space from the total Hilbert space and the associated quantum numbers.
ExactDiagonalization.EDCore.TargetSpace
— MethodTargetSpace(sector::Sector, sectors::Sector...)
Construct a target space from sectors.
ExactDiagonalization.EDCore.productable
— Functionproductable(sector₁::Sector, sector₂::Sector) -> Bool
Judge whether two sectors could be direct producted.
ExactDiagonalization.EDCore.sumable
— Functionsumable(sector₁::Sector, sector₂::Sector) -> Bool
Judge whether two sectors could be direct summed.
LinearAlgebra.eigen
— Methodeigen(ed::ED, sectors::Union{AbelianNumber, Sector}...; timer::TimerOutput=edtimer, kwargs...) -> EDEigen
eigen(ed::Algorithm{<:ED}, sectors::Union{AbelianNumber, Sector}...; kwargs...) -> EDEigen
Solve the eigen problem by the restarted Lanczos method provided by the Arpack package.
LinearAlgebra.eigen
— Methodeigen(m::EDMatrix; nev=6, which=:SR, tol=0.0, maxiter=300, sigma=nothing, v₀=dtype(m)[]) -> Eigen
Solve the eigen problem by the restarted Lanczos method provided by the Arpack package.
LinearAlgebra.eigen
— Methodeigen(ms::OperatorSum{<:EDMatrix}; nev::Int=1, tol::Real=0.0, maxiter::Int=300, v₀::Union{AbstractVector, Dict{<:Sector, <:AbstractVector}, Dict{<:AbelianNumber, <:AbstractVector}}=dtype(eltype(ms))[], timer::TimerOutput=edtimer)
Solve the eigen problem by the restarted Lanczos method provided by the Arpack package.
QuantumLattices.:⊕
— Method⊕(sector::Sector, sectors::Union{Sector, TargetSpace}...) -> TargetSpace
⊕(target::TargetSpace, sectors::Union{Sector, TargetSpace}...) -> TargetSpace
Get the direct sum of sectors and target spaces.
QuantumLattices.QuantumOperators.matrix
— Functionmatrix(ops::Operators, braket::NTuple{2, Sector}, table, dtype=valtype(eltype(ops))) -> SparseMatrixCSC{dtype, Int}
Get the CSC-formed sparse matrix representation of a set of operators.
Here, table
specifies the order of the operator ids.
QuantumLattices.QuantumOperators.matrix
— Methodmatrix(ed::ED, sectors::Union{AbelianNumber, Sector}...; timer::TimerOutput=edtimer, kwargs...) -> OperatorSum{<:EDMatrix}
matrix(ed::Algorithm{<:ED}, sectors::Union{AbelianNumber, Sector}...; kwargs...) -> OperatorSum{<:EDMatrix}
Get the sparse matrix representation of a quantum lattice system in the target space.