Hello everyone I trying to simulate bose hubbard model with cavity interaction for 2 label atom. I want to see the phase diagram using DMRG, for that my Hamiltonian is,

where,

using ITensors, ITensorMPS, LinearAlgebra
ITensors.space(::SiteType"polaritons") = 10
n = [0 0 0 0 0;
0 1 0 0 0;
0 0 2 0 0;
0 0 0 3 0;
0 0 0 0 4]
a = [0 1 0 0 0;
0 0 sqrt(2) 0 0;
0 0 0 sqrt(3) 0;
0 0 0 0 sqrt(4);
0 0 0 0 0]
a_t = [0 0 0 0 0;
1 0 0 0 0;
0 sqrt(2) 0 0 0;
0 0 sqrt(3) 0 0;
0 0 0 sqrt(4) 0]
I_atom = Diagonal(ones(2))
I_boson = Diagonal(ones(5))
Sz = [0.5 0;
0 -0.5]
sigma_x = [0 1;
1 0]
sigma_y = [0 -1.0*im;
1.0*im 0]
sigma_up = sigma_x + 1.0*im*sigma_y
sigma_dn = sigma_x - 1.0*im*sigma_y
ITensors.op(::OpName"A",::SiteType"polaritons") =
kron(a, I_atom)
ITensors.op(::OpName"Adag",::SiteType"polaritons") =
kron(a_t, I_atom)
ITensors.op(::OpName"Num",::SiteType"polaritons") =
kron(n, I_atom)
ITensors.op(::OpName"Sz",::SiteType"polaritons") =
kron(I_boson, Sz)
ITensors.op(::OpName"sig_up",::SiteType"polaritons") =
kron(I_boson, sigma_up)
ITensors.op(::OpName"sig_dn",::SiteType"polaritons") =
kron(I_boson, sigma_dn)
ITensors.op(::OpName"Id",::SiteType"polaritons") =
kron(I_boson, I_atom)
N = 10
s = Index(N, "polaritons")
a = op("A", s)
adag = op("Adag", s)
N = op("Num", s)
Sz = op("Sz", s)
sig_up = op("sig_up", s)
sig_dn = op("sig_dn", s)
Id = op("Id", s)
let
N = 10
t = 0.1
mu = 0.2
omega = 0.3
beta = 0.1
epsilon = 0.2
s = Index(N, "polaritons")
ampo = AutoMPO()
Sup = AutoMPO()
Sdn = AutoMPO()
S_z = AutoMPO()
for i in 1:N
Sup += (1.0, "sig_up", i)
end
for i in 1:N
Sdn += (1.0, "sig_dn", i)
end
for i in 1:N
Sdn += (1.0, "Sz", i)
end
for i in 1:N-1
ampo += -t, "adag", i, "a", i+1
ampo += -t, "a", i, "adag", i
end
for i in 1:N
ampo += -mu, "N", i
end
for i in 1:N
ampo += beta, "S_up", i, "a", i
ampo += beta, "S_dn", i, "adag", i
end
for i in 1:N
ampo += omega, "adag", i, "a", i
ampo += epsilon, "S_z", i
ampo += epsilon, 0.5*N, "Id", i
end
H = MPO(ampo, s)
return
end
ERROR: LoadError: MethodError: no method matching MPO(::Sum{Scaled{ComplexF64, Prod{Op}}}, ::Index{Int32})
The type `MPO` exists, but no method is defined for this combination of argument types when trying to construct it.
Closest candidates are:
MPO(::Any, ::Any, ::Any)
@ ITensors ~/.julia/packages/ITensors/AbUYR/src/lib/ITensorMPS/src/mpo.jl:13
MPO(::Vector{<:Index}, ::Any)
@ ITensors ~/.julia/packages/ITensors/AbUYR/src/lib/ITensorMPS/src/mpo.jl:97
MPO(::ITensors.LazyApply.Applied{typeof(sum), Tuple{Array{ITensors.LazyApply.Applied{typeof(*), Tuple{C, Prod{Op}}, @NamedTuple{}}, 1}}, @NamedTuple{}} where C, ::Vector{<:Index}; splitblocks, kwargs...)
@ ITensors ~/.julia/packages/ITensors/AbUYR/src/lib/ITensorMPS/src/opsum_to_mpo/opsum_to_mpo_generic.jl:289
...
Stacktrace:
[1] top-level scope
Please help