From 2fdcc04d948660ece5c27d0225add06c66f62a10 Mon Sep 17 00:00:00 2001 From: Yi-Te Huang Date: Tue, 18 Feb 2025 13:10:24 +0900 Subject: [PATCH] introduce `groundstate` --- CHANGELOG.md | 2 + docs/src/resources/api.md | 1 + .../QuantumObject/QuantumObject_functions.md | 1 + src/qobj/eigsolve.jl | 53 +++++++++++++++++-- test/core-test/eigenvalues_and_operators.jl | 10 ++++ 5 files changed, 62 insertions(+), 5 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index b2dcd5593..64e7a0979 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -10,6 +10,7 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 - Support for single `AbstractQuantumObject` in `sc_ops` for faster specific method in `ssesolve` and `smesolve`. ([#408]) - Change save callbacks from `PresetTimeCallback` to `FunctionCallingCallback`. ([#410]) - Align `eigenstates` and `eigenenergies` to QuTiP. ([#411]) +- Introduce `groundstate`. ([#412]) ## [v0.27.0] Release date: 2025-02-14 @@ -147,3 +148,4 @@ Release date: 2024-11-13 [#408]: https://github.com/qutip/QuantumToolbox.jl/issues/408 [#410]: https://github.com/qutip/QuantumToolbox.jl/issues/410 [#411]: https://github.com/qutip/QuantumToolbox.jl/issues/411 +[#412]: https://github.com/qutip/QuantumToolbox.jl/issues/412 diff --git a/docs/src/resources/api.md b/docs/src/resources/api.md index a6f5389a9..13a091d88 100644 --- a/docs/src/resources/api.md +++ b/docs/src/resources/api.md @@ -90,6 +90,7 @@ partial_transpose EigsolveResult eigenenergies eigenstates +groundstate LinearAlgebra.eigen LinearAlgebra.eigvals eigsolve diff --git a/docs/src/users_guide/QuantumObject/QuantumObject_functions.md b/docs/src/users_guide/QuantumObject/QuantumObject_functions.md index ac2c62698..0cbe1b0f4 100644 --- a/docs/src/users_guide/QuantumObject/QuantumObject_functions.md +++ b/docs/src/users_guide/QuantumObject/QuantumObject_functions.md @@ -43,6 +43,7 @@ Here is a table that summarizes all the supported linear algebra functions and a - [`eigenenergies`](@ref): return eigenenergies (eigenvalues) - [`eigenstates`](@ref): return [`EigsolveResult`](@ref) (contains eigenvalues and eigenvectors) +- [`groundstate`](@ref): return the ground state eigenvalue and corresponding eigenvector - [`eigvals`](@ref): return eigenvalues - [`eigen`](@ref): using dense eigen solver and return [`EigsolveResult`](@ref) (contains eigenvalues and eigenvectors) - [`eigsolve`](@ref): using sparse eigen solver and return [`EigsolveResult`](@ref) (contains eigenvalues and eigenvectors) diff --git a/src/qobj/eigsolve.jl b/src/qobj/eigsolve.jl index f180c0272..b6a1c2628 100644 --- a/src/qobj/eigsolve.jl +++ b/src/qobj/eigsolve.jl @@ -3,7 +3,7 @@ Eigen solvers and results for QuantumObject =# export EigsolveResult -export eigenenergies, eigenstates, eigsolve +export eigenenergies, eigenstates, groundstate, eigsolve export eigsolve_al @doc raw""" @@ -87,13 +87,19 @@ function Base.getproperty(res::EigsolveResult, key::Symbol) end end +_eigenvector_type(::Type{OperatorQuantumObject}) = Ket +_eigenvector_type(::Type{SuperOperatorQuantumObject}) = OperatorKet + Base.iterate(res::EigsolveResult) = (res.values, Val(:vector_list)) Base.iterate(res::EigsolveResult{T1,T2,Nothing}, ::Val{:vector_list}) where {T1,T2} = ([res.vectors[:, k] for k in 1:length(res.values)], Val(:vectors)) -Base.iterate(res::EigsolveResult{T1,T2,OperatorQuantumObject}, ::Val{:vector_list}) where {T1,T2} = - ([QuantumObject(res.vectors[:, k], Ket, res.dimensions) for k in 1:length(res.values)], Val(:vectors)) -Base.iterate(res::EigsolveResult{T1,T2,SuperOperatorQuantumObject}, ::Val{:vector_list}) where {T1,T2} = - ([QuantumObject(res.vectors[:, k], OperatorKet, res.dimensions) for k in 1:length(res.values)], Val(:vectors)) +Base.iterate( + res::EigsolveResult{T1,T2,OpType}, + ::Val{:vector_list}, +) where {T1,T2,OpType<:Union{OperatorQuantumObject,SuperOperatorQuantumObject}} = ( + [QuantumObject(res.vectors[:, k], _eigenvector_type(OpType), res.dimensions) for k in 1:length(res.values)], + Val(:vectors), +) Base.iterate(res::EigsolveResult, ::Val{:vectors}) = (res.vectors, Val(:done)) Base.iterate(res::EigsolveResult, ::Val{:done}) = nothing @@ -538,3 +544,40 @@ function eigenstates( return eigsolve(A; kwargs...) end end + +@doc raw""" + groundstate(A::QuantumObject; safe::Bool=true, tol::Real=1e-8, kwargs...) + +Calculate the ground state eigenvalue and corresponding eigenvector + +# Arguments +- `A::QuantumObject`: the [`QuantumObject`](@ref) to solve ground state eigenvalue and eigenvector +- `safe::Bool`: if `true` check for degenerate ground state. Default to `true`. +- `tol::Real`: the tolerance. Default is `1e-8`. +- `kwargs`: Additional keyword arguments passed to the solver. If `sparse=true`, the keyword arguments are passed to [`eigsolve`](@ref), otherwise to [`eigen`](@ref). + +# Returns +- `eigval::Number`: the ground state eigenvalue +- `eigvec::QuantumObject`: the ground state eigenvector +""" +function groundstate( + A::QuantumObject{OpType}; + safe::Bool = true, + tol::Real = 1e-8, + kwargs..., +) where {OpType<:Union{OperatorQuantumObject,SuperOperatorQuantumObject}} + # TODO: support for sparse eigsolve + #if !sparse + result = eigen(A; kwargs...) + #else + # eigvals = safe ? 2 : 1 # number of eigenvalues to calculate + # result = eigsolve(A; eigvals = eigvals, tol = tol, kwargs...) + #end + + # the tolarence should be less strick than the `tol` for the eigensolver + # so it's numerical errors are not seen as degenerate states. + evals = result.values + safe && (abs(evals[2] - evals[1]) < (10 * tol)) && @warn "Ground state may be degenerate." + + return evals[1], QuantumObject(result.vectors[:, 1], _eigenvector_type(OpType), result.dimensions) +end diff --git a/test/core-test/eigenvalues_and_operators.jl b/test/core-test/eigenvalues_and_operators.jl index 3faebf0dc..fdd0610f0 100644 --- a/test/core-test/eigenvalues_and_operators.jl +++ b/test/core-test/eigenvalues_and_operators.jl @@ -92,6 +92,14 @@ @test isapprox(vec2mat(vecs[1]).data * exp(-1im * angle(vecs[1][1])), vec2mat(vecs2[1]).data, atol = 1e-7) @test isapprox(vec2mat(vecs[1]).data * exp(-1im * angle(vecs[1][1])), vec2mat(state3[1]).data, atol = 1e-5) + # ground state # TODO: support for sparse eigsolve + U = rand_unitary(5) + M = U * Qobj(diagm([1, 1, 2, 3, 4])) * U' # degenerate ground state + gval_1, gvec_1 = @test_logs (:warn,) groundstate(M) + # gval_2, gvec_2 = @test_logs (:warn,) groundstate(M, sparse = true) + @test gval_1 ≈ 1# ≈ gval_2 + #@test isapprox(gvec_1, gvec_2, atol = 1e-6) + @testset "Type Inference (eigen)" begin N = 5 a = kron(destroy(N), qeye(N)) @@ -112,6 +120,8 @@ @inferred eigenstates(H, sparse = false) @inferred eigenstates(H, sparse = true) @inferred eigenstates(L, sparse = true) + @inferred groundstate(H)#, sparse = false) # TODO: support for sparse eigsolve + #@inferred groundstate(H, sparse = true) @inferred eigsolve_al(L, 1 \ (40 * κ), eigvals = 10) end end