# Bose Hubbard model with cavity interaction **URL:** \ **Category:** ITensor Julia Questions\ **Created:** [September 16, 2025, 6:16am UTC](https://itensor.discourse.group/t/bose-hubbard-model-with-cavity-interaction/2494 "2025-09-16T06:16:08Z")\ **Posts on this page:** 8\ **Page:** 1
**Author:** ![SOURITRA](https://avatars.discourse-cdn.com/v4/letter/s/ccd318/32.png) [@SOURITRA](https://itensor.discourse.group/u/SOURITRA)\ **Post date:** [September 16, 2025, 6:16am UTC](https://itensor.discourse.group/t/bose-hubbard-model-with-cavity-interaction/2494/1 "2025-09-16T06:16:08Z")
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, ![Screenshot from 2025-09-16 11-38-34](https://global.discourse-cdn.com/free1/uploads/itensor/original/1X/727d01a4796845581564ecf8b6c8986005cdd5d1.png) where, ![Screenshot from 2025-09-16 11-39-59](https://global.discourse-cdn.com/free1/uploads/itensor/original/1X/37c043dbf3a5194c354551846b13cdd065562aab.png) ```auto ``` ```auto 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 ``` ```auto 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 ```auto ``` ---
**Author:** ![mathew0036](https://yyz2.discourse-cdn.com/free1/user_avatar/itensor.discourse.group/mathew0036/32/214_2.png) [@mathew0036](https://itensor.discourse.group/u/mathew0036)\ **Post date:** [September 16, 2025, 7:50am UTC](https://itensor.discourse.group/t/bose-hubbard-model-with-cavity-interaction/2494/2 "2025-09-16T07:50:40Z")
Two problems I can point out: 1. > [@SOURITRA](#): > > `Sup = AutoMPO()` This should be `OpSum()` everywhere. AutoMPO() is deprecated code. 1. The error is probably coming from here > [@SOURITRA](#): > > `ampo += epsilon, 0.5*N, "Id", i` You can change this to `ampo += epsilon*0.5*N, "Id", i` I hope this helps! ---
**Author:** ![SOURITRA](https://avatars.discourse-cdn.com/v4/letter/s/ccd318/32.png) [@SOURITRA](https://itensor.discourse.group/u/SOURITRA)\ **Post date:** [September 16, 2025, 9:51am UTC](https://itensor.discourse.group/t/bose-hubbard-model-with-cavity-interaction/2494/3 "2025-09-16T09:51:05Z")
Still facing the same issue . The MPO operator showing the the error. ```auto 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 ``` ---
**Author:** ![mtfishman](https://yyz2.discourse-cdn.com/free1/user_avatar/itensor.discourse.group/mtfishman/32/11_2.png) [@mtfishman](https://itensor.discourse.group/u/mtfishman)\ **Post date:** [September 16, 2025, 4:11pm UTC](https://itensor.discourse.group/t/bose-hubbard-model-with-cavity-interaction/2494/4 "2025-09-16T16:11:15Z")
I see you are constructing the MPO like this: `H = MPO(ampo, s)`, however you define `s` as: `s = Index(N, "polaritons")` which is _single_ Index with a dimension `N`. Maybe you want `s = [Index(10, “polaritons”) for _ in 1:N]` to make a chain of `N` polaritons (each of dimension 10). ---
**Author:** ![SOURITRA](https://avatars.discourse-cdn.com/v4/letter/s/ccd318/32.png) [@SOURITRA](https://itensor.discourse.group/u/SOURITRA)\ **Post date:** [September 17, 2025, 4:46am UTC](https://itensor.discourse.group/t/bose-hubbard-model-with-cavity-interaction/2494/5 "2025-09-17T04:46:05Z")
Thank you for replying but It still showing error. ```auto ERROR: LoadError: ArgumentError: Overload of "op" or "op!" functions not found for operator name "N" and Index tags: ("polaritons",). Stacktrace: [1] op(name::String, s::Index{Int64}; adjoint::Bool, kwargs::@Kwargs{}) @ ITensors.SiteTypes ~/.julia/packages/ITensors/AbUYR/src/lib/SiteTypes/src/sitetype.jl:402 [2] op @ ~/.julia/packages/ITensors/AbUYR/src/lib/SiteTypes/src/sitetype.jl:236 [inlined] [3] op(s::Index{Int64}, opname::String) @ ITensors.SiteTypes ~/.julia/packages/ITensors/AbUYR/src/lib/SiteTypes/src/sitetype.jl:449 [4] computeSiteProd(sites::Vector{Index{Int64}}, ops::Prod{Op}) @ ITensors.ITensorMPS ~/.julia/packages/ITensors/AbUYR/src/lib/ITensorMPS/src/opsum_to_mpo/opsum_to_mpo_generic.jl:136 [5] svdMPO(ValType::Type{Float64}, os::Sum{Scaled{ComplexF64, Prod{Op}}}, sites::Vector{Index{Int64}}; mindim::Int64, maxdim::Int64, cutoff::Float64) @ ITensors.ITensorMPS ~/.julia/packages/ITensors/AbUYR/src/lib/ITensorMPS/src/opsum_to_mpo/opsum_to_mpo.jl:118 [6] svdMPO @ ~/.julia/packages/ITensors/AbUYR/src/lib/ITensorMPS/src/opsum_to_mpo/opsum_to_mpo.jl:5 [inlined] [7] #svdMPO#541 @ ~/.julia/packages/ITensors/AbUYR/src/lib/ITensorMPS/src/opsum_to_mpo/opsum_to_mpo.jl:147 [inlined] [8] svdMPO @ ~/.julia/packages/ITensors/AbUYR/src/lib/ITensorMPS/src/opsum_to_mpo/opsum_to_mpo.jl:144 [inlined] [9] MPO(os::Sum{Scaled{ComplexF64, Prod{Op}}}, sites::Vector{Index{Int64}}; splitblocks::Bool, kwargs::@Kwargs{}) @ ITensors.ITensorMPS ~/.julia/packages/ITensors/AbUYR/src/lib/ITensorMPS/src/opsum_to_mpo/opsum_to_mpo_generic.jl:299 [10] MPO(os::Sum{Scaled{ComplexF64, Prod{Op}}}, sites::Vector{Index{Int64}}) @ ITensors.ITensorMPS ~/.julia/packages/ITensors/AbUYR/src/lib/ITensorMPS/src/opsum_to_mpo/opsum_to_mpo_generic.jl:289 [11] top-level scope ``` ---
**Author:** ![SOURITRA](https://avatars.discourse-cdn.com/v4/letter/s/ccd318/32.png) [@SOURITRA](https://itensor.discourse.group/u/SOURITRA)\ **Post date:** [September 17, 2025, 7:24am UTC](https://itensor.discourse.group/t/bose-hubbard-model-with-cavity-interaction/2494/6 "2025-09-17T07:24:59Z")
I am still confused about how to put the Hamiltonian in MPO form what I posted earlier. It would be really helpful for me if anyone from this discussion group could help. ---
**Author:** ![mtfishman](https://yyz2.discourse-cdn.com/free1/user_avatar/itensor.discourse.group/mtfishman/32/11_2.png) [@mtfishman](https://itensor.discourse.group/u/mtfishman)\ **Post date:** [September 17, 2025, 12:39pm UTC](https://itensor.discourse.group/t/bose-hubbard-model-with-cavity-interaction/2494/7 "2025-09-17T12:39:45Z")
I see that you’re using an operator `"N"` when you are constructing your `OpSum` but you don’t define an operator `"N"` for your `"polaritons"` SiteType, maybe you’re missing that definition? ---
**Author:** ![SOURITRA](https://avatars.discourse-cdn.com/v4/letter/s/ccd318/32.png) [@SOURITRA](https://itensor.discourse.group/u/SOURITRA)\ **Post date:** [October 20, 2025, 9:53am UTC](https://itensor.discourse.group/t/bose-hubbard-model-with-cavity-interaction/2494/8 "2025-10-20T09:53:10Z")
Right now my code is working but the local density plot not showing good results. It shows that the system reaches superfluid state much earlier. I want to reproduce the PRL 99, 186401 (2007) results. Here I have shared the code please check where I am making mistakes. ```auto 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.5*[0 1; 1 0] sigma_y = 0.5*[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(I_atom, a) ITensors.op(::OpName"Adag",::SiteType"polaritons") = kron(I_atom, a_t) ITensors.op(::OpName"Num",::SiteType"polaritons") = kron(I_atom, n) + kron(sigma_up*sigma_dn, I_boson) ITensors.op(::OpName"Sz",::SiteType"polaritons") = kron(sz, I_boson) ITensors.op(::OpName"sig_up",::SiteType"polaritons") = kron(sigma_up, I_boson) ITensors.op(::OpName"sig_dn",::SiteType"polaritons") = kron(sigma_dn, I_boson) ITensors.op(::OpName"Id",::SiteType"polaritons") = kron(I_atom, I_boson) s = siteind("polaritons") A = op("A", s) Adag = op("Adag", s) Num = op("Num", s) Sz = op("Sz", s) sig_up = op("sig_up", s) sig_dn = op("sig_dn", s) Id = op("Id", s) let L = 5 #t = 0.24 mu = 0.5 omega = 0.5 beta = 0.1 epsilon = 0.5 Energy = Float64[] values = Float64[] s = siteinds("polaritons", L) for t in 0.0:0.012:0.24 println("t = ",t) state = fill(2, L) state[L] = 1 psi0 = productMPS(s,state) ampo = OpSum() for i in 1:L-1 ampo += -t, "Adag", i, "A", i+1 ampo += -t, "A", i, "Adag", i+1 end for i in 1:L ampo += -mu, "Num", i end for i in 1:L ampo += beta, "sig_up", i, "A", i ampo += beta, "sig_dn", i, "Adag", i end for i in 1:L ampo += epsilon, "Sz", i ampo += epsilon*0.5, "Id", i end for i in 1:L ampo += omega, "Adag", i, "A", i end H = MPO(ampo, s) nsweeps = 100 maxdim = [10,20,40,60,100,200] cutoff = [10E-10] energy,psi = dmrg(H, psi0; nsweeps, maxdim, cutoff) println(exp_val) push!(values, exp_val[2]) end println(values) return end ``` Please help me as early as possible. Thank you