From 3767d5d67a39ef064c0c1458ba57a53a69804054 Mon Sep 17 00:00:00 2001 From: Gertian Date: Mon, 19 Aug 2024 12:19:30 +0300 Subject: [PATCH 1/8] Define + for LocalOperators I defined additional for Localoperators --- src/operators/localoperator.jl | 14 ++++++++++++++ 1 file changed, 14 insertions(+) diff --git a/src/operators/localoperator.jl b/src/operators/localoperator.jl index ffe9e34a..b9473c17 100644 --- a/src/operators/localoperator.jl +++ b/src/operators/localoperator.jl @@ -57,3 +57,17 @@ function Base.repeat(O::LocalOperator, m::Int, n::Int) end return LocalOperator(lattice, terms...) end + +function Base.:+(O1::LocalOperator, O2::LocalOperator) + if O1.lattice != O2.lattice + throw("lattices should match") + end + terms = [] + for (inds, operator) in O1.terms + push!(terms, inds => operator) + end + for (inds, operator) in O2.terms + push!(terms, inds => operator) + end + return LocalOperator(O1.lattice, terms...) +end From 82b4f3be7f3d6a9f8ea87d5bc587f98239d10e4b Mon Sep 17 00:00:00 2001 From: Gertian Date: Tue, 20 Aug 2024 13:46:58 +0300 Subject: [PATCH 2/8] Update localoperator.jl 1) Overloaded checklattice for two LocalOperators 2)Changed LocalOperator+LocalOperator such that its more efficient when duplicate indices appear 3)Defined Number*LocalOperator 4)Defined - using + and * --- src/operators/localoperator.jl | 27 +++++++++++++++++++-------- 1 file changed, 19 insertions(+), 8 deletions(-) diff --git a/src/operators/localoperator.jl b/src/operators/localoperator.jl index b9473c17..272b94ca 100644 --- a/src/operators/localoperator.jl +++ b/src/operators/localoperator.jl @@ -27,6 +27,9 @@ while the second version throws an error if the lattices do not match. function checklattice(args...) return checklattice(Bool, args...) || throw(ArgumentError("Lattice mismatch.")) end +function checklattice(::Type{Bool}, H1::LocalOperator, H2::LocalOperator) + return H1.lattice == H2.lattice +end function checklattice(::Type{Bool}, peps::InfinitePEPS, O::LocalOperator) return size(peps) == size(O.lattice) end @@ -58,16 +61,24 @@ function Base.repeat(O::LocalOperator, m::Int, n::Int) return LocalOperator(lattice, terms...) end +function Base.:*(α::Number, O2::LocalOperator) + return LocalOperator(O2.lattice, map(t -> (t[1] => α * t[2]), O2.terms)...) +end + function Base.:+(O1::LocalOperator, O2::LocalOperator) - if O1.lattice != O2.lattice - throw("lattices should match") - end - terms = [] - for (inds, operator) in O1.terms - push!(terms, inds => operator) - end + checklattice(O1, O2) + terms = [O1.terms...] for (inds, operator) in O2.terms - push!(terms, inds => operator) + i = findfirst(existing_inds -> existing_inds == inds, map(first, terms)) + if !isnothing(i) + terms[i] = (inds => terms[i][2] + operator) + else + push!(terms, inds => operator) + end end return LocalOperator(O1.lattice, terms...) end + +function Base.:-(O1::LocalOperator, O2::LocalOperator) + return O1 + (-1) * O2 +end From 3049fc415d68a913bbbd935ea094c7178a631c90 Mon Sep 17 00:00:00 2001 From: Gertian Date: Tue, 20 Aug 2024 17:53:50 +0300 Subject: [PATCH 3/8] + - *for localOperators 1) Extended LocalOperator constructor so that it eliminates operators with norm < tolerance. 2) Used this to make + - * implementations more streamlined --- src/operators/localoperator.jl | 30 ++++++++++++++++++------------ 1 file changed, 18 insertions(+), 12 deletions(-) diff --git a/src/operators/localoperator.jl b/src/operators/localoperator.jl index 272b94ca..f438513f 100644 --- a/src/operators/localoperator.jl +++ b/src/operators/localoperator.jl @@ -5,16 +5,31 @@ struct LocalOperator{T<:Tuple,S} lattice::Matrix{S} terms::T end -function LocalOperator(lattice::Matrix{S}, terms::Pair...) where {S} +function LocalOperator(lattice::Matrix{S}, terms::Pair... ; min_norm_operators = eps()^(3/4)) where {S} lattice′ = PeriodicArray(lattice) + relevant_terms = [] for (inds, operator) in terms + # Check if the indices of the operator are valid with themselves and the lattice @assert operator isa AbstractTensorMap @assert numout(operator) == numin(operator) == length(inds) for i in 1:length(inds) @assert space(operator, i) == lattice′[inds[i]] + end + # Check if we already have an operator acting on the coordinates + i = findfirst(existing_inds -> existing_inds == inds, map(first, relevant_terms)) + if !isnothing(i) # We are adding to an existing operator + new_operator = relevant_terms[i][2] + operator + if norm(new_operator) > min_norm_operators + relevant_terms[i] = (inds => new_operator) + else + deleteat!(relevant_terms, i) + end + else # It's a new operator, add it if its norm is large enough + norm(operator) > min_norm_operators && push!(relevant_terms, inds => operator) end end - return LocalOperator{typeof(terms),S}(lattice, terms) + relevant_terms = Tuple(relevant_terms) + return LocalOperator{typeof(relevant_terms), S}(lattice, relevant_terms) end """ @@ -67,16 +82,7 @@ end function Base.:+(O1::LocalOperator, O2::LocalOperator) checklattice(O1, O2) - terms = [O1.terms...] - for (inds, operator) in O2.terms - i = findfirst(existing_inds -> existing_inds == inds, map(first, terms)) - if !isnothing(i) - terms[i] = (inds => terms[i][2] + operator) - else - push!(terms, inds => operator) - end - end - return LocalOperator(O1.lattice, terms...) + return LocalOperator(O1.lattice, O1.terms..., O2.terms...) end function Base.:-(O1::LocalOperator, O2::LocalOperator) From 84b638f1ceb827bef2c54270538f9b89cf4f34f3 Mon Sep 17 00:00:00 2001 From: lkdvos Date: Wed, 21 Aug 2024 10:19:19 +0200 Subject: [PATCH 4/8] Add docstring LocalOperator --- src/operators/localoperator.jl | 23 +++++++++++++++++++++++ 1 file changed, 23 insertions(+) diff --git a/src/operators/localoperator.jl b/src/operators/localoperator.jl index f438513f..e5bad1b5 100644 --- a/src/operators/localoperator.jl +++ b/src/operators/localoperator.jl @@ -1,6 +1,29 @@ # Hamiltonian consisting of local terms # ------------------------------------- +""" + struct LocalOperator{T<:Tuple,S} + +A sum of local operators acting on a lattice. The lattice is stored as a matrix of vector spaces, +and the terms are stored as a tuple of pairs of indices and operators. + +# Fields + +- `lattice::Matrix{S}`: The lattice on which the operator acts. +- `terms::T`: The terms of the operator, stored as a tuple of pairs of indices and operators. + +# Constructors + + LocalOperator(lattice::Matrix{S}, terms::Pair...) + LocalOperator{T,S}(lattice::Matrix{S}, terms::T) where {T,S} # expert mode + +# Examples + +```julia +lattice = fill(ℂ^2, 1, 1) # single-site unitcell +O1 = LocalOperator(lattice, ((1, 1),) => σx, ((1, 1), (1, 2)) => σx ⊗ σx, ((1, 1), (2, 1)) => σx ⊗ σx) +``` +""" struct LocalOperator{T<:Tuple,S} lattice::Matrix{S} terms::T From 2c9abb9db9283765601d19816673a85d53ffdf25 Mon Sep 17 00:00:00 2001 From: lkdvos Date: Wed, 21 Aug 2024 10:19:53 +0200 Subject: [PATCH 5/8] Update LocalOperator constructor --- src/operators/localoperator.jl | 47 ++++++++++++++++++---------------- 1 file changed, 25 insertions(+), 22 deletions(-) diff --git a/src/operators/localoperator.jl b/src/operators/localoperator.jl index e5bad1b5..4b2fd29f 100644 --- a/src/operators/localoperator.jl +++ b/src/operators/localoperator.jl @@ -27,32 +27,35 @@ O1 = LocalOperator(lattice, ((1, 1),) => σx, ((1, 1), (1, 2)) => σx ⊗ σx, ( struct LocalOperator{T<:Tuple,S} lattice::Matrix{S} terms::T -end -function LocalOperator(lattice::Matrix{S}, terms::Pair... ; min_norm_operators = eps()^(3/4)) where {S} - lattice′ = PeriodicArray(lattice) - relevant_terms = [] - for (inds, operator) in terms + function LocalOperator{T,S}(lattice::Matrix{S}, terms::T) where {T,S} + plattice = PeriodicArray(lattice) # Check if the indices of the operator are valid with themselves and the lattice - @assert operator isa AbstractTensorMap - @assert numout(operator) == numin(operator) == length(inds) - for i in 1:length(inds) - @assert space(operator, i) == lattice′[inds[i]] - end - # Check if we already have an operator acting on the coordinates - i = findfirst(existing_inds -> existing_inds == inds, map(first, relevant_terms)) - if !isnothing(i) # We are adding to an existing operator - new_operator = relevant_terms[i][2] + operator - if norm(new_operator) > min_norm_operators - relevant_terms[i] = (inds => new_operator) - else - deleteat!(relevant_terms, i) + for (inds, operator) in terms + @assert operator isa AbstractTensorMap + @assert numout(operator) == numin(operator) == length(inds) + @assert spacetype(operator) == S + + for i in 1:length(inds) + @assert space(operator, i) == plattice[inds[i]] end - else # It's a new operator, add it if its norm is large enough - norm(operator) > min_norm_operators && push!(relevant_terms, inds => operator) end + return new{T,S}(lattice, terms) + end +end +function LocalOperator( + lattice::Matrix, terms::Pair...; atol=maximum(x -> eps(scalartype(x[2]))^(3 / 4), terms) +) + allinds = getindex.(terms, 1) + alloperators = getindex.(terms, 2) + + relevant_terms = [] + for inds in unique(allinds) + operator = sum(alloperators[findall(==(inds), allinds)]) + norm(operator) > atol && push!(relevant_terms, inds => operator) end - relevant_terms = Tuple(relevant_terms) - return LocalOperator{typeof(relevant_terms), S}(lattice, relevant_terms) + + terms_tuple = Tuple(relevant_terms) + return LocalOperator{typeof(terms_tuple),eltype(lattice)}(lattice, terms_tuple) end """ From 44bfe069fae90db58a824b9ed27a5282ee2a7eb9 Mon Sep 17 00:00:00 2001 From: lkdvos Date: Wed, 21 Aug 2024 10:20:02 +0200 Subject: [PATCH 6/8] Formatter --- src/operators/localoperator.jl | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/src/operators/localoperator.jl b/src/operators/localoperator.jl index 4b2fd29f..32931611 100644 --- a/src/operators/localoperator.jl +++ b/src/operators/localoperator.jl @@ -1,4 +1,3 @@ - # Hamiltonian consisting of local terms # ------------------------------------- """ @@ -107,7 +106,7 @@ function Base.:*(α::Number, O2::LocalOperator) end function Base.:+(O1::LocalOperator, O2::LocalOperator) - checklattice(O1, O2) + checklattice(O1, O2) return LocalOperator(O1.lattice, O1.terms..., O2.terms...) end From 7d8a839c59d861c5e7f658a9a662e97ce58fc585 Mon Sep 17 00:00:00 2001 From: lkdvos Date: Wed, 21 Aug 2024 10:25:56 +0200 Subject: [PATCH 7/8] Clean-up linear algebra localoperator --- src/operators/localoperator.jl | 16 +++++++++++----- 1 file changed, 11 insertions(+), 5 deletions(-) diff --git a/src/operators/localoperator.jl b/src/operators/localoperator.jl index 32931611..8ebbe8a5 100644 --- a/src/operators/localoperator.jl +++ b/src/operators/localoperator.jl @@ -101,15 +101,21 @@ function Base.repeat(O::LocalOperator, m::Int, n::Int) return LocalOperator(lattice, terms...) end -function Base.:*(α::Number, O2::LocalOperator) - return LocalOperator(O2.lattice, map(t -> (t[1] => α * t[2]), O2.terms)...) +# Linear Algebra +# -------------- +function Base.:*(α::Number, O::LocalOperator) + scaled_terms = map(((inds, operator),) -> (inds => α * operator), O.terms) + return LocalOperator{typeof(scaled_terms),eltype(O.lattice)}(O.lattice, scaled_terms) end +Base.:*(O::LocalOperator, α::Number) = α * O + +Base.:/(O::LocalOperator, α::Number) = O * inv(α) +Base.:\(α::Number, O::LocalOperator) = inv(α) * O function Base.:+(O1::LocalOperator, O2::LocalOperator) checklattice(O1, O2) return LocalOperator(O1.lattice, O1.terms..., O2.terms...) end -function Base.:-(O1::LocalOperator, O2::LocalOperator) - return O1 + (-1) * O2 -end +Base.:-(O::LocalOperator) = -1 * O +Base.:-(O1::LocalOperator, O2::LocalOperator) = O1 + (-O2) From 51e4b96aca8885048fc355c4e51eb352ef4649bc Mon Sep 17 00:00:00 2001 From: lkdvos Date: Wed, 21 Aug 2024 10:56:20 +0200 Subject: [PATCH 8/8] Fix eps(real(scalartype))) error --- src/operators/localoperator.jl | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/src/operators/localoperator.jl b/src/operators/localoperator.jl index 8ebbe8a5..65796b88 100644 --- a/src/operators/localoperator.jl +++ b/src/operators/localoperator.jl @@ -42,7 +42,9 @@ struct LocalOperator{T<:Tuple,S} end end function LocalOperator( - lattice::Matrix, terms::Pair...; atol=maximum(x -> eps(scalartype(x[2]))^(3 / 4), terms) + lattice::Matrix, + terms::Pair...; + atol=maximum(x -> eps(real(scalartype(x[2])))^(3 / 4), terms), ) allinds = getindex.(terms, 1) alloperators = getindex.(terms, 2)