Exact Diagonalization Preview

Exact diagonalization, or ED, solves a declared finite many-body eigenproblem without introducing a physical approximation beyond the finite representation itself. One chooses a finite geometry, basis, local cutoff, boundary condition, and symmetry sector; represents the Hamiltonian on that space; computes all or selected eigenpairs; and evaluates observables from the resulting finite-system states.

A finite calculation retains finite-size, boundary, and local-cutoff errors:

$$\begin{gathered} \text{exact for the represented finite problem} \\ \neq \\ \text{exact for the thermodynamic system}. \end{gathered}$$

An eigensolver residual can be near machine precision while the physical result still has finite-size, boundary, cutoff, or model error. Conversely, a modest cluster can be uniquely valuable because it gives unbiased access to every state in the represented space, supplies sharp unit tests, and benchmarks methods that reach larger systems through additional approximations.

Required background. Computational method choice distinguishes numerical convergence from physical limits; matrix diagonalization introduces eigenpairs and residuals.

Complete and targeted finite diagonalization

The literature uses ED in two related senses.

Complete diagonalization

A dense or structure-aware solver returns every eigenvalue and, when requested, a complete orthonormal eigenbasis:

$$H = \sum_{n=1}^{D} E_n |n\rangle \langle n|.$$

This form supports exact finite-temperature traces, arbitrary unitary evolution in the finite space, full level statistics, and complete Lehmann sums. Its memory cost is at least quadratic if a dense Hamiltonian or all eigenvectors are stored.

Targeted finite diagonalization

In many-body practice, ED can also mean constructing the exact finite Hamiltonian but using a sparse iterative solver to obtain only a few eigenpairs. The finite operator is still represented without a variational ansatz, but the output is incomplete:

$$\left\{ E_0,\ldots,E_{k-1} \right\}, \qquad k\ll D.$$

This route is appropriate for ground states and low excitations. It does not automatically provide a complete thermal trace or all spectral weight. Lanczos develops the connection between the finite operator, its Krylov subspace, and the computed observables.

Whenever “ED” is reported, state which meaning is intended. Specify the finite geometry, boundary condition, Hilbert space and any local cutoff before constructing the matrix. These define the problem being solved.

Constructing the Basis

Let

$$\mathcal S_L = \operatorname{span} \left\{ |a\rangle \right\}_{a=0}^{D-1}$$

be an orthonormal finite basis. A computational basis needs both a physical definition and an indexing convention:

$$|a\rangle \longleftrightarrow \text{configuration} \longleftrightarrow \text{integer index}.$$

The map must be one-to-one on the represented sector. Duplicate configurations, omitted states, or inconsistent ordering can produce a Hermitian matrix with the wrong physics, so basis validation precedes diagonalization.

Spin-one-half bitstrings

For $L$ spin-one-half sites, choose

$$b_i = \begin{cases} 1, & |\uparrow\rangle_i,\\ 0, & |\downarrow\rangle_i. \end{cases}$$

A product state is

$$|b_{L-1}\cdots b_1b_0\rangle,$$

and one possible integer encoding is

$$a = \sum_{i=0}^{L-1} b_i2^i.$$

The site-to-bit convention must be stated. Reversing endianness does not change a correctly transformed spectrum, but it changes every operator mask, basis label, and reshaping convention. A one-site product state is an inexpensive test of the convention.

Use dimensionless spin operators $\mathbf S=\boldsymbol\sigma/2$ throughout this page; physical angular momentum is $\hbar\mathbf S$. The spin projection is

$$S_i^z |a\rangle = \left( b_i-\frac12 \right) |a\rangle.$$

If total $S^z$ is conserved, a fixed-$N_\uparrow$ sector contains only configurations satisfying

$$\sum_i b_i = N_\uparrow.$$

Its dimension is

$$D_{N_\uparrow} = \binom{L}{N_\uparrow}.$$

The counting derivation belongs to Scaling of Hilbert Space. Here the point is operational: the sector basis must contain exactly that many unique states, and every Hamiltonian action must remain inside it.

Fermionic occupation strings

For spinless fermionic modes ordered as $0,1,\ldots,M-1$, an occupation basis is

$$|\mathbf n\rangle = |n_{M-1}\cdots n_1n_0\rangle, \qquad n_i\in\{0,1\}.$$

The annihilation action includes the parity of occupied modes preceding $j$:

$$c_j|\mathbf n\rangle = (-1)^{p_j} n_j |\ldots,0_j,\ldots\rangle,$$

with

$$p_j = \sum_{i<j}n_i.$$

Changing the mode ordering changes intermediate signs and basis labels, although consistently transformed observables agree. On the full two-mode Fock space, check the anticommutator before assembling an interacting Hamiltonian:

$$\left( c_i c_j^\dagger + c_j^\dagger c_i \right) |\mathbf n\rangle = \delta_{ij} |\mathbf n\rangle.$$

The creation-string convention is $|\mathbf n\rangle=(c_0^\dagger)^{n_0}\cdots(c_{M-1}^\dagger)^{n_{M-1}}|0\rangle$, despite displaying the most significant bit first. Fermionic operator actions explain these signs. In a fixed-number sector, the individual $c_j$ maps to a different sector; projecting every factor back into the original sector would incorrectly make it vanish.

Truncated bosonic occupations

Bosonic local spaces are infinite, but finitely many modes at fixed total number already give a finite sector, with no additional local approximation. When number is not fixed, a finite representation requires a cutoff. A local cutoff uses

$$0 \leq n_i \leq n_{\max},$$

whereas a total-number cutoff uses

$$\sum_i n_i \leq N_{\max}.$$

These define different subspaces. The same symbol “cutoff 5” is meaningless without saying which rule is used. Observables sensitive to high occupation, such as onsite number variance, must be converged under cutoff enlargement.

The ladder factors

$$\begin{aligned} a_i|n_i\rangle &= \sqrt{n_i} |n_i-1\rangle, \\ a_i^\dagger|n_i\rangle &= \sqrt{n_i+1} |n_i+1\rangle \end{aligned}$$

must be combined with a boundary rule at $n_i=n_{\max}$. The projected ladder operator obeys $[a,a^\dagger]=I-(n_{\max}+1)|n_{\max}\rangle\langle n_{\max}|$, rather than the full oscillator commutator. Construct the desired Hamiltonian matrix elements before imposing the retained basis: generally $PABP\ne(PAP)(PBP)$ because the intermediate state can leave the basis and return. For example, $Paa^\dagger P$ acts as $n_{\max}+1$ on the top state, while the product of truncated ladders gives zero.

From Local Terms to a Sparse Hamiltonian

Write a local lattice Hamiltonian as

$$H = \sum_{X\subset\Lambda_L} h_X,$$

where each $h_X$ acts on a small set of sites or modes. Matrix elements are

$$H_{ab} = \langle a|H|b\rangle.$$

The useful computational object is often the action

$$H|b\rangle = \sum_a H_{ab}|a\rangle,$$

not a dense table of all $D^2$ possible entries. A local term changes only a few occupation labels, so each basis state connects to a small fraction of the full basis.

Define

$$\operatorname{nnz}(H) = \#\left\{ (a,b): H_{ab}\neq0 \right\}.$$

For many local models,

$$\operatorname{nnz}(H) \sim \operatorname{poly}(L)D,$$

even though $D$ itself grows exponentially. This makes sparse or matrix-free Hamiltonian application possible; it does not remove the exponential vector length.

Assembly as a configuration graph

It is useful to view the represented Hamiltonian as a weighted graph:

• basis states are vertices;
• diagonal terms give onsite vertex weights;
• off-diagonal terms connect configurations;
• matrix elements are edge weights;
• conserved labels diagonal in this configuration basis split the graph into disconnected components; translation or total-spin symmetry can instead require linear combinations of vertices.

This picture exposes several tests. Every transition must land on a valid basis state. Reverse edges must carry conjugate amplitudes. A claimed conserved sector must be closed under the full Hamiltonian. Termwise closure is a convenient sufficient condition; separate terms can also cancel leakage.

Hermiticity

For an orthonormal basis of a closed system,

$$H_{ab} = H_{ba}^\ast.$$

Common assembly errors include:

• adding a hopping or spin-flip transition in only one direction;
• double-counting an undirected bond;
• applying inconsistent complex phases on reverse links;
• using different basis-index conventions in diagonal and off-diagonal terms;
• forgetting a boundary bond or adding it twice.

Hermiticity is necessary but not sufficient. Two identically wrong off-diagonal entries can remain Hermitian.

Local Example: A Heisenberg Bond

For two spin-one-half sites,

$$H_{12} = J\, \mathbf S_1 \cdot \mathbf S_2.$$

Using dimensionless spin operators, $\mathbf S=\boldsymbol\sigma/2$, write

$$\mathbf S_1 \cdot \mathbf S_2 = S_1^zS_2^z + \frac12 \left( S_1^+S_2^- + S_1^-S_2^+ \right).$$

The local action is

$$\begin{aligned} \mathbf S_1\cdot\mathbf S_2 |\uparrow\uparrow\rangle &= \frac14 |\uparrow\uparrow\rangle, \\ \mathbf S_1\cdot\mathbf S_2 |\downarrow\downarrow\rangle &= \frac14 |\downarrow\downarrow\rangle, \end{aligned}$$

and

$$\begin{aligned} \mathbf S_1\cdot\mathbf S_2 |\uparrow\downarrow\rangle &= -\frac14 |\uparrow\downarrow\rangle + \frac12 |\downarrow\uparrow\rangle, \\ \mathbf S_1\cdot\mathbf S_2 |\downarrow\uparrow\rangle &= -\frac14 |\downarrow\uparrow\rangle + \frac12 |\uparrow\downarrow\rangle. \end{aligned}$$

These four rules are enough to assemble every Heisenberg bond in a product basis. The Heisenberg model develops the many-site physics.

The complete two-spin spectrum

Order the basis as

$$\left( |\uparrow\uparrow\rangle, |\uparrow\downarrow\rangle, |\downarrow\uparrow\rangle, |\downarrow\downarrow\rangle \right).$$

The matrix is

$$H_{12} = \frac{J}{4} \begin{pmatrix} 1 & 0 & 0 & 0\\ 0 & -1 & 2 & 0\\ 0 & 2 & -1 & 0\\ 0 & 0 & 0 & 1 \end{pmatrix}.$$

The parallel states are already eigenstates. The antiparallel block has symmetric and antisymmetric combinations:

$$\begin{aligned} |T_0\rangle &= \frac{ |\uparrow\downarrow\rangle + |\downarrow\uparrow\rangle }{\sqrt2}, & E_T &= \frac{J}{4}, \\ |S\rangle &= \frac{ |\uparrow\downarrow\rangle - |\downarrow\uparrow\rangle }{\sqrt2}, & E_S &= -\frac{3J}{4}. \end{aligned}$$

Together with $|\uparrow\uparrow\rangle$ and $|\downarrow\downarrow\rangle$, the symmetric states form a triplet at $J/4$ and the antisymmetric state is a singlet at $-3J/4$. Two Spin-One-Half Particles derives the angular-momentum addition. ED recovers it from basis actions.

This tiny problem supports a strong unit-test suite:

$$\operatorname{Tr}H_{12} = 0,$$ $$\operatorname{Tr}H_{12}^2 = \frac{3J^2}{4},$$

and

$$\left[ H_{12}, S_1^z+S_2^z \right] = 0.$$

A many-site implementation should reproduce the same local matrix elements before any large sparse calculation is trusted.

Dense, Sparse, and Matrix-Free Representations

The represented operator and the eigensolver are separate choices.

Dense complete route

Storing a general complex dense Hamiltonian requires approximately

$$M_{\mathrm{dense}} \simeq 16D^2 \quad \text{bytes}$$

in complex double precision, before workspace and eigenvectors. A complete dense Hermitian eigensolve typically requires order-$D^3$ arithmetic, with algorithm-dependent constants and storage reductions.

This route is valuable when $D$ is small enough because completeness enables:

• every finite-system energy level;
• exact traces over the represented space;
• arbitrary operator matrices in the eigenbasis;
• full finite-system propagators;
• complete level and degeneracy comparisons.

Sparse stored route

A sparse representation stores nonzero values and their indices. Its memory scales with $\operatorname{nnz}(H)$ rather than $D^2$, plus indexing overhead. It is useful when the operator is reused for many matrix-vector products or several observables.

Sparse storage does not imply a sparse eigenvector. Generic interacting eigenstates can have nonzero amplitude on essentially every basis state.

Matrix-free route

A matrix-free implementation computes

$$|y\rangle = H|x\rangle$$

by applying local terms directly. It avoids storing the sparse matrix and can reduce memory, but it is harder to inspect. Tiny-system matrices, randomized Hermiticity checks, and exact action tests become especially important.

Storage formats and operation counts are explained in Sparse Matrices.

Memory needed for vectors and matrices

A single complex state vector requires

$$M_{\mathrm{vec}} = 16D \quad \text{bytes}.$$

For an unconstrained $L=20$ spin-one-half system,

$$D = 2^{20} = 1{,}048{,}576.$$

Thus one vector is about $16$ MiB, while a dense complex matrix is about $16$ TiB:

$$\begin{aligned} M_{\mathrm{vec}} &\simeq 16\ \mathrm{MiB}, \\ M_{\mathrm{dense}} &\simeq 16\ \mathrm{TiB}. \end{aligned}$$

The fixed-$N_\uparrow=10$ sector has

$$D_{10} = \binom{20}{10} = 184{,}756.$$

One vector then needs about $2.82$ MiB, but a dense complex matrix still needs about $509$ GiB. Symmetry reduction can make sparse vector methods practical long after dense storage has failed.

Real calculations need several work vectors, basis-index data, sparse indices, observables, and solver workspace. Quoting only one-vector memory is a lower bound.

What a Complete Eigensystem Provides

Let

$$H|n\rangle = E_n|n\rangle, \qquad \langle m|n\rangle = \delta_{mn}.$$

Eigenstate observables

For an operator $O$,

$$\langle O\rangle_n = \langle n|O|n\rangle.$$

Off-diagonal matrix elements

$$O_{mn} = \langle m|O|n\rangle$$

control transitions, response, and dynamics. Operator conventions must match those used in the Hamiltonian basis.

Correlation functions

Equal-time correlations follow directly:

$$C_{ij}^{(n)} = \langle n| O_iO_j |n\rangle.$$

Connected correlations subtract one-point products. Their physical definitions and normalizations are explained in Connected Correlation Functions.

Exact finite-temperature traces

A complete spectrum gives

$$Z_L(\beta) = \sum_{n=1}^{D} e^{-\beta E_n},$$

and

$$\langle O\rangle_{\beta,L} = \frac{ \sum_n e^{-\beta E_n} \langle n|O|n\rangle }{ Z_L(\beta) }.$$

This is exact for the represented finite Hilbert space. It is not necessarily exact for an untruncated bosonic space or the thermodynamic system.

A low-energy subset can approximate low-temperature thermodynamics only when omitted weight is controlled. If $E_{\mathrm{cut}}$ is the first omitted energy, the suppression scale is

$$e^{-\beta \left( E_{\mathrm{cut}}-E_0 \right)},$$

but degeneracy and the number of omitted states also matter. “Low temperature” must be compared with both gaps and entropy.

An error bound for a partial thermal trace

For $\beta>0$, suppose the retained states are exact eigenstates of the finite Hamiltonian, with partition sum $Z_{\mathrm r}>0$. Suppose there are at most $R$ omitted states and every omitted energy is at least a certified $E_{\mathrm{cut}}$. Set

$$B=R e^{-\beta E_{\mathrm{cut}}}, \qquad Z_{\mathrm o}\le B, \qquad \delta=\frac{Z_{\mathrm o}}{Z_{\mathrm r}+Z_{\mathrm o}} \le\frac{B}{Z_{\mathrm r}+B}.$$

If no states are omitted, the error is zero. Otherwise, normalize the Gibbs weights separately on the retained and omitted subspaces to obtain $\rho_{\mathrm r}$ and $\rho_{\mathrm o}$. Their supports are orthogonal, and the full state is

$$\rho=(1-\delta)\rho_{\mathrm r}+\delta\rho_{\mathrm o}, \qquad \frac12\|\rho-\rho_{\mathrm r}\|_1=\delta.$$

Thus $\delta$ is exactly the trace distance of this spectral truncation. For a Hermitian observable with spectrum in $[o_{\min},o_{\max}]$,

$$|\langle O\rangle_\beta-\operatorname{Tr}(\rho_{\mathrm r}O)| \le(o_{\max}-o_{\min})\frac{B}{Z_{\mathrm r}+B}.$$

No commutation of $O$ with $H$ is required: each normalized state’s expectation lies in the same spectral interval. For a probability, use $0\le O\le I$ and interval width one. An operator-norm bound gives the less sharp but convenient width $2\|O\|$.

The retained free energy $F_{\mathrm r}=-\beta^{-1}\ln Z_{\mathrm r}$ also bounds the full free energy:

$$F_{\mathrm r}-\frac1\beta\ln\!\left(1+\frac{B}{Z_{\mathrm r}}\right) \le F\le F_{\mathrm r}.$$

These bounds quantify omitted thermal weight, separately from errors in computed eigenpairs or matrix elements. A few converged Ritz values do not by themselves certify that no omitted state lies below $E_{\mathrm{cut}}$. An energy shift rescales $B$ and $Z_{\mathrm r}$ equally, so their ratio and the observable bound remain unchanged.

Dynamical spectral sums

At zero temperature, a finite-system spectral measure for $O$ has the form

$$S_O(\omega) = \sum_n \left| \langle n|O|0\rangle \right|^2 \delta \left( \omega-\omega_{n0} \right),$$

where

$$\omega_{n0} = \frac{ E_n-E_0 }{\hbar}.$$

The exact finite result is a set of delta peaks. Any broadening used in a plot is an added resolution prescription. Spectral Functions defines the spectral conventions. Dynamical Correlation Functions Numerically shows how to benchmark this direct line list against Lanczos and real-time estimators with a matched kernel.

Entanglement from eigenvectors

First embed a state computed in a constrained or symmetry-adapted basis into the product basis (or assemble the coefficient matrix directly from its configurations). A fixed-number basis does not factor into subsystem bases. For a spatial bipartition $A\cup B$, the resulting coefficients can be arranged as

$$|\psi\rangle = \sum_{a,b} C_{ab} |a\rangle_A |b\rangle_B.$$

Then

$$\rho_A = CC^\dagger.$$

The arrangement depends on site ordering, bit convention, and subsystem choice. For fermionic modes, also fix the ordered Fock-space tensor convention and include the parity signs when reordering modes into the chosen grouping. Product states and Bell pairs are essential tests. Entanglement Spectrum develops this calculation and its interpretation.

Degeneracy and Basis Dependence

Inside an exactly degenerate eigenspace, the solver may return any orthonormal basis. Individual vectors can rotate under tiny perturbations or changes of numerical library while the invariant subspace remains correct.

If $\mathcal D$ is a degenerate subspace, define its projector:

$$P_{\mathcal D} = \sum_{n\in\mathcal D} |n\rangle\langle n|.$$

The projector is basis independent. For a symmetry-resolved interpretation, diagonalize the commuting symmetry operator within $\mathcal D$ or construct the sector before solving. Do not infer broken symmetry from an arbitrary numerical combination of degenerate states.

Eigenvector phases are also arbitrary:

$$|n\rangle \longrightarrow e^{i\phi_n}|n\rangle.$$

Expectation values are unchanged, but raw vector components and off-diagonal phases are not. Regression tests should compare phase-invariant quantities or align phases by an explicit convention.

Checking the eigensystem and observables

For every normalized computed eigenpair,

$$|r_n\rangle = \left( H-E_n \right) |n\rangle.$$

Report absolute or scaled residuals and check

$$\langle m|n\rangle \simeq \delta_{mn}.$$

For a complete eigensystem,

$$\sum_n E_n = \operatorname{Tr}H,$$

and

$$\sum_n E_n^2 = \operatorname{Tr}H^2.$$

These identities detect missing eigenvalues, sector mismatches, and some assembly errors. Higher moments can be useful when independently available.

Observable checks

• exact sum rules;
• positive probabilities and nonnegative variances;
• known limits at zero coupling or decoupled sites;
• symmetry-forbidden matrix elements, including transitions to other sectors;
• derivative identities such as Hellmann–Feynman;
• agreement between direct time evolution and spectral reconstruction on tiny systems.

For a nondegenerate eigenstate,

$$\frac{dE_n}{d\lambda} = \left\langle n\left| \frac{\partial H} {\partial\lambda} \right|n\right\rangle.$$

This tests the Hamiltonian parameterization, state, observable, and finite-difference procedure together.

Small-System Limits

ED has a controlled finite-system exactness and an uncontrolled reach problem. Several cautions follow.

Discrete spectra are not continua

Finite systems have discrete energy levels. A dense cluster of peaks can approximate a continuum only with a declared joint limit in size and resolution.

Symmetry breaking is subtle

Finite eigenstates can preserve an exact symmetry even where the thermodynamic system breaks it. Order parameters may vanish while squared order parameters, structure factors, susceptibilities, or quasi-degenerate state families reveal ordering tendencies.

Shape can matter as much as site count

Two clusters with the same $L$ can have different shortest loops, aspect ratios, coordination defects, or momentum grids. A sequence chosen only by increasing $L$ can oscillate between incompatible geometries.

Boundary conditions can select sectors

Periodic, open, twisted, and antiperiodic boundaries change allowed momenta and finite-size corrections. Boundary averaging can reduce some shell effects but does not erase the need to report each boundary condition.

Accessible sizes can precede the asymptotic regime

A smooth trend over three clusters is not proof of asymptotic scaling. Finite-Size Effects explains shape, shell, gap, boundary, and resolution effects; Thermodynamic Limit defines the limiting target.

Exercises

2. Estimate memory before allocation

For an unconstrained $L=24$ spin-one-half system, estimate the memory for one complex double-precision vector and for a dense complex matrix. Use binary units.

Solution

The dimension is

$$D = 2^{24}.$$

At $16$ bytes per complex number,

$$M_{\mathrm{vec}} = 16\cdot2^{24} \ \text{bytes} = 256\ \text{MiB}.$$

The dense matrix requires

$$\begin{aligned} M_{\mathrm{dense}} &= 16 \left( 2^{24} \right)^2 \ \text{bytes} \\ &= 16\cdot2^{48} \ \text{bytes} \\ &= 4\ \text{PiB}. \end{aligned}$$

This excludes eigensolver workspace and eigenvector storage. Sparse or matrix-free action is mandatory, but even that still requires several $256$ MiB vectors.

3. Test a claimed conserved sector

A Hamiltonian is assembled in a fixed-$N_\uparrow$ spin basis. Give two independent numerical tests that the sector is closed.

Solution

First, apply the full Hamiltonian to every sector basis state on a tiny system and verify that every nonzero output has the same $N_\uparrow$ and maps to a valid sector index. For a decomposition whose terms individually conserve number, testing each term is a useful stronger diagnostic.

Second, construct the number operator

$$N_\uparrow = \sum_i \left( S_i^z+\frac12 \right)$$

on the full tiny-system Hilbert space, before restricting to the sector, and verify

$$\left[ H,N_\uparrow \right] = 0$$

to numerical precision. The transition test exercises the basis map; the full-space commutator tests the original operator. Inside a fixed-$N_\uparrow$ block, $N_\uparrow$ is a scalar multiple of the identity, so its commutator is zero even for an incorrectly projected Hamiltonian. Equivalently, test the leakage $(I-P)HP$ in the full space, where $P$ projects onto the proposed block.

4. Bound a partial thermal trace

At positive inverse temperature $\beta>0$, suppose all omitted states have energy at least $E_{\mathrm{cut}}$ and there are at most $R$ omitted states. Give an upper bound on their partition-function contribution.

As a concrete test, retain only a nondegenerate ground state with $E_0=0$ and omit $R=1000$ states, all at energy $\Delta>0$. Compute the omitted probability at $\beta\Delta=8$. How large must $\beta\Delta$ be for every effect $0\le O\le I$ to have absolute expectation error at most $1\%$?

Solution

Each omitted Boltzmann factor satisfies

$$e^{-\beta E_n} \leq e^{-\beta E_{\mathrm{cut}}}.$$

Therefore the omitted contribution obeys

$$Z_{\mathrm{omit}} \leq R e^{-\beta E_{\mathrm{cut}}}.$$

Relative to the ground-state contribution,

$$\frac{ Z_{\mathrm{omit}} }{ e^{-\beta E_0} } \leq R e^{-\beta \left( E_{\mathrm{cut}}-E_0 \right)}.$$

The factor $R$ shows why a gap alone is not enough: the number of omitted states can offset Boltzmann suppression.

Here $Z_{\mathrm r}=1$ and $Z_{\mathrm o}=1000e^{-\beta\Delta}$, so the bound is attained:

$$\delta=\frac{1000e^{-8}}{1+1000e^{-8}} \simeq0.251196.$$

Each omitted state’s Boltzmann factor is small, but together the omitted states carry about $25.1\%$ of the probability. The projector onto the excited subspace attains this error. Requiring $\delta\le0.01$ gives

$$1000e^{-\beta\Delta}\le\frac{0.01}{0.99}, \qquad \beta\Delta\ge\ln(99000)\simeq11.5029.$$

The required inverse temperature includes the logarithm of the omitted degeneracy. A gap alone cannot set the thermal truncation tolerance.

5. Compare degenerate eigenspaces

Two diagonalization libraries return different orthonormal vectors for a twofold-degenerate eigenvalue. How should the results be compared?

Solution

Compare the projectors onto the degenerate subspaces:

$$P = |1\rangle\langle1| + |2\rangle\langle2|.$$

The individual vectors can be related by any $2\times2$ unitary rotation and arbitrary phases, so component-wise comparison is not meaningful. One can compare projector norms, principal angles between subspaces, traces of observables over the subspace, and symmetry-resolved vectors obtained by diagonalizing a commuting operator within it.

References

• P. Weinberg and M. Bukov, “QuSpin: A Python Package for Dynamics and Exact Diagonalisation of Quantum Many Body Systems, Part I: Spin Chains,” SciPost Physics 2, 003 (2017), doi:10.21468/SciPostPhys.2.1.003.
• E. Dagotto, “Correlated Electrons in High-Temperature Superconductors,” Reviews of Modern Physics 66, 763–840 (1994), doi:10.1103/RevModPhys.66.763.
• J. M. Zhang and R. X. Dong, “Exact Diagonalization: The Bose–Hubbard Model as an Example,” European Journal of Physics 31, 591–602 (2010), doi:10.1088/0143-0807/31/3/016.
• H. Fehske, R. Schneider, and A. Weiße, eds., Computational Many-Particle Physics, Springer, 2008, doi:10.1007/978-3-540-74686-7.
• C. Lanczos, “An Iteration Method for the Solution of the Eigenvalue Problem of Linear Differential and Integral Operators,” Journal of Research of the National Bureau of Standards 45, 255–282 (1950), doi:10.6028/jres.045.026.
• H. Q. Lin, “Exact Diagonalization of Quantum-Spin Models,” Physical Review B 42, 6561–6567 (1990), doi:10.1103/PhysRevB.42.6561.
• G. H. Golub and C. F. Van Loan, Matrix Computations, 4th ed., Johns Hopkins University Press, 2013.
• J. Jaklič and P. Prelovšek, “Finite-Temperature Properties of Doped Antiferromagnets,” Advances in Physics 49, 1–92 (2000), doi:10.1080/000187300243381.
• Y. Saad, Numerical Methods for Large Eigenvalue Problems, 2nd ed., SIAM, 2011, doi:10.1137/1.9781611970739.
• A. W. Sandvik, “Computational Studies of Quantum Spin Systems,” in AIP Conference Proceedings 1297, 135–338 (2010), doi:10.1063/1.3518900.