# Multiorbital electron mixed Hamiltonian

**URL:** <https://itensor.discourse.group/t/multiorbital-electron-mixed-hamiltonian/1444>\
**Category:** DMRG and Numerical Methods\
**Tags:** julia\
**Created:** [January 24, 2024, 10:05pm UTC](https://itensor.discourse.group/t/multiorbital-electron-mixed-hamiltonian/1444 "2024-01-24T22:05:35Z")\
**Posts on this page:** 4\
**Page:** 1

<div class="post-metadata">

**Author:** ![fegarcia](https://avatars.discourse-cdn.com/v4/letter/f/977dab/32.png) [@fegarcia](https://itensor.discourse.group/u/fegarcia)\
**Post date:** [January 24, 2024, 10:05pm UTC](https://itensor.discourse.group/t/multiorbital-electron-mixed-hamiltonian/1444/1 "2024-01-24T22:05:35Z")

</div>

Hi, thanks for the development of this library. I’ve been trying to create a set up similar to the one described in this paper ([Phys. Rev. Lett. 126, 197202 (2021) - Quantum Spin Torque Driven Transmutation of an Antiferromagnetic Mott Insulator](https://journals.aps.org/prl/abstract/10.1103/PhysRevLett.126.197202)) were, in a 1D chain, all sites are of the ‘Electron’ type, except for a central region where on top of the electron sites, the authors add another orbital for electrons to interact.  
I’ve been trying to use a mixed site setup as described in other examples. Thus, I’ve created a custom sites for the 2 electron orbital central region. I also want to conserve the total number of fermions. When I try to run DMRG on the constructed Hamiltonian I get the following error: Block Block(1, 1, 4) has flux QN((“Nf”,1,-1),(“Nt”,-1,-1)) that is inconsistent with the desired flux QN(“Nf”,0,-1).  
The error appears when I consider the hopping elements between the central region and the two other single electron regions.  
Here is a minimal code:

First I define the new site type for the central region including the symmetry of total number of fermions.

```julia
using ITensors

ITensors.enable_debug_checks()

function ITensors.space(::SiteType"2E";
                        conserve_qns=false, conserve_sz=false)
    if conserve_qns
        return [QN("Nt", 0, -1) => 1, QN("Nt", 1, -1) => 4,
                QN("Nt", 2, -1) => 6, QN("Nt", 3, -1) => 4,
                QN("Nt", 4, -1) => 1]
    end
    return 16
end

```

Then I label the 16 states

```julia
ITensors.val(::ValName"AEmpBEmp", ::SiteType"2E") = 1
ITensors.val(::ValName"AupBEmp", ::SiteType"2E") = 2
ITensors.val(::ValName"AdnBEmp", ::SiteType"2E") = 3
ITensors.val(::ValName"AEmpBup", ::SiteType"2E") = 4
ITensors.val(::ValName"AEmpBdn", ::SiteType"2E") = 5
ITensors.val(::ValName"AupdnBEmp", ::SiteType"2E") = 6
ITensors.val(::ValName"AEmpBupdn", ::SiteType"2E") = 7
ITensors.val(::ValName"AupBup", ::SiteType"2E") = 8
ITensors.val(::ValName"AdnBdn", ::SiteType"2E") = 9
ITensors.val(::ValName"AupBdn", ::SiteType"2E") = 10
ITensors.val(::ValName"AdnBup", ::SiteType"2E") = 11
ITensors.val(::ValName"AupdnBup", ::SiteType"2E") = 12
ITensors.val(::ValName"AupdnBdn", ::SiteType"2E") = 13
ITensors.val(::ValName"AupBupdn", ::SiteType"2E") = 14
ITensors.val(::ValName"AdnBupdn", ::SiteType"2E") = 15
ITensors.val(::ValName"AupdnBupdn", ::SiteType"2E") = 16

ITensors.state(::StateName"AEmpBEmp", ::SiteType"2E") = [1.0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]
ITensors.state(::StateName"AupBEmp", ::SiteType"2E") = [0, 1.0, 0, 0, 0, 0, 0, 0, 0 ,0, 0, 0, 0, 0, 0, 0]
ITensors.state(::StateName"AdnBEmp", ::SiteType"2E") = [0, 0, 1.0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]
ITensors.state(::StateName"AEmpBup", ::SiteType"2E") = [0, 0, 0, 1.0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]
ITensors.state(::StateName"AEmpBdn", ::SiteType"2E") = [0, 0, 0, 0, 1.0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]
ITensors.state(::StateName"AupdnBEmp", ::SiteType"2E") = [0, 0, 0, 0, 0, 1.0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]
ITensors.state(::StateName"AEmpBupdn", ::SiteType"2E") = [0, 0, 0, 0, 0, 0, 1.0, 0, 0, 0 ,0, 0, 0, 0, 0, 0]
ITensors.state(::StateName"AupBup", ::SiteType"2E") = [0, 0, 0, 0, 0, 0, 0, 1.0, 0, 0, 0, 0, 0, 0, 0, 0]
ITensors.state(::StateName"AupBdn", ::SiteType"2E") = [0, 0, 0, 0, 0, 0, 0, 0, 0, 1.0, 0, 0, 0, 0, 0, 0]
ITensors.state(::StateName"AdnBup", ::SiteType"2E") = [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1.0, 0, 0, 0, 0, 0]
ITensors.state(::StateName"AupdnBup", ::SiteType"2E") = [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1.0, 0, 0, 0, 0]
ITensors.state(::StateName"AupdnBdn", ::SiteType"2E") = [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1.0, 0, 0, 0]
ITensors.state(::StateName"AupBupdn", ::SiteType"2E") = [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1.0, 0, 0]
ITensors.state(::StateName"AdnBupdn", ::SiteType"2E") = [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1.0, 0]
ITensors.state(::StateName"AupdnBupdn", ::SiteType"2E") = [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1.0]

```

Then I create some basic operators. For the example I’m including only the hopping operators of the first orbital, which would be the orbital later connected to the single electron regions. I’m also adding the F fermonic operator

```julia
function ITensors.op!(Op::ITensor, ::OpName"Nt", ::SiteType"2E", s::Index)
    Op[s' => 2, s => 2] = 1
    Op[s' => 3, s => 3] = 1
    Op[s' => 4, s => 4] = 1
    Op[s' => 5, s => 5] = 1
    Op[s' => 6, s => 6] = 2
    Op[s' => 7, s => 7] = 2
    Op[s' => 8, s => 8] = 2
    Op[s' => 9, s => 9] = 2
    Op[s' => 10, s => 10] = 2
    Op[s' => 11, s => 11] = 2
    Op[s' => 12, s => 12] = 3
    Op[s' => 13, s => 13] = 3
    Op[s' => 14, s => 14] = 3
    Op[s' => 15, s => 15] = 3
    return Op[s' => 16, s => 16] = 4
end

function ITensors.op!(Op::ITensor, ::OpName"Caup", ::SiteType"2E", s::Index)
    Op[s' => 1, s => 2] = 1
    Op[s' => 3, s => 6] = -1
    Op[s' => 4, s => 8] = -1
    Op[s' => 5, s => 10] = -1
    Op[s' => 7, s => 14] = 1
    Op[s' => 9, s => 13] = 1
    Op[s' => 11, s => 12] = 1
    return Op[s' => 15, s => 16] = -1
end

function ITensors.op!(Op::ITensor, ::OpName"Cadagup", ::SiteType"2E", s::Index)
    Op[s' => 2, s => 1] = 1
    Op[s' => 6, s => 3] = -1
    Op[s' => 8, s => 4] = -1
    Op[s' => 10, s => 5] = -1
    Op[s' => 14, s => 7] = 1
    Op[s' => 13, s => 9] = 1
    Op[s' => 12, s => 11] = 1
    return Op[s' => 16, s => 15] = -1
end

function ITensors.op!(Op::ITensor, ::OpName"F", ::SiteType"2E", s::Index)
    Op[s' => 1, s => 1] = 1
    Op[s' => 2, s => 2] = -1
    Op[s' => 3, s => 3] = -1
    Op[s' => 4, s => 4] = -1
    Op[s' => 5, s => 5] = -1
    Op[s' => 6, s => 6] = 1
    Op[s' => 7, s => 7] = 1
    Op[s' => 8, s => 8] = 1
    Op[s' => 9, s => 9] = 1
    Op[s' => 10, s => 10] = 1
    Op[s' => 11, s => 11] = 1
    Op[s' => 12, s => 12] = -1
    Op[s' => 13, s => 13] = -1
    Op[s' => 14, s => 14] = -1
    Op[s' => 15, s => 15] = -1
    return Op[s' => 16, s => 16] = 1
end

```

I make sure to include the has\_fermion\_string for the necessary operators

```Julia
ITensors.has_fermion_string(::OpName"Cadagup", ::SiteType"2E") = true
ITensors.has_fermion_string(::OpName"Caup", ::SiteType"2E") = true

```

Then I construct the Hamiltonian (for the example I only include one ‘Electron’ site to each side of the central ‘2E’ site)

```julia
N = 3
sites = siteinds(n -> n!=2 ? "Electron" : "2E", N; conserve_qns=true, conserve_sz=false)
println(sites)
os = OpSum()
os += -1, "Cdagup", 1, "Caup", 2
os += -1, "Cadagup", 2, "Cup", 1
os += -1, "Cadagup", 2, "Cup", 1
os += -1, "Cdagup", 1, "Caup", 2

H = MPO(os, sites)
state = ["Up", "AupBup", "Emp"]
ψ = randomMPS(sites, state, 2)
sweeps = Sweeps(30)
maxdim!(sweeps, 20, 20, 30, 50, 100, 100, 100, 200, 200, 200, 200, 200, 400)
noise!(sweeps, 1E-6, 1E-7, 1E-9)
 cutoff!(sweeps, 1E-9)
E, ϕ = dmrg(H, ψ, sweeps)

```

Thank you.

---

<div class="post-metadata">

**Author:** ![miles](https://yyz2.discourse-cdn.com/free1/user_avatar/itensor.discourse.group/miles/32/6_2.png) [@miles](https://itensor.discourse.group/u/miles)\
**Post date:** [February 4, 2024, 12:18am UTC](https://itensor.discourse.group/t/multiorbital-electron-mixed-hamiltonian/1444/2 "2024-02-04T00:18:44Z")

</div>

Hi, thanks for your patience with our slow reply. It looks like you are using the library correctly and way to go about delving into some of the more advanced features. My best guess about the error you are getting is that your Hamiltonian may not actually conserve the set of quantum numbers you are wishing to either due to an error in typing it in, or a slight conceptual error somewhere.

One thing you can do is call the functions `checkflux` and `flux` on your initial wavefunction to check that it has a well-defined flux and it is the one you want.

Another thing to definitely do is to put

```julia
ITensor.enable_debug_checks()

```

at the top of your code which will turn on additional, possibly more informative error messages.

If none of these reveal the problem, please send me an email with a short, fully working example code that reproduces the issue and I can probably find it for you.

---

<div class="post-metadata">

**Author:** ![fegarcia](https://avatars.discourse-cdn.com/v4/letter/f/977dab/32.png) [@fegarcia](https://itensor.discourse.group/u/fegarcia)\
**Post date:** [February 13, 2024, 1:32am UTC](https://itensor.discourse.group/t/multiorbital-electron-mixed-hamiltonian/1444/3 "2024-02-13T01:32:58Z")

</div>

Hello, Miles. Thank you for the answer. Yes, the problem was spotted when I used the flux option. The error was very simple, namely I had named the QN for the new site different than the one for the electron site and thus two QN were defined and the Hamiltonian did not conserve them.

---

<div class="post-metadata">

**Author:** ![miles](https://yyz2.discourse-cdn.com/free1/user_avatar/itensor.discourse.group/miles/32/6_2.png) [@miles](https://itensor.discourse.group/u/miles)\
**Post date:** [February 13, 2024, 2:51am UTC](https://itensor.discourse.group/t/multiorbital-electron-mixed-hamiltonian/1444/4 "2024-02-13T02:51:25Z")

</div>

Glad to hear you found the issue!
