From 705fd9e0e8edfc24169e1b6cc2020911082b9393 Mon Sep 17 00:00:00 2001 From: lkdvos Date: Wed, 11 Sep 2024 11:28:34 +0200 Subject: [PATCH] Add potts model --- src/MPSKitModels.jl | 1 + src/models/hamiltonians.jl | 39 +++++++++++++++++++++- src/operators/spinoperators.jl | 61 ++++++++++++++++++++++++++++++++++ test/potts.jl | 21 ++++++++++++ 4 files changed, 121 insertions(+), 1 deletion(-) create mode 100644 test/potts.jl diff --git a/src/MPSKitModels.jl b/src/MPSKitModels.jl index 2426ab8..32bc1bf 100644 --- a/src/MPSKitModels.jl +++ b/src/MPSKitModels.jl @@ -37,6 +37,7 @@ export e⁺, e⁻, e⁺⁺, e⁻⁻, e⁺⁻, e⁻⁺ export transverse_field_ising export kitaev_model +export quantum_potts export heisenberg_XXX, heisenberg_XXZ, heisenberg_XYZ export bilinear_biquadratic_model export hubbard_model, bose_hubbard_model diff --git a/src/models/hamiltonians.jl b/src/models/hamiltonians.jl index d7d5152..5b2a37d 100644 --- a/src/models/hamiltonians.jl +++ b/src/models/hamiltonians.jl @@ -1,7 +1,6 @@ #=========================================================================================== Ising model ===========================================================================================# -ising_kwargs = (; J=1.0, g=1.0) """ transverse_field_ising([elt::Type{<:Number}], [symmetry::Type{<:Sector}], @@ -234,6 +233,44 @@ function bilinear_biquadratic_model(elt::Type{<:Number}=ComplexF64, end end +#=========================================================================================== + Potts models +===========================================================================================# +""" + quantum_potts([elt::Type{<:Number}], [symmetry::Type{<:Sector}], + [lattice::AbstractLattice]; q=3, J=1.0, g=1.0) + +MPO for the hamiltonian of the quantum Potts model, as defined by +```math +H = - J \\sum_{\\langle i,j \\rangle} Z_i^\\dagger Z_j + Z_i Z_j^\\dagger +- g \\sum_i (X_i + X_i^\\dagger) +``` +where the operators ``Z`` and ``X`` are the ``q``-rotation operators satisfying +``Z^q = X^q = 1`` and ``ZX = \\omega XZ`` where ``\\omega = e^{2πi/q}``. +""" +function quantum_potts end +function quantum_potts(lattice::AbstractLattice; kwargs...) + return quantum_potts(ComplexF64, Trivial, lattice; kwargs...) +end +function quantum_potts(symmetry::Type{<:Sector}, + lattice::AbstractLattice=InfiniteChain(1); kwargs...) + return quantum_potts(ComplexF64, symmetry, lattice; kwargs...) +end +function quantum_potts(elt::Type{<:Number}, lattice::AbstractLattice; + kwargs...) + return quantum_potts(elt, Trivial, lattice; kwargs...) +end +function quantum_potts(elt::Type{<:Number}=ComplexF64, + symmetry::Type{<:Sector}=Trivial, + lattice::AbstractLattice=InfiniteChain(1); + q=3, J=1.0, g=1.0) + return @mpoham sum(nearest_neighbours(lattice)) do (i, j) + return -J * potts_exchange(elt, symmetry; q){i,j} + end - sum(vertices(lattice)) do i + return g * potts_field(elt, symmetry; q){i} + end +end + #=========================================================================================== Hubbard models ===========================================================================================# diff --git a/src/operators/spinoperators.jl b/src/operators/spinoperators.jl index 30a3c68..1d2743d 100644 --- a/src/operators/spinoperators.jl +++ b/src/operators/spinoperators.jl @@ -444,3 +444,64 @@ const SS = S_exchange """Pauli exchange operator.""" σσ(args...; kwargs...) = 4 * S_exchange(args...; kwargs...) + +""" + potts_exchange([eltype::Type{<:Number}], [symmetry::Type{<:Sector}]; q=3) + +The Potts exchange operator ``Z ⊗ Z' + Z' ⊗ Z``, where ``Z^q = 1``. +""" +function potts_exchange end +potts_exchange(; kwargs...) = potts_exchange(ComplexF64, Trivial; kwargs...) +potts_exchange(elt::Type{<:Number}; kwargs...) = potts_exchange(elt, Trivial; kwargs...) +function potts_exchange(symmetry::Type{<:Sector}; kwargs...) + return potts_exchange(ComplexF64, symmetry; kwargs...) +end + +function potts_exchange(elt::Type{<:Number}, ::Type{Trivial}; q=3) + pspace = ComplexSpace(q) + Z = TensorMap(zeros, elt, pspace ← pspace) + for i in 1:q + Z[i, i] = cis(2π * (i - 1) / q) + end + return Z ⊗ Z' + Z' ⊗ Z +end +function potts_exchange(elt::Type{<:Number}, ::Type{ZNIrrep{Q}}; q=Q) where {Q} + @assert q == Q "q must match the irrep charge" + pspace = Vect[ZNIrrep{q}](i => 1 for i in 0:(q - 1)) + aspace = Vect[ZNIrrep{q}](1 => 1, -1 => 1) + Z_left = TensorMap(ones, elt, pspace ← pspace ⊗ aspace) + Z_right = TensorMap(ones, elt, aspace ⊗ pspace ← pspace) + return contract_twosite(Z_left, Z_right) +end + +""" + potts_field([eltype::Type{<:Number}], [symmetry::Type{<:Sector}]; q=3) + +The Potts field operator ``X + X'``, where ``X^q = 1``. +""" +function potts_field end +potts_field(; kwargs...) = potts_field(ComplexF64, Trivial; kwargs...) +potts_field(elt::Type{<:Number}; kwargs...) = potts_field(elt, Trivial; kwargs...) +function potts_field(symmetry::Type{<:Sector}; kwargs...) + return potts_field(ComplexF64, symmetry; kwargs...) +end + +function potts_field(elt::Type{<:Number}, ::Type{Trivial}; q=3) + pspace = ComplexSpace(q) + X = TensorMap(zeros, elt, pspace ← pspace) + for i in 1:q + X[mod1(i - 1, q), i] = one(elt) + end + return X + X' +end +# TODO: generalize to arbitrary q +function potts_field(elt::Type{<:Number}, ::Type{ZNIrrep{Q}}; q=Q) where {Q} + @assert q == Q "q must match the irrep charge" + @assert q == 3 "only q = 3 is implemented" + pspace = Vect[ZNIrrep{q}](i => 1 for i in 0:(q - 1)) + X = TensorMap(zeros, elt, pspace ← pspace) + for (c, b) in blocks(X) + b .= isone(c) ? 2 : -1 + end + return X +end diff --git a/test/potts.jl b/test/potts.jl new file mode 100644 index 0000000..a54998e --- /dev/null +++ b/test/potts.jl @@ -0,0 +1,21 @@ +using MPSKit +using TensorKit + +alg = VUMPS(; maxiter=25, verbosity=0) +E₀ = -(4 / 3 + 2sqrt(3) / π) + +@testset "Trivial" begin + H = quantum_potts(; q=3) + ψ = InfiniteMPS(3, 32) + @test imag(expectation_value(ψ, H)) ≈ 0 atol = 1e-10 + ψ, envs, δ = find_groundstate(ψ, H, alg) + @test E₀ ≈ expectation_value(ψ, H, envs) atol = 1e-2 +end + +@testset "Z3Irrep" begin + H = quantum_potts(Z3Irrep; q=3) + ψ = InfiniteMPS(Rep[ℤ₃](i => 1 for i in 0:2), Rep[ℤ₃](i => 12 for i in 0:2)) + @test imag(expectation_value(ψ, H)) ≈ 0 atol = 1e-10 + ψ, envs, δ = find_groundstate(ψ, H, alg) + @test E₀ ≈ expectation_value(ψ, H, envs) atol = 1e-2 +end