The Ising CFT spectrum

This tutorial is meant to show the finite size CFT spectrum for the quantum Ising model. We do this by first employing an exact diagonalization technique, and then extending the analysis to larger system sizes through the use of MPS techniques.

using MPSKit, MPSKitModels, TensorKit, Plots, KrylovKit
using LinearAlgebra: eigen, diagm, Hermitian

The hamiltonian is defined on a finite lattice with periodic boundary conditions, which can be implemented as follows:

L = 12
H = periodic_boundary_conditions(transverse_field_ising(), L);

Exact diagonalisation

In MPSKit, there is support for exact diagonalisation by leveraging the fact that applying the hamiltonian to an untruncated MPS will result in an effective hamiltonian on the center site which implements the action of the entire hamiltonian. Thus, optimizing the middle tensor is equivalent to optimixing a state in the entire Hilbert space, as all other tensors are just unitary matrices that mix the basis.

energies, states = exact_diagonalization(H; num=18, alg=Lanczos(; krylovdim=200));
plot(real.(energies);
     seriestype=:scatter,
     legend=false,
     ylabel="energy",
     xlabel="#eigenvalue")

Extracting momentum

Given a state, it is possible to assign a momentum label through the use of the translation operator. This operator can be defined in MPO language either diagramatically as

or in the code as:

id = complex(isomorphism(ℂ^2, ℂ^2))
@tensor O[-1 -2; -3 -4] := id[-1, -3] * id[-2, -4]
T = periodic_boundary_conditions(DenseMPO(O), L);

We can then calculate the momentum of the groundstate as the expectation value of this operator. However, there is a subtlety because of the degeneracies in the energy eigenvalues. The eigensolver will find an orthonormal basis within each energy subspace, but this basis is not necessarily a basis of eigenstates of the translation operator. In order to fix this, we diagonalize the translation operator within each energy subspace.

momentum(ψᵢ, ψⱼ=ψᵢ) = angle(dot(ψᵢ, T * ψⱼ))

function fix_degeneracies(basis)
    N = zeros(ComplexF64, length(basis), length(basis))
    M = zeros(ComplexF64, length(basis), length(basis))
    for i in eachindex(basis), j in eachindex(basis)
        N[i, j] = dot(basis[i], basis[j])
        M[i, j] = momentum(basis[i], basis[j])
    end

    vals, vecs = eigen(Hermitian(N))
    M = (vecs' * M * vecs)
    M /= diagm(vals)

    vals, vecs = eigen(M)
    return angle.(vals)
end

momenta = Float64[]
append!(momenta, fix_degeneracies(states[1:1]))
append!(momenta, fix_degeneracies(states[2:2]))
append!(momenta, fix_degeneracies(states[3:3]))
append!(momenta, fix_degeneracies(states[4:5]))
append!(momenta, fix_degeneracies(states[6:9]))
append!(momenta, fix_degeneracies(states[10:11]))
append!(momenta, fix_degeneracies(states[12:12]))
append!(momenta, fix_degeneracies(states[13:16]))
append!(momenta, fix_degeneracies(states[17:18]))

plot(momenta,
     real.(energies[1:18]);
     seriestype=:scatter,
     xlabel="momentum",
     ylabel="energy",
     legend=false)

Finite bond dimension

If we limit the maximum bond dimension of the MPS, we get an approximate solution, but we can reach higher system sizes.

L_mps = 20
H_mps = periodic_boundary_conditions(transverse_field_ising(), L_mps)
D = 64
ψ, envs, δ = find_groundstate(FiniteMPS(L_mps, ℂ^2, ℂ^D), H_mps, DMRG());

[ Info: DMRG init:	obj = -1.986175369583e+01	err = 7.3962e-02
[ Info: DMRG   1:	obj = -2.549098908186e+01	err = 9.7854336373e-03	time = 2.85 sec
[ Info: DMRG   2:	obj = -2.549098968632e+01	err = 1.1440955438e-06	time = 0.65 sec
[ Info: DMRG   3:	obj = -2.549098968636e+01	err = 1.1899350273e-07	time = 0.61 sec
[ Info: DMRG   4:	obj = -2.549098968636e+01	err = 1.3300005913e-08	time = 0.34 sec
[ Info: DMRG   5:	obj = -2.549098968636e+01	err = 5.5087464909e-09	time = 0.31 sec
[ Info: DMRG   6:	obj = -2.549098968636e+01	err = 3.2413738015e-09	time = 0.30 sec
[ Info: DMRG   7:	obj = -2.549098968636e+01	err = 2.4243792973e-09	time = 0.31 sec
[ Info: DMRG   8:	obj = -2.549098968636e+01	err = 2.1231366707e-09	time = 0.56 sec
[ Info: DMRG   9:	obj = -2.549098968636e+01	err = 1.9891538232e-09	time = 0.44 sec
[ Info: DMRG  10:	obj = -2.549098968636e+01	err = 1.8864377638e-09	time = 0.32 sec
[ Info: DMRG  11:	obj = -2.549098968636e+01	err = 1.7682141746e-09	time = 0.33 sec
[ Info: DMRG  12:	obj = -2.549098968636e+01	err = 1.6262519428e-09	time = 0.30 sec
[ Info: DMRG  13:	obj = -2.549098968636e+01	err = 1.4678381585e-09	time = 0.30 sec
[ Info: DMRG  14:	obj = -2.549098968636e+01	err = 1.3037732764e-09	time = 0.30 sec
[ Info: DMRG  15:	obj = -2.549098968636e+01	err = 1.1436368174e-09	time = 0.30 sec
[ Info: DMRG  16:	obj = -2.549098968636e+01	err = 9.9399752806e-10	time = 0.30 sec
[ Info: DMRG  17:	obj = -2.549098968636e+01	err = 8.5842868469e-10	time = 0.29 sec
[ Info: DMRG  18:	obj = -2.549098968636e+01	err = 7.3821724717e-10	time = 0.29 sec
[ Info: DMRG  19:	obj = -2.549098968636e+01	err = 6.3317095217e-10	time = 0.29 sec
[ Info: DMRG  20:	obj = -2.549098968636e+01	err = 5.4224930858e-10	time = 0.47 sec
[ Info: DMRG  21:	obj = -2.549098968636e+01	err = 4.6407320953e-10	time = 0.30 sec
[ Info: DMRG  22:	obj = -2.549098968636e+01	err = 3.9709950959e-10	time = 0.28 sec
[ Info: DMRG  23:	obj = -2.549098968636e+01	err = 3.3985399711e-10	time = 0.29 sec
[ Info: DMRG  24:	obj = -2.549098968636e+01	err = 2.9097790850e-10	time = 0.28 sec
[ Info: DMRG  25:	obj = -2.549098968636e+01	err = 2.4926442816e-10	time = 0.28 sec
[ Info: DMRG  26:	obj = -2.549098968636e+01	err = 2.1365944710e-10	time = 0.28 sec
[ Info: DMRG  27:	obj = -2.549098968636e+01	err = 1.8325685774e-10	time = 0.29 sec
[ Info: DMRG  28:	obj = -2.549098968636e+01	err = 1.5728170079e-10	time = 0.27 sec
[ Info: DMRG  29:	obj = -2.549098968636e+01	err = 1.3507433129e-10	time = 0.27 sec
[ Info: DMRG  30:	obj = -2.549098968636e+01	err = 1.1601025066e-10	time = 0.32 sec
[ Info: DMRG conv 31:	obj = -2.549098968636e+01	err = 9.9471043716e-11	time = 13.06 sec

Excitations on top of the groundstate can be found through the use of the quasiparticle ansatz. This returns quasiparticle states, which can be converted to regular FiniteMPS objects.

E_ex, qps = excitations(H, QuasiparticleAnsatz(), ψ, envs; num=16)
states_mps = vcat(ψ, map(qp -> convert(FiniteMPS, qp), qps))
E_mps = map(x -> expectation_value(x, H_mps), states_mps)

T_mps = periodic_boundary_conditions(DenseMPO(O), L_mps)
momenta_mps = Float64[]
append!(momenta_mps, fix_degeneracies(states[1:1]))
append!(momenta_mps, fix_degeneracies(states[2:2]))
append!(momenta_mps, fix_degeneracies(states[3:3]))
append!(momenta_mps, fix_degeneracies(states[4:5]))
append!(momenta_mps, fix_degeneracies(states[6:9]))
append!(momenta_mps, fix_degeneracies(states[10:11]))
append!(momenta_mps, fix_degeneracies(states[12:12]))
append!(momenta_mps, fix_degeneracies(states[13:16]))

plot(momenta_mps,
     real.(energies[1:16]);
     seriestype=:scatter,
     xlabel="momentum",
     ylabel="energy",
     legend=false)

This page was generated using Literate.jl.