Quantum Numbers

Overview

In this notebook we will introduce the idea of quantum number conservation and how to enforce it in Renormalizer.

Quantum number conservation

In the context of Shuai group’s research, the most common case of quantum number conservation is the conservation of electron/exciton in the Holstein model.

\[\hat H = \sum_{ij} J_{ij} \hat a^\dagger_i \hat a_{i+1} + \sum_{ik}\frac{1}{2} (\hat p_{ik}^2 + \omega_k^2 \hat x_{ik}^2) + \sum_{ik} \hat a^\dagger_i \hat a_i \hat x_{ik}\]

Here, each electronic degree of freedom, described by \(\hat a^\dagger_i\) and \(\hat a\), is coupled with a set of harmonic oscillators indexed by \(k\).

For the demonstration purposes, we may neglect the vibrations and focus only on the electronic part, which leads us to the Hückel model (or tight binding model)

\[\hat H = \sum_{ij} J_{ij} \hat a^\dagger_i \hat a_{i+1}\]

This is a model that conserves the total particle number. Formally speaking, this means \(\hat H\) commutes with the total particle number operator \(\hat N = \sum_i \hat a^\dagger_i \hat a_i\). More intuitively, for any state \(|\psi\rangle\) with a particular particle number, \(\hat H |\psi\rangle\) has the same particle number. For example, suppose there are only two sites, and

\[|\psi\rangle = |01\rangle+|10\rangle\]

Thus,

\[\hat H |\psi\rangle = J_{12}|10\rangle+J_{21}|01\rangle\]

The two states with different quantum numbers, \(|00\rangle\) and \(|11\rangle\), will not present in \(\hat H |\psi\rangle\). In this model, the total particle number, is called a (good) quantum number.

There are two main reasons why we should care about quantum number - For many chemical applications the quantum number can be considered as an input to the model, just as the coefficients \(J_{ij}\). For example, the total number of exciton in a excitonic coupling model should be 1. If the conservation of quantum number is not enforced, numerical error may eventually lead us to incorrect solution. - Another advantage of using quantum number is that it saves memory and accelerates calculation. By looking at quantum numbers we can assert some of the matrix/vector elements must be zero. The sparsity can then be exploited for more efficient tensor manipulation, such as SVD.

Setting quantum number for states

In renormalizer, the first place to set the quantum number is basis sets. Most BasisSet classes have the sigmaqn argument that determines the quantum number of each basis.

from renormalizer import BasisSimpleElectron

2025-12-20 11:02:56,620[INFO] Use NumPy as backend
2025-12-20 11:02:56,621[INFO] numpy random seed is 9012
2025-12-20 11:02:56,621[INFO] random seed is 1092
2025-12-20 11:02:56,631[INFO] Git Commit Hash: 06430a9cbb7af930317f368a36eff512a8d5ce70
2025-12-20 11:02:56,632[INFO] use 64 bits

from renormalizer.utils.log import init_log, INFO
# set log level to INFO
init_log(INFO)

Here we set up two basis sets and sets the quantum number for \(|0\rangle\) to 0 and \(|1\rangle\) to 1.

b1 = BasisSimpleElectron(0, sigmaqn=[0, 1])
b2 = BasisSimpleElectron(1, sigmaqn=[0, 1])
basis = [b1, b2]

We next build a random MPS based on the basis sets and see its effect

from renormalizer import Mps, Model

model = Model(basis, ham_terms=[])

mps = Mps.random(model, qntot=1, m_max=2)

mps[0].array

array([[[1., 0.],
        [0., 1.]]])

mps[1].array

array([[[ 0.        ],
        [-0.7089592 ]],

       [[-0.70524949],
        [ 0.        ]]])

mps.todense()

array([ 0.        , -0.7089592 , -0.70524949,  0.        ])

We can see that although the MPS is random, half of the matrix elements are zero, due to particle number conservation. As a result, the overall dense state vector has a well-defined particle number of 1.

For comparison, in the following the MPS when quantum number is not activated is shown

model2 = Model([BasisSimpleElectron(i, sigmaqn=[0, 0]) for i in range(2)], ham_terms=[])
mps2 = Mps.random(model2, qntot=0, m_max=2)

mps2[0].array

array([[[-0.27201201, -0.96229386],
        [ 0.96229386, -0.27201201]]])

mps2[1].array

array([[[-0.56196169],
        [ 0.61166976]],

       [[ 0.34668986],
        [ 0.43573537]]])

mps2.todense()

array([-0.18075719, -0.58568699, -0.63507609,  0.47008079])

Setting quantum number for operators

Just as states, operators are also associated with quantum number. The quantum number of an operator shows the change of quantum number if the operator is applied to a state. For example, the quantum number of \(\hat a^\dagger\) is 1, since

\[\hat a^\dagger |0\rangle = |1\rangle\]

Similarly, the quantum number of \(\hat a\) is -1.

In Renormalizer, to fully take advantage of quatnum number conservation, we need to set the quantum number for the operators. The Op class accepts the qn argument for the quantum number of each elementary operator.

from renormalizer import Op, Mpo

Apply the creation operator to the MPS. The total quantum number of the MPS increases from 1 to 2, and the resulting state is \(|11\rangle\).

mps3 = Mpo(model, Op(r"a^\dagger", 0, qn=1)) @ mps

mps3.qntot

array([2])

mps3.todense()

array([ 0.       ,  0.       ,  0.       , -0.7089592])

For complex symbols the quantum number for each elementary symbol should be specified.

ham_terms = Op(r"a^\dagger a", [0, 1], qn=[1, -1]) + Op(r"a^\dagger a", [1, 0], qn=[1, -1])

To summarize, three places are relevant to the setup of quantum number - The basis sets - The operators - Total quantum number in MPS

Multiple Quantum Numbers

Renormalizer supports the conservation of multiple quantum numbers. For example, in ab initio electronic structure calculations, both the number of alpha spin electrons and the number of beta spin electrons are conserved. In such cases, quantum numbers should be set to a numpy array of integers.

In the following, we assume there are two basis sets "up" and "down" for spin-up electron and spin-down electron, respectively. Each basis set is a BasisHalfSpin instance.

We assign the quantum number 0 to the spin up state [1, 0] and the quantum number 1 to the spin down state [0, 1].

The creation operator \(a^\dagger_{\rm{up}}\) creates a spin-up electron and does not affect spin-down electron, and thus its quantum number of [1, 0].

Note that the quantum number of [1, 0] is not directly related to the quantum state [1, 0].

import numpy as np

# spin up creation operator
Op(r"a^\dagger", "up", qn=np.array([1, 0]))

Op('a^\\dagger', ['up'], 1.0, [[1, 0]])

Similarly \(a_{\rm{up}}\) annihilates a spin-up electron and does not affect spin-down electron, and thus its quantum number of [-1, 0].

# spin up annihilation operator
Op(r"a", "up", qn=np.array([-1, 0]))

Op('a', ['up'], 1.0, [[-1, 0]])

The creation and annihilation operators for the spin-down electron have the quantum number of [0, 1] and [0, -1] respectively.

# spin down creation operator
Op(r"a^\dagger", "down", qn=np.array([0, 1]))

Op('a^\\dagger', ['down'], 1.0, [[0, 1]])

# spin down annihilation operator
Op(r"a", "down", qn=np.array([0, -1]))

Op('a', ['down'], 1.0, [[0, -1]])

Disable Quantum Numbers

Setting quantum numbers is not necessary. Usually, quantum number conservation is important for lattice models. However, there are also numerous cases where quantum number is not conserved. Examples include the spin-boson model and many excited state dynamics models. In these cases quantum numbers should be disabled, otherwise the result could be catastropic.

Consider a toy model \(H=t|0\rangle\langle1| + t|1\rangle\langle0| + g|1\rangle\langle1|\). Quantum number is not conserved. If we set the quantum number incorrectly, we may get into trouble.

from renormalizer import BasisMultiElectron, optimize_mps

# Parameters
t = 1.0
g = 0.5

# Hamiltonian terms
# by default quantum numbers for "a^\dagger" and "a" are set to 1 and -1 respectively.
ham_terms = [
    Op(r"a^\dagger a", [0, 1], t),  # Hopping
    Op(r"a^\dagger a", [1, 0], t),  # Hopping
    Op(r"a^\dagger a", [1, 1], g),  # Site energy
]

# Define basis with quantum number conservation (WRONG for this system)
# |0> has quantum number 0, |1> has quantum number 1
basis_wrong = [
    BasisMultiElectron([0, 1], sigmaqn=[0, 1]),
]

model_wrong = Model(basis_wrong, ham_terms)
mpo_wrong = Mpo(model_wrong)

# DMRG optimization in the qn=0 quantum space
mps_wrong = Mps.random(model_wrong, qntot=0, m_max=10)
mps_wrong.optimize_config.method = "1site"
energies_wrong, _ = optimize_mps(mps_wrong, mpo_wrong)

2025-12-20 11:02:56,809[INFO] optimization method: 1site
2025-12-20 11:02:56,810[INFO] e_rtol: 1e-06
2025-12-20 11:02:56,810[INFO] e_atol: 1e-08
2025-12-20 11:02:56,812[INFO] procedure: [[10, 0.4], [20, 0.2], [30, 0.1], [40, 0], [40, 0]]
2025-12-20 11:02:56,817[INFO] DMRG has converged!
2025-12-20 11:02:56,818[INFO] mps current size: 16.0B, Matrix product bond dim:[1, 1]

The solution is limited to the state [1, 0] because of quantum number enforcement.

print(energies_wrong[-1])
print(mps_wrong.todense())

0.0
[1. 0.]

Next, we disable quantum number and perform the optimization.

In Renormalizer, quantum numbers for “a^:nbsphinx-math:dagger” and “a” are set to 1 and -1 respectively when initializing the Op instance. They must be explicitly set to 0 or set to None to disable their quantum numbers.

# disable quantum number
ham_terms = [
    Op(r"a^\dagger a", [0, 1], t, qn=[0, 0]),  # Hopping
    Op(r"a^\dagger a", [1, 0], t, qn=None),  # Hopping, equivalent
    Op(r"a^\dagger a", [1, 1], g, qn=[0, 0]),  # Site energy
]
# disable quantum number
basis_correct = [
    BasisMultiElectron([0, 1], sigmaqn=[0, 0]),
]

# perform the optimization as usual
model_correct = Model(basis_correct, ham_terms)
mpo_correct = Mpo(model_correct)
mps_correct = Mps.random(model_correct, qntot=0, m_max=10)
mps_correct.optimize_config.method = "1site"
energies_correct, _ = optimize_mps(mps_correct, mpo_correct)

2025-12-20 11:02:56,832[INFO] optimization method: 1site
2025-12-20 11:02:56,833[INFO] e_rtol: 1e-06
2025-12-20 11:02:56,834[INFO] e_atol: 1e-08
2025-12-20 11:02:56,834[INFO] procedure: [[10, 0.4], [20, 0.2], [30, 0.1], [40, 0], [40, 0]]
2025-12-20 11:02:56,840[INFO] DMRG has converged!
2025-12-20 11:02:56,841[INFO] mps current size: 16.0B, Matrix product bond dim:[1, 1]

We can obtain the correct solution in the full Hilbert space

print(energies_correct[-1])
print(mps_correct.todense())

-0.7807764064044151
[ 0.78820544 -0.61541221]

# reference value by exact diagonolization
np.linalg.eigvals(mpo_wrong.todense())

array([-0.78077641,  1.28077641])