From 70dd248656250b63490808d4dba951b1ac4ad8d5 Mon Sep 17 00:00:00 2001 From: Boris De Vos Date: Fri, 26 Sep 2025 16:26:54 +0200 Subject: [PATCH 01/15] introduce `local_expectation_value` for infinite case --- src/MPSKit.jl | 2 +- src/algorithms/expval.jl | 26 ++++++++++++++++++++------ 2 files changed, 21 insertions(+), 7 deletions(-) diff --git a/src/MPSKit.jl b/src/MPSKit.jl index 51ee82e32..c2220b080 100644 --- a/src/MPSKit.jl +++ b/src/MPSKit.jl @@ -42,7 +42,7 @@ export DynamicalDMRG, NaiveInvert, Jeckelmann export exact_diagonalization, fidelity_susceptibility # toolbox: -export expectation_value, correlator, variance +export expectation_value, local_expectation_value, correlator, variance export correlation_length, marek_gap, transfer_spectrum export entropy, entanglement_spectrum export open_boundary_conditions, periodic_boundary_conditions diff --git a/src/algorithms/expval.jl b/src/algorithms/expval.jl index e8ac68738..d9ec0d89d 100644 --- a/src/algorithms/expval.jl +++ b/src/algorithms/expval.jl @@ -106,12 +106,19 @@ function expectation_value( envs::AbstractMPSEnvironments = environments(ψ, H) ) return sum(1:length(ψ)) do i - return contract_mpo_expval( - ψ.AC[i], envs.GLs[i], H[i][:, 1, 1, end], envs.GRs[i][end] - ) + return local_expectation_value(ψ, H, envs, i) end end +function local_expectation_value( + ψ::InfiniteMPS, H::InfiniteMPOHamiltonian, + envs::AbstractMPSEnvironments = environments(ψ, H), site::Int=1 + ) + return contract_mpo_expval( + ψ.AC[site], envs.GLs[site], H[site][:, 1, 1, end], envs.GRs[site][end] + ) +end + # DenseMPO # -------- function expectation_value(ψ::FiniteMPS, mpo::FiniteMPO) @@ -128,11 +135,18 @@ function expectation_value( envs::MultilineEnvironments = environments(ψ, O) ) return prod(product(1:size(ψ, 1), 1:size(ψ, 2))) do (i, j) - GL = envs[i].GLs[j] - GR = envs[i].GRs[j] - return contract_mpo_expval(ψ.AC[i, j], GL, O[i, j], GR, ψ.AC[i + 1, j]) + return local_expectation_value(ψ, O, envs, (i, j)) end end +function local_expectation_value( + ψ::MultilineMPS, O::MultilineMPO{<:InfiniteMPO}, + envs::MultilineEnvironments = environments(ψ, O), + (i, j)::Tuple{Int, Int} = size(ψ) + ) + GL = envs[i].GLs[j] + GR = envs[i].GRs[j] + return contract_mpo_expval(ψ.AC[i, j], GL, O[i, j], GR, ψ.AC[i + 1, j]) +end function expectation_value(ψ::MultilineMPS, mpo::MultilineMPO, envs...) # TODO: fix environments return prod(x -> expectation_value(x...), zip(parent(ψ), parent(mpo))) From a7a0ef2a82fc77e95c02611259073363493fb717 Mon Sep 17 00:00:00 2001 From: Boris De Vos Date: Fri, 26 Sep 2025 16:27:52 +0200 Subject: [PATCH 02/15] `(local_)expectation_value` for MPO tensor w/ infinite mps --- src/algorithms/expval.jl | 20 ++++++++++++++++++++ 1 file changed, 20 insertions(+) diff --git a/src/algorithms/expval.jl b/src/algorithms/expval.jl index d9ec0d89d..498e27e73 100644 --- a/src/algorithms/expval.jl +++ b/src/algorithms/expval.jl @@ -119,6 +119,26 @@ function local_expectation_value( ) end +# MPO tensor +#----------- +function expectation_value( + ψ::InfiniteMPS, O::MPOTensor, + envs::AbstractMPSEnvironments + ) + return sum(1:length(ψ)) do site + return local_expectation_value(ψ, O, envs, site) + end +end + +function local_expectation_value( + ψ::InfiniteMPS, O::MPOTensor, + envs::AbstractMPSEnvironments, site::Int=1 + ) + return contract_mpo_expval( + ψ.AC[site], envs.GLs[site], O, envs.GRs[site] + ) +end + # DenseMPO # -------- function expectation_value(ψ::FiniteMPS, mpo::FiniteMPO) From 7d22ffe67b190a66b22360f0e931d9a1e1b063f7 Mon Sep 17 00:00:00 2001 From: Boris De Vos Date: Mon, 29 Sep 2025 12:51:05 +0200 Subject: [PATCH 03/15] test for hamiltonians --- test/operators.jl | 3 +++ 1 file changed, 3 insertions(+) diff --git a/test/operators.jl b/test/operators.jl index 3ea1ae003..a817c3c35 100644 --- a/test/operators.jl +++ b/test/operators.jl @@ -245,8 +245,11 @@ module TestOperators e1 = expectation_value(ψ2, H1) e3 = expectation_value(ψ2, H3) + e3_1 = local_expectation_value(ψ2, H3, environments(ψ2, H3), 1) + e3_2 = local_expectation_value(ψ2, H3, environments(ψ2, H3), 2) @test e1 + e3 ≈ expectation_value(ψ2, H1 + H3) atol = 1.0e-10 + @test e3 ≈ e3_1 + e3_2 atol = 1.0e-10 H4 = H1 + H3 h4 = H4 * H4 From 3f04c193eb60fa95e38da8168ecbba91fcb7edc5 Mon Sep 17 00:00:00 2001 From: Boris De Vos Date: Mon, 29 Sep 2025 12:51:37 +0200 Subject: [PATCH 04/15] test for MPO tensors --- test/algorithms.jl | 27 +++++++++++++++++++++++++++ test/setup.jl | 32 +++++++++++++++++++++++++++++--- 2 files changed, 56 insertions(+), 3 deletions(-) diff --git a/test/algorithms.jl b/test/algorithms.jl index b009c121b..95ee8e271 100644 --- a/test/algorithms.jl +++ b/test/algorithms.jl @@ -462,6 +462,33 @@ module TestAlgorithms end end + @testset "2D infinite partition functions" verbose = true begin + beta = 0.5 # ferromagnetic phase + f_th = -2.0515856253898357 + m_th = 0.911319377877496 + e_th = -1.7455645753125533 + + alg = VOMPS(; tol=1e-8, verbosity=verbosity_conv) + O_mpo = classical_ising(; β = beta) + ψ₀ = InfiniteMPS(ℂ^2, ℂ^10) + ψ, envs = leading_boundary(ψ₀, O_mpo, alg) + + λ = expectation_value(ψ, O_mpo, envs) + f = -log(λ) / beta + @test f ≈ f_th atol = 1.0e-10 + + O = ising_bulk_tensor(beta) + denom = expectation_value(ψ, O, envs) + + M = ising_magnetization_tensor(beta) + m = local_expectation_value(ψ, M, envs) / denom + @test m ≈ m_th atol = 1.0e-10 + + E = ising_energy_tensor(beta) + e = local_expectation_value(ψ, E, envs) / denom + @test e ≈ e_th atol = 1e-2 + end + @testset "excitations" verbose = true begin @testset "infinite (ham)" begin H = repeat(force_planar(heisenberg_XXX()), 2) diff --git a/test/setup.jl b/test/setup.jl index a6b07bbc0..ca9de4df2 100644 --- a/test/setup.jl +++ b/test/setup.jl @@ -18,7 +18,8 @@ export force_planar export symm_mul_mpo export transverse_field_ising, heisenberg_XXX, bilinear_biquadratic_model, XY_model, kitaev_model -export classical_ising, finite_classical_ising, sixvertex +export classical_ising, sixvertex +export ising_bond_tensor, ising_magnetisation_tensor, ising_energy_tensor # using TensorOperations @@ -300,15 +301,40 @@ function ising_bond_tensor(β) return nt end -function classical_ising(; β = log(1 + sqrt(2)) / 2, L = Inf) +function ising_bulk_tensor(β) nt = ising_bond_tensor(β) - δbulk = zeros(ComplexF64, (2, 2, 2, 2)) δbulk[1, 1, 1, 1] = 1 δbulk[2, 2, 2, 2] = 1 @tensor obulk[-1 -2; -3 -4] := δbulk[1 2; 3 4] * nt[-1; 1] * nt[-2; 2] * nt[-3; 3] * nt[-4; 4] Obulk = TensorMap(obulk, ℂ^2 * ℂ^2, ℂ^2 * ℂ^2) + return Obulk +end + +function ising_magnetisation_tensor(β) + nt = ising_bond_tensor(β) + M = zeros(2, 2, 2, 2) + M[1, 1, 1, 1] = 1 + M[2, 2, 2, 2] = -1 + @tensor m[-1 -2; -3 -4] := M[1 2; 3 4] * nt[-1; 1] * nt[-2; 2] * nt[-3; 3] * nt[-4; 4] + return TensorMap(m, ℂ^2 * ℂ^2 ← ℂ^2 * ℂ^2) +end + +function ising_energy_tensor(β) + nt = ising_bond_tensor(β) + J = 1.0 + e = ComplexF64[-J J; J -J] .* nt + @tensor e_hor[-1 -2; -3 -4] := + O[1 2; 3 4] * nt[-1; 1] * nt[-2; 2] * nt[-3; 3] * e[-4; 4] + @tensor e_vert[-1 -2; -3 -4] := + O[1 2; 3 4] * nt[-1; 1] * nt[-2; 2] * e[-3; 3] * nt[-4; 4] + e = e_hor + e_vert + return TensorMap(e, ℂ^2 * ℂ^2 ← ℂ^2 * ℂ^2) +end + +function classical_ising(; β = log(1 + sqrt(2)) / 2, L = Inf) + Obulk = ising_bulk_tensor(β) L == Inf && return InfiniteMPO([Obulk]) From 372fab0ae1356aa3ce713a50e74ef7878d3f93d1 Mon Sep 17 00:00:00 2001 From: Boris De Vos Date: Mon, 29 Sep 2025 13:24:59 +0200 Subject: [PATCH 05/15] export for tests --- test/setup.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/test/setup.jl b/test/setup.jl index ca9de4df2..0949972b0 100644 --- a/test/setup.jl +++ b/test/setup.jl @@ -19,7 +19,7 @@ export symm_mul_mpo export transverse_field_ising, heisenberg_XXX, bilinear_biquadratic_model, XY_model, kitaev_model export classical_ising, sixvertex -export ising_bond_tensor, ising_magnetisation_tensor, ising_energy_tensor +export ising_bond_tensor, ising_bulk_tensor, ising_magnetisation_tensor, ising_energy_tensor # using TensorOperations From c07bd7c4d39e4190569176942082aa0793893125 Mon Sep 17 00:00:00 2001 From: Boris De Vos Date: Mon, 29 Sep 2025 13:50:56 +0200 Subject: [PATCH 06/15] typo --- test/algorithms.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/test/algorithms.jl b/test/algorithms.jl index 95ee8e271..590744eca 100644 --- a/test/algorithms.jl +++ b/test/algorithms.jl @@ -480,7 +480,7 @@ module TestAlgorithms O = ising_bulk_tensor(beta) denom = expectation_value(ψ, O, envs) - M = ising_magnetization_tensor(beta) + M = ising_magnetisation_tensor(beta) m = local_expectation_value(ψ, M, envs) / denom @test m ≈ m_th atol = 1.0e-10 From 5d210a745121234dd1bda332cc1c16fa3f139e60 Mon Sep 17 00:00:00 2001 From: Boris De Vos Date: Mon, 29 Sep 2025 14:16:40 +0200 Subject: [PATCH 07/15] another one --- test/setup.jl | 1 + 1 file changed, 1 insertion(+) diff --git a/test/setup.jl b/test/setup.jl index 0949972b0..eb4a9e11d 100644 --- a/test/setup.jl +++ b/test/setup.jl @@ -323,6 +323,7 @@ end function ising_energy_tensor(β) nt = ising_bond_tensor(β) + O = ising_bulk_tensor(β) J = 1.0 e = ComplexF64[-J J; J -J] .* nt @tensor e_hor[-1 -2; -3 -4] := From 36456db9e8505a8644be75d0f673e2fdeeadd2eb Mon Sep 17 00:00:00 2001 From: Boris De Vos Date: Mon, 29 Sep 2025 14:51:44 +0200 Subject: [PATCH 08/15] and another one + format --- src/algorithms/expval.jl | 4 ++-- test/algorithms.jl | 6 +++--- test/setup.jl | 2 +- 3 files changed, 6 insertions(+), 6 deletions(-) diff --git a/src/algorithms/expval.jl b/src/algorithms/expval.jl index 498e27e73..332a4de42 100644 --- a/src/algorithms/expval.jl +++ b/src/algorithms/expval.jl @@ -112,7 +112,7 @@ end function local_expectation_value( ψ::InfiniteMPS, H::InfiniteMPOHamiltonian, - envs::AbstractMPSEnvironments = environments(ψ, H), site::Int=1 + envs::AbstractMPSEnvironments = environments(ψ, H), site::Int = 1 ) return contract_mpo_expval( ψ.AC[site], envs.GLs[site], H[site][:, 1, 1, end], envs.GRs[site][end] @@ -132,7 +132,7 @@ end function local_expectation_value( ψ::InfiniteMPS, O::MPOTensor, - envs::AbstractMPSEnvironments, site::Int=1 + envs::AbstractMPSEnvironments, site::Int = 1 ) return contract_mpo_expval( ψ.AC[site], envs.GLs[site], O, envs.GRs[site] diff --git a/test/algorithms.jl b/test/algorithms.jl index 590744eca..78ae0c87c 100644 --- a/test/algorithms.jl +++ b/test/algorithms.jl @@ -468,7 +468,7 @@ module TestAlgorithms m_th = 0.911319377877496 e_th = -1.7455645753125533 - alg = VOMPS(; tol=1e-8, verbosity=verbosity_conv) + alg = VOMPS(; tol = 1.0e-8, verbosity = verbosity_conv) O_mpo = classical_ising(; β = beta) ψ₀ = InfiniteMPS(ℂ^2, ℂ^10) ψ, envs = leading_boundary(ψ₀, O_mpo, alg) @@ -479,14 +479,14 @@ module TestAlgorithms O = ising_bulk_tensor(beta) denom = expectation_value(ψ, O, envs) - + M = ising_magnetisation_tensor(beta) m = local_expectation_value(ψ, M, envs) / denom @test m ≈ m_th atol = 1.0e-10 E = ising_energy_tensor(beta) e = local_expectation_value(ψ, E, envs) / denom - @test e ≈ e_th atol = 1e-2 + @test e ≈ e_th atol = 1.0e-2 end @testset "excitations" verbose = true begin diff --git a/test/setup.jl b/test/setup.jl index eb4a9e11d..6e7e3b487 100644 --- a/test/setup.jl +++ b/test/setup.jl @@ -323,7 +323,7 @@ end function ising_energy_tensor(β) nt = ising_bond_tensor(β) - O = ising_bulk_tensor(β) + O = ising_bulk_tensor(β).data J = 1.0 e = ComplexF64[-J J; J -J] .* nt @tensor e_hor[-1 -2; -3 -4] := From eec5a46ffa7472ef2bd20782d612222edc4178a4 Mon Sep 17 00:00:00 2001 From: Boris De Vos Date: Tue, 30 Sep 2025 08:33:12 +0200 Subject: [PATCH 09/15] change `expectation_value` for mpo tensors --- src/MPSKit.jl | 2 +- src/algorithms/expval.jl | 22 +++++++--------------- test/algorithms.jl | 6 +++--- test/operators.jl | 3 --- 4 files changed, 11 insertions(+), 22 deletions(-) diff --git a/src/MPSKit.jl b/src/MPSKit.jl index c2220b080..51ee82e32 100644 --- a/src/MPSKit.jl +++ b/src/MPSKit.jl @@ -42,7 +42,7 @@ export DynamicalDMRG, NaiveInvert, Jeckelmann export exact_diagonalization, fidelity_susceptibility # toolbox: -export expectation_value, local_expectation_value, correlator, variance +export expectation_value, correlator, variance export correlation_length, marek_gap, transfer_spectrum export entropy, entanglement_spectrum export open_boundary_conditions, periodic_boundary_conditions diff --git a/src/algorithms/expval.jl b/src/algorithms/expval.jl index 332a4de42..60239aa98 100644 --- a/src/algorithms/expval.jl +++ b/src/algorithms/expval.jl @@ -106,11 +106,11 @@ function expectation_value( envs::AbstractMPSEnvironments = environments(ψ, H) ) return sum(1:length(ψ)) do i - return local_expectation_value(ψ, H, envs, i) + return _local_expectation_value(ψ, H, envs, i) end end -function local_expectation_value( +function _local_expectation_value( ψ::InfiniteMPS, H::InfiniteMPOHamiltonian, envs::AbstractMPSEnvironments = environments(ψ, H), site::Int = 1 ) @@ -122,18 +122,10 @@ end # MPO tensor #----------- function expectation_value( - ψ::InfiniteMPS, O::MPOTensor, - envs::AbstractMPSEnvironments - ) - return sum(1:length(ψ)) do site - return local_expectation_value(ψ, O, envs, site) - end -end - -function local_expectation_value( - ψ::InfiniteMPS, O::MPOTensor, - envs::AbstractMPSEnvironments, site::Int = 1 + ψ::InfiniteMPS, mpo_pair::Tuple{InfiniteMPO, Pair}, + envs::AbstractMPSEnvironments = environments(ψ, first(mpo_pair)) ) + site, O = mpo_pair[2] return contract_mpo_expval( ψ.AC[site], envs.GLs[site], O, envs.GRs[site] ) @@ -155,10 +147,10 @@ function expectation_value( envs::MultilineEnvironments = environments(ψ, O) ) return prod(product(1:size(ψ, 1), 1:size(ψ, 2))) do (i, j) - return local_expectation_value(ψ, O, envs, (i, j)) + return _local_expectation_value(ψ, O, envs, (i, j)) end end -function local_expectation_value( +function _local_expectation_value( ψ::MultilineMPS, O::MultilineMPO{<:InfiniteMPO}, envs::MultilineEnvironments = environments(ψ, O), (i, j)::Tuple{Int, Int} = size(ψ) diff --git a/test/algorithms.jl b/test/algorithms.jl index 78ae0c87c..13e505435 100644 --- a/test/algorithms.jl +++ b/test/algorithms.jl @@ -481,11 +481,11 @@ module TestAlgorithms denom = expectation_value(ψ, O, envs) M = ising_magnetisation_tensor(beta) - m = local_expectation_value(ψ, M, envs) / denom - @test m ≈ m_th atol = 1.0e-10 + m = expectation_value(ψ, (O_mpo, 1 => M)) / denom + @test abs(m) ≈ m_th atol = 1.0e-10 # account for spin flip E = ising_energy_tensor(beta) - e = local_expectation_value(ψ, E, envs) / denom + e = expectation_value(ψ, (O_mpo, 1 => E)) / denom @test e ≈ e_th atol = 1.0e-2 end diff --git a/test/operators.jl b/test/operators.jl index a817c3c35..3ea1ae003 100644 --- a/test/operators.jl +++ b/test/operators.jl @@ -245,11 +245,8 @@ module TestOperators e1 = expectation_value(ψ2, H1) e3 = expectation_value(ψ2, H3) - e3_1 = local_expectation_value(ψ2, H3, environments(ψ2, H3), 1) - e3_2 = local_expectation_value(ψ2, H3, environments(ψ2, H3), 2) @test e1 + e3 ≈ expectation_value(ψ2, H1 + H3) atol = 1.0e-10 - @test e3 ≈ e3_1 + e3_2 atol = 1.0e-10 H4 = H1 + H3 h4 = H4 * H4 From 46abd10bb0ef43622bdeceda991aedad7bdaab84 Mon Sep 17 00:00:00 2001 From: Boris De Vos Date: Tue, 30 Sep 2025 08:43:49 +0200 Subject: [PATCH 10/15] docstring --- src/algorithms/expval.jl | 16 ++++++++++------ 1 file changed, 10 insertions(+), 6 deletions(-) diff --git a/src/algorithms/expval.jl b/src/algorithms/expval.jl index 60239aa98..baf32d13f 100644 --- a/src/algorithms/expval.jl +++ b/src/algorithms/expval.jl @@ -1,22 +1,26 @@ """ expectation_value(ψ, O, [environments]) expectation_value(ψ, inds => O) + expectation_value(ψ, (mpo, site => O), [environments]) Compute the expectation value of an operator `O` on a state `ψ`. Optionally, it is possible to make the computations more efficient by also passing in previously calculated `environments`. In general, the operator `O` may consist of an arbitrary MPO `O <: AbstractMPO` that acts -on all sites, or a local operator `O = inds => operator` acting on a subset of sites. -In the latter case, `inds` is a tuple of indices that specify the sites on which the operator -acts, while the operator is either a `AbstractTensorMap` or a `FiniteMPO`. +on all sites, a local operator `O = inds => operator` acting on a subset of sites, or a local MPO tensor +acting on a site within a network whose environment is determined by another MPO `mpo`. +In the second case, `inds` is a tuple of indices that specify the sites on which the operator +acts, while the operator is either a `AbstractTensorMap` or a `FiniteMPO`. In the latter case, +the operator is a `AbstractTensorMap` that acts on the physical space of a single site. # Arguments * `ψ::AbstractMPS` : the state on which to compute the expectation value -* `O::Union{AbstractMPO,Pair}` : the operator to compute the expectation value of. - This can either be an `AbstractMPO`, or a pair of indices and local operator.. +* `O::Union{AbstractMPO,Pair,AbstractTensorMap}` : the operator to compute the expectation value of. + This can either be an `AbstractMPO`, a pair of indices and local operator, or a local MPO tensor + represented as a `AbstractTensorMap`. * `environments::AbstractMPSEnvironments` : the environments to use for the calculation. If not given, they will be calculated. - + Depending on the type of `O`, these will be the environments of the operator `O` or the MPO `mpo`. # Examples ```jldoctest julia> ψ = FiniteMPS(ones(Float64, (ℂ^2)^4)); From 52d0d03e8fbd968c9aa3b7ab4f713c3ffc26e780 Mon Sep 17 00:00:00 2001 From: Boris De Vos Date: Tue, 30 Sep 2025 09:03:11 +0200 Subject: [PATCH 11/15] fix test --- test/algorithms.jl | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/test/algorithms.jl b/test/algorithms.jl index 13e505435..8a78cf44d 100644 --- a/test/algorithms.jl +++ b/test/algorithms.jl @@ -478,7 +478,9 @@ module TestAlgorithms @test f ≈ f_th atol = 1.0e-10 O = ising_bulk_tensor(beta) - denom = expectation_value(ψ, O, envs) + denom = expectation_value(ψ, O_mpo, envs) + denom2 = expectation_value(mps, (O_mpo, 1 => O), envs) + @test denom ≈ denom2 atol = 1.0e-10 M = ising_magnetisation_tensor(beta) m = expectation_value(ψ, (O_mpo, 1 => M)) / denom From a2c8fcee65e3a8928d0192c2c79afc8c9ace6a30 Mon Sep 17 00:00:00 2001 From: Boris De Vos Date: Tue, 30 Sep 2025 09:30:14 +0200 Subject: [PATCH 12/15] typo --- test/algorithms.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/test/algorithms.jl b/test/algorithms.jl index 8a78cf44d..a14a2a113 100644 --- a/test/algorithms.jl +++ b/test/algorithms.jl @@ -479,7 +479,7 @@ module TestAlgorithms O = ising_bulk_tensor(beta) denom = expectation_value(ψ, O_mpo, envs) - denom2 = expectation_value(mps, (O_mpo, 1 => O), envs) + denom2 = expectation_value(ψ, (O_mpo, 1 => O), envs) @test denom ≈ denom2 atol = 1.0e-10 M = ising_magnetisation_tensor(beta) From 4744604645b74421a9e622538cfd4ac9ef6efb8b Mon Sep 17 00:00:00 2001 From: Boris De Vos Date: Tue, 30 Sep 2025 10:05:35 +0200 Subject: [PATCH 13/15] fix and relax test --- test/algorithms.jl | 2 +- test/setup.jl | 1 + 2 files changed, 2 insertions(+), 1 deletion(-) diff --git a/test/algorithms.jl b/test/algorithms.jl index a14a2a113..ad25d39bc 100644 --- a/test/algorithms.jl +++ b/test/algorithms.jl @@ -484,7 +484,7 @@ module TestAlgorithms M = ising_magnetisation_tensor(beta) m = expectation_value(ψ, (O_mpo, 1 => M)) / denom - @test abs(m) ≈ m_th atol = 1.0e-10 # account for spin flip + @test abs(m) ≈ m_th atol = 1.0e-8 # account for spin flip E = ising_energy_tensor(beta) e = expectation_value(ψ, (O_mpo, 1 => E)) / denom diff --git a/test/setup.jl b/test/setup.jl index 6e7e3b487..f7c232b33 100644 --- a/test/setup.jl +++ b/test/setup.jl @@ -324,6 +324,7 @@ end function ising_energy_tensor(β) nt = ising_bond_tensor(β) O = ising_bulk_tensor(β).data + O = reshape(O, 2, 2, 2, 2) J = 1.0 e = ComplexF64[-J J; J -J] .* nt @tensor e_hor[-1 -2; -3 -4] := From a165a9f33230ea1a7899aad8df0f7e5d6c26bed7 Mon Sep 17 00:00:00 2001 From: Boris De Vos Date: Tue, 30 Sep 2025 10:33:26 +0200 Subject: [PATCH 14/15] simplify to stop my confusion --- test/algorithms.jl | 5 ++--- test/setup.jl | 49 ++++++++++++++++++++-------------------------- 2 files changed, 23 insertions(+), 31 deletions(-) diff --git a/test/algorithms.jl b/test/algorithms.jl index ad25d39bc..56c2d1d58 100644 --- a/test/algorithms.jl +++ b/test/algorithms.jl @@ -477,16 +477,15 @@ module TestAlgorithms f = -log(λ) / beta @test f ≈ f_th atol = 1.0e-10 - O = ising_bulk_tensor(beta) + O, M, E = classical_ising_tensors(beta) + denom = expectation_value(ψ, O_mpo, envs) denom2 = expectation_value(ψ, (O_mpo, 1 => O), envs) @test denom ≈ denom2 atol = 1.0e-10 - M = ising_magnetisation_tensor(beta) m = expectation_value(ψ, (O_mpo, 1 => M)) / denom @test abs(m) ≈ m_th atol = 1.0e-8 # account for spin flip - E = ising_energy_tensor(beta) e = expectation_value(ψ, (O_mpo, 1 => E)) / denom @test e ≈ e_th atol = 1.0e-2 end diff --git a/test/setup.jl b/test/setup.jl index f7c232b33..d4bddeeed 100644 --- a/test/setup.jl +++ b/test/setup.jl @@ -18,8 +18,7 @@ export force_planar export symm_mul_mpo export transverse_field_ising, heisenberg_XXX, bilinear_biquadratic_model, XY_model, kitaev_model -export classical_ising, sixvertex -export ising_bond_tensor, ising_bulk_tensor, ising_magnetisation_tensor, ising_energy_tensor +export classical_ising_tensors, classical_ising, sixvertex # using TensorOperations @@ -294,49 +293,43 @@ function kitaev_model(; t = 1.0, mu = 1.0, Delta = 1.0, L = Inf) end end -function ising_bond_tensor(β) - t = [exp(β) exp(-β); exp(-β) exp(β)] +function classical_ising_tensors(beta) + J = 1.0 + K = beta * J + + # Boltzmann weights + t = ComplexF64[exp(K) exp(-K); exp(-K) exp(K)] r = eigen(t) nt = r.vectors * sqrt(Diagonal(r.values)) * r.vectors - return nt -end -function ising_bulk_tensor(β) - nt = ising_bond_tensor(β) + # local partition function tensor δbulk = zeros(ComplexF64, (2, 2, 2, 2)) δbulk[1, 1, 1, 1] = 1 δbulk[2, 2, 2, 2] = 1 @tensor obulk[-1 -2; -3 -4] := δbulk[1 2; 3 4] * nt[-1; 1] * nt[-2; 2] * nt[-3; 3] * nt[-4; 4] - Obulk = TensorMap(obulk, ℂ^2 * ℂ^2, ℂ^2 * ℂ^2) - return Obulk -end -function ising_magnetisation_tensor(β) - nt = ising_bond_tensor(β) - M = zeros(2, 2, 2, 2) - M[1, 1, 1, 1] = 1 - M[2, 2, 2, 2] = -1 + # magnetization tensor + M = copy(δbulk) + M[2, 2, 2, 2] *= -1 @tensor m[-1 -2; -3 -4] := M[1 2; 3 4] * nt[-1; 1] * nt[-2; 2] * nt[-3; 3] * nt[-4; 4] - return TensorMap(m, ℂ^2 * ℂ^2 ← ℂ^2 * ℂ^2) -end -function ising_energy_tensor(β) - nt = ising_bond_tensor(β) - O = ising_bulk_tensor(β).data - O = reshape(O, 2, 2, 2, 2) - J = 1.0 + # bond interaction tensor and energy-per-site tensor e = ComplexF64[-J J; J -J] .* nt @tensor e_hor[-1 -2; -3 -4] := - O[1 2; 3 4] * nt[-1; 1] * nt[-2; 2] * nt[-3; 3] * e[-4; 4] + δbulk[1 2; 3 4] * nt[-1; 1] * nt[-2; 2] * nt[-3; 3] * e[-4; 4] @tensor e_vert[-1 -2; -3 -4] := - O[1 2; 3 4] * nt[-1; 1] * nt[-2; 2] * e[-3; 3] * nt[-4; 4] + δbulk[1 2; 3 4] * nt[-1; 1] * nt[-2; 2] * e[-3; 3] * nt[-4; 4] e = e_hor + e_vert - return TensorMap(e, ℂ^2 * ℂ^2 ← ℂ^2 * ℂ^2) -end + + # fixed tensor map space for all three + TMS = ℂ^2 ⊗ ℂ^2 ← ℂ^2 ⊗ ℂ^2 + + return TensorMap(obulk, TMS), TensorMap(m, TMS), TensorMap(e, TMS) +end; function classical_ising(; β = log(1 + sqrt(2)) / 2, L = Inf) - Obulk = ising_bulk_tensor(β) + Obulk, _, _ = classical_ising_tensors(β) L == Inf && return InfiniteMPO([Obulk]) From 60af5a4d0746d9c5a97a2fb9b048fc6ec4f96bf9 Mon Sep 17 00:00:00 2001 From: Boris De Vos Date: Tue, 30 Sep 2025 12:28:13 +0200 Subject: [PATCH 15/15] removed too much --- test/algorithms.jl | 2 +- test/setup.jl | 11 +++++++++-- 2 files changed, 10 insertions(+), 3 deletions(-) diff --git a/test/algorithms.jl b/test/algorithms.jl index 56c2d1d58..64b75a47a 100644 --- a/test/algorithms.jl +++ b/test/algorithms.jl @@ -478,7 +478,7 @@ module TestAlgorithms @test f ≈ f_th atol = 1.0e-10 O, M, E = classical_ising_tensors(beta) - + denom = expectation_value(ψ, O_mpo, envs) denom2 = expectation_value(ψ, (O_mpo, 1 => O), envs) @test denom ≈ denom2 atol = 1.0e-10 diff --git a/test/setup.jl b/test/setup.jl index d4bddeeed..8c5cf92e2 100644 --- a/test/setup.jl +++ b/test/setup.jl @@ -293,14 +293,20 @@ function kitaev_model(; t = 1.0, mu = 1.0, Delta = 1.0, L = Inf) end end -function classical_ising_tensors(beta) +function ising_bond_tensor(β) J = 1.0 - K = beta * J + K = β * J # Boltzmann weights t = ComplexF64[exp(K) exp(-K); exp(-K) exp(K)] r = eigen(t) nt = r.vectors * sqrt(Diagonal(r.values)) * r.vectors + return nt +end + +function classical_ising_tensors(beta) + nt = ising_bond_tensor(beta) + J = 1.0 # local partition function tensor δbulk = zeros(ComplexF64, (2, 2, 2, 2)) @@ -333,6 +339,7 @@ function classical_ising(; β = log(1 + sqrt(2)) / 2, L = Inf) L == Inf && return InfiniteMPO([Obulk]) + nt = ising_bond_tensor(β) δleft = zeros(ComplexF64, (1, 2, 2, 2)) δleft[1, 1, 1, 1] = 1 δleft[1, 2, 2, 2] = 1