QAtlas.jl

_MAX_ED_SITES = 12

Hard cap on chain length for dense-ED helpers. 2^12 = 4096 → 4096×4096 Hermitian eigendecomposition takes ~1 s on a laptop; 2^14 ~ 16384 is already several seconds and 64 GB+ of workspace, which is past the "cheap reference data" regime this helper targets.

const AbstractModel = AbstractQAtlasModel

Backward-compatible alias. Existing downstream code dispatches on ::AbstractModel; new code should use ::AbstractQAtlasModel directly or — preferably — a concrete model struct.

AbstractQAtlasModel

Abstract parent type for every QAtlas model. Concrete subtypes carry their physics parameters as typed fields (e.g. struct TFIM <: AbstractQAtlasModel; J::Float64; h::Float64; end).

The older Model{S}(params::Dict) phantom-typed wrapper is still an AbstractQAtlasModel (via the deprecated alias below) but new models must use concrete structs.

AbstractQuantity

Abstract parent type for quantities. New code defines concrete structs (e.g. struct MagnetizationX <: AbstractQuantity end) so dispatch is static and naming is explicit (axis, entropy variant, …). The older Quantity{S} phantom-type wrapper is retained for legacy symbol-based dispatch; see the Quantity(::Symbol) shim below.

BoundaryCondition

Abstract parent type. The three concrete subtypes carry system-size information where applicable, so fetch can read it from the BC instead of kwargs.

• Infinite — thermodynamic limit; no size.
• PBC(N::Int) — periodic boundary conditions at finite N.
• OBC(N::Int) — open boundary conditions at finite N.

For backward compatibility, the zero-argument constructors PBC() and OBC() exist and set N = 0, which signals "caller will pass N via kwargs" — legacy fetch methods still look at kwargs[:N]. New fetch methods read bc.N directly.

CentralCharge() <: AbstractQuantity

Central charge c of the emergent CFT. For 1D critical systems extracted from the Calabrese–Cardy entanglement formula; universality pages return literature values.

CriticalExponents() <: AbstractQuantity

Standard set of equilibrium critical exponents {α, β, γ, δ, ν, η} of a universality class. Returns a NamedTuple.

For exact values: fields are Rational{Int}. For numerical estimates: fields are Float64 with corresponding _err fields (e.g., β_err) giving the uncertainty.

CriticalTemperature() <: AbstractQuantity

Exact critical temperature of a classical model.

E8() <: AbstractQAtlasModel

The E8 integrable field theory — the low-energy continuum theory of the 2D Ising model perturbed by a small longitudinal magnetic field at T = T_c (Zamolodchikov 1989). The only physics parameter is the underlying magnetic field perturbation strength, but QAtlas currently only exposes the universal mass-ratio spectrum, so this struct carries no fields.

E8Spectrum() <: AbstractQuantity

Zamolodchikov E8 mass spectrum (8 stable particles). Concrete implementation lives in src/universalities/E8.jl; the type is defined here so src/core/alias.jl can reference it without circular loads.

Energy() <: AbstractQuantity

Ground-state / thermal energy expectation. Per-site or total depends on the model's convention (documented on the model's fetch method).

EnergyLocal() <: AbstractQuantity

Bond-resolved energy density vector, length N_bulk − 1 for a bond Hamiltonian Σ_b h_b.

ExactSpectrum

Dispatch tag for the full sorted eigenvalue spectrum of a finite model.

FermiVelocity() <: AbstractQuantity

Fermi velocity v_F = ∂ε/∂k |_{k_F}. Meaningful for non-interacting / mean-field fermionic band structures (tight-binding lattices, Bogoliubov-de Gennes diagonalisations). In QAtlas this is the type returned by models like Honeycomb (at the Dirac cones), the other tight-binding lattices, and the TFIM Majorana mode at the critical field.

FidelitySusceptibility() <: AbstractQuantity

Fidelity susceptibility χ_F(λ) = −∂²⟨ψ(λ)|ψ(λ + δλ)⟩/∂δλ².

FreeEnergy() <: AbstractQuantity

Helmholtz free energy per site, f = -β⁻¹ log Z / N.

const Graphene = Honeycomb

Backward-compatible alias for Honeycomb. New code should prefer the lattice-named model.

GroundStateEnergyDensity

Dispatch tag for the ground-state energy per site in the thermodynamic limit (N → ∞).

GrowthExponents() <: AbstractQuantity

KPZ-type growth / roughness / dynamic exponents. Returns (β_growth, α_rough, z) instead of the equilibrium set.

Heisenberg1D

Dispatch tag for the spin-1/2 antiferromagnetic Heisenberg model on a 1D chain (or more generally any finite spin-1/2 cluster). Hamiltonian:

H = J Σ_{⟨i,j⟩} S_i · S_j,   spin-1/2, J > 0 antiferromagnetic

Honeycomb(; t::Real = 1.0, Lx::Int = 0, Ly::Int = 0) <: AbstractQAtlasModel

Nearest-neighbor tight-binding model on the honeycomb lattice (spinless, single-orbital) — historically known as graphene.

Hamiltonian: H = -t Σ_{⟨i,j⟩} (c†_i c_j + h.c.).

Fields:

• t — nearest-neighbor hopping amplitude (default 1.0).
• Lx, Ly — unit-cell counts in the two lattice directions. 0 is a sentinel meaning "caller passes via kwargs at fetch time".

For backward compatibility with the pre-v0.13 call form fetch(Graphene(), TightBindingSpectrum(); Lx, Ly, t=…), the fetch method still accepts Lx, Ly, t as kwargs (they override the struct fields when supplied). The legacy name Graphene is retained as a type alias in src/deprecate/legacy_honeycomb.jl.

Infinite()

Thermodynamic-limit boundary condition — no finite size.

Ising2D <: AbstractQAtlasModel

Dispatch tag for the two-dimensional Ising universality class. Acts as a convenience wrapper over Universality(:Ising) with the dimension pinned to d = 2, so call sites do not have to repeat d = 2 on every fetch:

QAtlas.fetch(Ising2D(), CriticalExponents())
# ≡ QAtlas.fetch(Universality(:Ising), CriticalExponents(); d = 2)

The struct is kept as a distinct nominal type (rather than const Ising2D = Universality{:Ising}) because the d kwarg is fixed at 2 here — a type-level const alias would defeat that fix. Subtyped to AbstractQAtlasModel so it composes with the canonical fetch(::AbstractQAtlasModel, ::AbstractQuantity, ::BoundaryCondition) signature.

See also: Universality, CriticalExponents.

IsingSquare(; J::Real = 1.0, Lx::Int = 0, Ly::Int = 0) <: AbstractQAtlasModel

Classical 2D Ising model on a square lattice with periodic boundary conditions (PBC) in both directions.

Hamiltonian: H = -J Σ_{⟨i,j⟩} σᵢ σⱼ, σᵢ ∈ {-1, +1}.

Physics parameters are carried as typed struct fields:

• J — Ising coupling (default 1.0, J > 0 ferromagnetic).
• Lx, Ly — lattice extents. 0 is a legacy sentinel meaning "thermodynamic limit / unspecified"; finite-size quantities like PartitionFunction require both to be positive.

For backward compatibility the fetch methods accept Lx, Ly, J, β as kwargs too — kwargs override the struct fields when both are supplied.

See also: PartitionFunction, CriticalTemperature, SpontaneousMagnetization.

KPZ1D <: AbstractQAtlasModel

Dispatch tag for the 1+1-dimensional KPZ (Kardar–Parisi–Zhang) universality class. Acts as a convenience wrapper over Universality(:KPZ) with the dimension pinned to d = 1:

QAtlas.fetch(KPZ1D(), CriticalExponents())
# ≡ QAtlas.fetch(Universality(:KPZ), GrowthExponents(); d = 1)

The struct is kept as a distinct nominal type (rather than const KPZ1D = Universality{:KPZ}) because the dimension is fixed here — a type-level alias would defeat that fix. Subtyped to AbstractQAtlasModel so dispatch composes with the canonical fetch(::AbstractQAtlasModel, ::AbstractQuantity, ::BoundaryCondition) signature.

See also: Universality, GrowthExponents.

Kagome

Dispatch tag for the nearest-neighbor tight-binding model on the kagome lattice (spinless, single-orbital, three sublattices).

Hamiltonian: H = -t Σ{⟨i,j⟩} (c†i c_j + h.c.)

The spectrum exhibits a dispersionless flat band at E = +2t touching the lower dispersive bands at the Γ-point.

KitaevHoneycomb(; Kx=1.0, Ky=1.0, Kz=1.0) <: AbstractQAtlasModel

Kitaev honeycomb model with coupling amplitudes on the three bond types. Kx = Ky = Kz is the isotropic gapless point. See module header for the full Hamiltonian and solution outline.

Lieb

Dispatch tag for the nearest-neighbor tight-binding model on the Lieb (line-centred square) lattice (spinless, single-orbital, three sublattices). Hamiltonian:

H = -t Σ_{⟨i,j⟩} (c†_i c_j + h.c.)

The spectrum contains a dispersionless flat band at E = 0 for every momentum, with an additional three-fold band touching at the M-point.

LuttingerParameter() <: AbstractQuantity

Luttinger liquid parameter K. Meaningful for critical 1D models with U(1) symmetry (e.g. XXZ in the critical regime |Δ| < 1).

LuttingerVelocity() <: AbstractQuantity

Luttinger-liquid / bosonisation velocity u (a.k.a. v_{LL}) of the low-energy linear-dispersion mode in a 1D critical interacting system. Used by models like XXZ1D in the Luttinger regime |Δ| < 1, the Heisenberg chain at the SU(2) point, and any other bosonised 1D critical theory.

For a free-fermion model this coincides with FermiVelocity; for interacting systems u includes the Luttinger renormalisation.

MagnetizationX() <: AbstractQuantity

Bulk-averaged ⟨σˣ⟩ in Pauli convention (= 2 ⟨Sˣ⟩ in spin-1/2 units). For a spin-1/2 chain H = -J ΣSᶻSᶻ - h ΣSˣ this is the transverse magnetization; the axis-explicit name avoids the "transverse" / "longitudinal" ambiguity that depends on the model's Hamiltonian choice.

MagnetizationXLocal() <: AbstractQuantity

Site-resolved ⟨σˣ_i⟩ vector of length N_bulk.

MagnetizationY() <: AbstractQuantity

Bulk-averaged ⟨σʸ⟩.

MagnetizationZ() <: AbstractQuantity

Bulk-averaged ⟨σᶻ⟩. For Z₂-symmetric phases on an infinite system this is the order parameter at low temperature; finite-system fetch methods may return the absolute value / the ordered-phase limit as documented.

MagnetizationZLocal() <: AbstractQuantity

Site-resolved ⟨σᶻ_i⟩ vector of length N_bulk.

MassGap() <: AbstractQuantity

Energy gap between the ground state and the first excited state.

MeanField() <: AbstractQAtlasModel

Mean-field (Landau) universality class. Exact for d ≥ d_c (upper critical dimension). Kept as a top-level alias so that existing fetch(MeanField(), ...) callers keep working after the v0.13 API redesign; the canonical form is Universality(:MeanField).

Model{M} <: AbstractQAtlasModel  (deprecated)

Phantom-typed Dict wrapper kept for backward compatibility. The Model(:TFIM; J=1.0, h=1.0) constructor below still works but is routed through the Symbol-dispatch deprecation shim in src/deprecate/legacy_fetch.jl. Prefer concrete model structs for new code.

OBC(N::Int)
OBC(; N::Int = 0)

Open boundary condition. N is the chain length. N = 0 is a legacy sentinel meaning "size unspecified — caller passes it via kwargs"; fetch methods that accept OBC(0) must look up kwargs[:N].

PBC(N::Int)
PBC(; N::Int = 0)

Periodic boundary condition. See OBC for the N = 0 sentinel.

PartitionFunction() <: AbstractQuantity

Canonical partition function Z = Σ_σ exp(-β H(σ)).

Quantity{Q} <: AbstractQuantity  (deprecated)

Phantom-typed wrapper kept for the legacy symbol API. New code should use concrete quantity structs such as Energy(), MagnetizationX(), ZZCorrelation(; mode=:static).

RenyiEntropy(α) <: AbstractQuantity

Rényi entropy of order α, S_α = (1 − α)⁻¹ log Tr ρ_A^α.

• α = 1 recovers VonNeumannEntropy (implementations may dispatch accordingly).
• α = 2 is the second Rényi entropy, frequently measured experimentally.
• α > 0, α ≠ 1 are the supported generic cases.

The inner constructor rejects α ≤ 0 and α = 1 (use VonNeumannEntropy() explicitly) — this is intentional, to force the call site to be explicit about which entropy it wants.

SpecificHeat() <: AbstractQuantity

Specific heat per site, c_v(β) = β² (⟨H²⟩ − ⟨H⟩²) / N.

const SpinWaveVelocity = LuttingerVelocity

Spin-chain community alias for LuttingerVelocity. The "spin wave velocity" (e.g. in the Haldane / Affleck literature on the AFM Heisenberg chain) is the same quantity as the Luttinger velocity once bosonised; both dispatch through the same fetch method via the type identity.

SpontaneousMagnetization() <: AbstractQuantity

Spontaneous magnetization of a classical model as a function of temperature. Retained under this name for backward compatibility with the classical-Ising literature; new code may prefer the axis-explicit MagnetizationZ together with a T < T_c context.

SusceptibilityXX() <: AbstractQuantity

Static transverse susceptibility, χ_xx(β) = β · (⟨M_x²⟩ − ⟨M_x⟩²) / N.

SusceptibilityYY() <: AbstractQuantity

Analogue for the y-axis.

SusceptibilityZZ() <: AbstractQuantity

Uniform longitudinal susceptibility, χ_zz(β) = β · (⟨M_z²⟩ − ⟨M_z⟩²) / N.

TFIM(; J = 1.0, h = 1.0) <: AbstractQAtlasModel

The 1D transverse field Ising model with Hamiltonian

H = -J Σ_i σᶻ_i σᶻ_{i+1} - h Σ_i σˣ_i

J > 0 is ferromagnetic, h is the transverse field. The critical point sits at h = J.

ThermalEntropy() <: AbstractQuantity

Thermal / thermodynamic entropy per site, s(β) = −∂f/∂T where f is the free energy per site. Real-valued, non-negative, monotone in T.

TightBindingSpectrum() <: AbstractQuantity

Single-particle Bloch spectrum of a tight-binding model. Returned as a sorted Vector{Float64} of length n_orbitals · Lx · Ly.

Triangular

Dispatch tag for the nearest-neighbor tight-binding model on the triangular lattice (spinless, single-orbital, one sublattice).

Hamiltonian: H = -t Σ{⟨i,j⟩} (c†i c_j + h.c.)

The single band ranges from -6t (at the Γ-point) to +3t (at the K-points), reflecting geometric frustration: the absence of particle-hole (chiral) symmetry on a non-bipartite lattice.

Universality{C}

Parametric dispatch tag for universality classes. C is a Symbol identifying the class (:Ising, :XY, :Heisenberg, :Potts3, :Potts4, :Percolation, :KPZ, etc.).

Use with CriticalExponents (equilibrium) or GrowthExponents (KPZ-type) and a d keyword to select the spatial dimension:

QAtlas.fetch(Universality(:Ising), CriticalExponents(); d=2)   # exact Rational
QAtlas.fetch(Universality(:Ising), CriticalExponents(); d=3)   # numerical + _err
QAtlas.fetch(Universality(:Ising), CriticalExponents(); d=4)   # mean-field

VonNeumannEntropy() <: AbstractQuantity

Von Neumann entanglement entropy of a reduced density matrix: S_vN = −Tr ρ_A log ρ_A. Requires a subsystem specification through the model's fetch kwargs (e.g. ℓ, the subsystem length).

XXCorrelation{M}() <: AbstractQuantity
XXCorrelation(; mode::Symbol = :static)

Real-space 2-point ⟨σˣ_i σˣ_j⟩ correlator. See ZZCorrelation for the mode semantics.

XXStructureFactor() <: AbstractQuantity

Fourier-space equivalent of XXCorrelation.

XXZ1D(; J::Real = 1.0, Δ::Real = 0.0) <: AbstractQAtlasModel

Spin-1/2 XXZ chain

H = J Σ_i [ S^x_i S^x_{i+1} + S^y_i S^y_{i+1} + Δ S^z_i S^z_{i+1} ]

Convention: J > 0 is antiferromagnetic. Δ = 1 is the isotropic Heisenberg AF point, Δ = 0 is the XX (free-fermion) point, Δ = -1 is the isotropic ferromagnet. For |Δ| < 1 the chain is critical (Luttinger liquid, central charge c = 1).

YYCorrelation{M}() <: AbstractQuantity
YYCorrelation(; mode::Symbol = :static)

Real-space 2-point ⟨σʸ_i σʸ_j⟩ correlator.

YYStructureFactor() <: AbstractQuantity

Fourier-space equivalent of YYCorrelation.

ZZCorrelation{M}() <: AbstractQuantity
ZZCorrelation(; mode::Symbol = :static)

Real-space 2-point correlator ⟨σᶻ_i σᶻ_j⟩. The mode M::Symbol is a phantom type parameter so dispatch can specialise on it.

Supported mode values (by convention; individual models need only implement the ones they support):

• :static — equal-time, thermal or zero-temperature value
• :connected — ⟨σᶻ_i σᶻ_j⟩ − ⟨σᶻ_i⟩⟨σᶻ_j⟩
• :dynamic — retarded real-time correlator ⟨σᶻ_i(t) σᶻ_j(0)⟩
• :lightcone — space-time spreading ⟨σᶻ_i(t) σᶻ_j(0)⟩ as a matrix over (site, time)

The companion type for Fourier-space structure factors is ZZStructureFactor, kept separate because it carries (q, ω) arguments instead of (i, j, t).

ZZStructureFactor() <: AbstractQuantity

Fourier-space structure factor S_zz(q, ω) = ∫ dt e^{iωt} (1/N) Σ_{ij} e^{iq·(i-j)} ⟨σᶻ_i(t)σᶻ_j(0)⟩ (or its static limit, depending on the model's fetch signature).

Kept as a separate type from ZZCorrelation because the argument domain is (q, ω) instead of (i, j, t) and because existing users already expect a dedicated StructureFactor quantity.

_bc_size(bc::BoundaryCondition, kwargs) -> Int

Return the effective system size for bc. Prefers bc.N when it is positive; otherwise looks up kwargs[:N]; otherwise throws. Legacy fetch methods can use this helper to accept both OBC(N=24) and OBC(); N=24 call forms.

_bc_with_legacy_N(bc::OBC, m::Model{:TFIM}) -> OBC

If bc.N == 0 (zero-arg OBC() call) but the legacy Model Dict has an :N entry, promote it to a sized OBC(N) so the concrete fetch can read bc.N uniformly.

_build_wick_matrix(idx_t, idx_0, Σ, R, RΣ) -> Matrix{ComplexF64}

Assemble the Wick contraction matrix F for a product

γ_{idx_t[1]}(t) γ_{idx_t[2]}(t) … γ_{idx_0[1]}(0) γ_{idx_0[2]}(0) …

from the ground-state Majorana covariance Σ, the time evolution matrix R = exp(h t), and the precomputed product RΣ = R * Σ. F is complex antisymmetric; its Pfaffian gives ⟨…⟩_GS.

_ed_thermal_energy(H::AbstractMatrix, β::Real) -> Float64

⟨H⟩_β = Tr(H exp(-βH)) / Tr(exp(-βH)). H must be Hermitian.

_ed_thermal_expectation(evals, evecs, O_diag, β) -> Float64

Return ⟨O⟩_β = Σ_n O_nn exp(-β eₙ) / Σ_n exp(-β eₙ) given a pre-computed eigendecomposition (evals, evecs) of H and the diagonal of O in the H-eigenbasis, O_diag[n] = ⟨n|O|n⟩.

Uses a log-sum-exp shift (eₘᵢₙ) so β up to ~10^2 / bandwidth is numerically stable.

_ising_sq_transfer_matrix(Ly, β, J) -> Matrix

Return the 2^Ly × 2^Ly symmetric transfer matrix for a row of Ly Ising spins with PBC along the row direction.

The symmetric (Hermitian) form is defined by:

T_{σ,σ'} = exp(βJ/2 · Eₕ(σ)) · exp(βJ · Eᵥ(σ,σ')) · exp(βJ/2 · Eₕ(σ'))

where Eₕ(σ) = Σⱼ₌₁^Ly σⱼ σ_{(j mod Ly)+1} (horizontal bonds within row, PBC) Eᵥ(σ,σ') = Σⱼ₌₁^Ly σⱼ σ'ⱼ (vertical bonds between rows)

Note: for Ly = 2, each horizontal bond is counted twice by the PBC sum (σ₁σ₂ + σ₂σ₁ = 2σ₁σ₂). The same applies along the transfer direction: tr(T^Lx) double-counts the vertical bonds when Lx = 2, because the cyclic product T[σ,σ'] T[σ',σ] weights the single pair of rows with the vertical Boltzmann factor twice. This doubling is an intrinsic property of the PBC transfer-matrix formula at L = 2 and is used consistently by the brute-force cross-check in test/util/classical_partition.jl (its bond list is constructed to match this convention by enumerating (i, j) ↔ (i, (j mod Ly) + 1) and (i, j) ↔ ((i mod Lx) + 1, j) for every site).

The function is generic in β and J so that automatic-differentiation dual numbers (e.g. ForwardDiff.Dual) propagate through the matrix elements. The element type of the returned matrix is inferred from typeof(exp(β * J)).

_kitaev_fk_abs²(m, θ₁, θ₂) -> Float64

|f(k)|² = Kₓ² + Kᵧ² + K_z² + 2Kₓ Kᵧ cos θ₁ + 2Kᵧ K_z cos θ₂ + 2K_z Kₓ cos(θ₁ − θ₂).

_kitaev_pbc_sector_energy(m, Lx, Ly, νx, νy) -> Float64

Per-site ground state energy on the Lx × Ly torus within a fixed topological flux sector (νx, νy) ∈ {0, 1/2}². νx = 1/2 flips every bond that crosses the m1 = Lx → m1 = 1 seam (anti-periodic fermion boundary condition along the x-direction); νy = 1/2 does the same in the y-direction. The four sectors correspond to the four inequivalent choices (W_x, W_y) ∈ {±1}² of the two non-contractible Wilson-loop operators.

Built as a singular-value decomposition of the Majorana hopping matrix M with entries ±K_γ, where the sign flips whenever a bond of type γ crosses one of the selected seams. Per-site energy = −Σ σ_k / (2·Lx·Ly).

_majorana_covariance_gs(h::AbstractMatrix) -> Matrix{Float64}

Ground-state Majorana covariance Σ of a quadratic Hamiltonian whose Majorana matrix is h. Σ is real antisymmetric, with ⟨γ_a γ_b⟩_GS = δ_{ab} + i Σ_{ab}. The formula Σ = -i · sign(i h) is used, evaluated through an eigendecomposition of the Hermitian matrix i h.

_majorana_evolution(h::AbstractMatrix, t::Real) -> Matrix{Float64}

exp(h * t) for a real antisymmetric h, returned as a real orthogonal matrix. Tiny imaginary noise from the matrix exponential is projected away.

_majorana_ham(N, J, h) -> Matrix{Float64}

Return the 2N × 2N real antisymmetric Majorana Hamiltonian matrix h of the OBC TFIM such that H = (i/4) Σ h_{ab} γ_a γ_b.

_majorana_thermal_covariance(h::AbstractMatrix, β::Real) -> Matrix{Float64}

Thermal-equilibrium Majorana covariance at inverse temperature β. The T = 0 ground-state formula Σ_GS = -i sign(i h) generalises to

Σ(β) = -i tanh((β/2) · i h)

so that ⟨γ_a γ_b⟩_β = δ_{ab} + i Σ(β)_{ab}. In the canonical 2×2 BdG block this gives the off-diagonal entry tanh(βΛ/2), recovering the Fermi–Dirac population of each quasiparticle and reducing to ±1 as β → ∞.

β = Inf falls back to _majorana_covariance_gs.

_obc_hopping_matrix(model, Lx, Ly) -> Matrix{Float64}

Lx·Ly × Lx·Ly bipartite hopping matrix for the flux-free-sector Majorana problem with open boundary conditions in both lattice directions. M[a, b] = K_γ whenever the A-sublattice site of unit cell a connects to the B-sublattice site of unit cell b via a bond of type γ (matching Lattice2D's Honeycomb topology conventions: type1 ↔ Kz within the same cell, type2 ↔ Kx crossing (+a₁, −a₂), type3 ↔ Ky crossing (−a₂)). Bonds that would leave the Lx × Ly strip are dropped.

Row indexing: a = (m1 − 1) Ly + m2 for (m1, m2) ∈ 1..Lx × 1..Ly.

_pauli_string(N::Int, site_ops::Pair{Int,Matrix{ComplexF64}}...) -> Matrix{ComplexF64}

Return the 2^N × 2^N tensor product that places each σ_k at its listed site and the identity elsewhere. site_ops is an arbitrary sequence of (site => matrix) pairs; missing sites get _σ0.

# σˣ at site 3 in an N=5 chain:
_pauli_string(5, 3 => _σx)

# σˣᵢ σᶻⱼ for (i,j) = (2,4):
_pauli_string(5, 2 => _σx, 4 => _σz)

_peschel_mode_entropy(ν::Float64) -> Float64

Per-mode von Neumann entropy

s(ν) = -[(1-ν)/2 · ln((1-ν)/2) + (1+ν)/2 · ln((1+ν)/2)]

as a function of the Majorana covariance singular value ν ∈ [0, 1]. Numerically stable at ν → 0 (s → ln 2) and ν → 1 (s → 0).

_sx_expect(Σ, i) -> Float64

⟨σˣ_i⟩_GS = Σ[2i-1, 2i], since σˣ_i = -i γ_{2i-1} γ_{2i}.

_sx_sx_corr(N, J, h, i, j, t; β = Inf) -> ComplexF64

⟨σˣ_i(t) σˣ_j(0)⟩_β for the OBC TFIM. Reduces to a 4×4 Pfaffian since σˣ_k = -i γ_{2k-1} γ_{2k}. β = Inf ⇒ ground state.

_sz_sz_corr(N, J, h, i, j, t; β = Inf) -> ComplexF64

⟨σᶻ_i(t) σᶻ_j(0)⟩_β for the OBC TFIM at inverse temperature β (default Inf ⇒ ground state).

Implementation: Wick / Pfaffian over the (2i-1) + (2j-1) = 2(i+j)-2 Majorana operators constituting the product. The overall phase is (-i)^{i+j-2}. The thermal generalisation enters only through the Majorana 2-point function

⟨γ_a γ_b⟩_β = δ_{ab} + i Σ(β)_{ab},   Σ(β) = -i tanh((β/2) i h),

so the body of the Pfaffian assembly is unchanged.

_sz_sz_spreading(N, J, h, center, times) -> Matrix{ComplexF64}

Return C[it, ix] = ⟨σᶻ_ix(t_it) σᶻ_center(0)⟩_GS for ix ∈ 1:N, t_it ∈ 1:length(times). Reuses the Majorana Hamiltonian, covariance, and (per time-step) evolution matrix, so the cost is O(length(times) · N · M³) with M = 2(center + N) - 2.

_sz_sz_static_thermal(N, J, h, β; i = nothing, j = nothing) -> Matrix{Float64}

Static ⟨σᶻ_i σᶻ_j⟩_β for the OBC TFIM at inverse temperature β. If both i and j are given, returns a single value (wrapped in a 1×1 matrix); otherwise returns the full N×N matrix of equal-time correlators.

_tfim_bdg_spectrum(N, J, h) -> Vector{Float64}

Return the N positive BdG quasiparticle energies Λₙ > 0 for the OBC TFIM with N sites, Ising coupling J, and transverse field h.

The 2N×2N BdG matrix is: H_BdG = [[A, B]; [-B, -A]] where A (tridiagonal, symmetric) encodes hopping + onsite energy, and B (antisymmetric) encodes the pairing terms from JW transformation.

A_{ii}   = 2h
A_{i,i±1} = -J
B_{i,i+1} = +J,  B_{i+1,i} = -J

_tfim_dispersion(k, J, h) -> Float64

Single-particle BdG quasiparticle energy at momentum k for the TFIM with couplings J (Z-Z coupling) and h (transverse field):

Λ(k) = 2 √(J² + h² - 2 J h cos k).

_tfim_from_legacy_model(m::Model{:TFIM}) -> TFIM

Extract J, h from the legacy Model{:TFIM}(params) Dict and build the concrete struct.

_tfim_thermo_infinite(quantity, J, h, β) -> Float64

Compute one of the per-site thermodynamic potentials of the infinite TFIM at inverse temperature β. quantity is a Symbol from (:free_energy, :entropy, :specific_heat, :transverse_magnetization, :transverse_susceptibility).

The integrals are evaluated via adaptive Gauss-Kronrod quadrature.

_tfim_thermo_obc(quantity, N, J, h, β) -> Float64

Per-site thermodynamic potential for the OBC finite-N TFIM, computed by summing the contribution of each BdG quasiparticle mode.

The transverse magnetisation and its susceptibility require the full single-particle Bogoliubov coefficients, not just the spectrum, so this routine diagonalises the BdG matrix internally to obtain them.

_tfim_transverse_obc(quantity, N, J, h, β) -> Float64

Compute m_x or χ_xx per site for OBC finite N by direct site-resolved BdG expectation. Uses the matrix exponential formula

Σ(β) = -i tanh(β/2 · i h_BdG)

(see TFIM_dynamics.jl) and identifies

⟨σˣ_i⟩ = Σ[2i-1, 2i].

The transverse susceptibility is obtained by central numerical differentiation of the magnetisation with respect to h.

_xxz1d_from_legacy_model(m::Model{:XXZ1D}) -> XXZ1D

_xxz1d_hamiltonian_matrix(model::XXZ1D, N::Int) -> Matrix{ComplexF64}

Assemble the 2^N × 2^N OBC Hamiltonian

H = J Σᵢ [ Sˣ_i Sˣ_{i+1} + Sʸ_i Sʸ_{i+1} + Δ Sᶻ_i Sᶻ_{i+1} ]
  = (J/4) Σᵢ [ σˣ σˣ + σʸ σʸ + Δ σᶻ σᶻ ]

via explicit tensor products built from the Pauli primitives in src/core/dense_ed.jl. Capped by _MAX_ED_SITES.

_zz_static_structure_factor(N, J, h, β, q) -> Float64

S_zz(q, β) = (1/N) Σ_{i,j} e^{-i q (i-j)} ⟨σᶻ_i σᶻ_j⟩_β evaluated by direct double sum from the thermal Pfaffian correlator. For OBC the lattice lacks translation invariance so this is the boundary-aware definition; in the bulk of a long enough chain it converges to the translation-invariant value.

_zz_uniform_susceptibility(N, J, h, β) -> Float64

Uniform (q = 0) longitudinal susceptibility per site,

χ_zz(β) = (β/N) Σ_{i,j} ⟨σᶻ_i σᶻ_j⟩_β,

obtained from a finite-temperature direct double sum (assumes ⟨σᶻ⟩β = 0 in the OBC ground manifold of the TFIM, which holds for any h ≠ 0 finite N). This is the static isothermal susceptibility — the fluctuation-dissipation form `χ = β · ⟨M²⟩c / Nfor the magnetisationM = Σᵢ σᶻᵢ`.

fetch(::E8, ::E8Spectrum, ::Infinite) -> Vector{Float64}

Returns the analytical E8 mass spectrum [m₁, m₂, …, m₈] normalised by m₁ = 1. m₂/m₁ = 2 cos(π/5) = ϕ (golden ratio) is the hallmark of the E8 symmetry.

fetch(::Heisenberg1D, ::Energy, ::OBC; beta, J=1.0) -> Float64

Per-site thermal energy ⟨H⟩_β / N for the spin-½ antiferromagnetic Heisenberg OBC chain at finite N (the isotropic point Δ = 1 of XXZ1D). Routes through fetch(::XXZ1D, ::Energy, ::OBC).

Since Heisenberg1D currently carries no J field, callers must pass J as a kwarg (default J = 1.0). Downstream bridges (e.g. ITensorModels to_qatlas(::Heisenberg1D)) lose J on conversion; use XXZ1D(; J, Δ=1) directly if you need a non-unit coupling.

fetch(::Heisenberg1D, ::ExactSpectrum; N, J=1.0) -> Vector{Float64}

Return the sorted exact spectrum of the spin-1/2 Heisenberg Hamiltonian on an N-site chain or ring with boundary condition bc.

Supported cases

• N=2, bc=:OBC (dimer): [-3J/4, J/4, J/4, J/4] (singlet Es = -3J/4, triplet Et = J/4, three-fold degenerate).
• N=4, bc=:PBC (4-site ring): [-2J, -J×3, 0×7, +J×5] Ground state E₀ = -2J (unique singlet). The ferromagnetic quintet sits at E = +J. The full degeneracy structure is: 1 singlet + 1 triplet + (1 singlet + 2 triplets at E=0) + 1 quintet.

Arguments

• N::Int: number of spin-1/2 sites
• J::Real: Heisenberg coupling constant (default 1.0; J > 0 AFM)
• bc::Symbol: boundary condition, :OBC (default) or :PBC

References

A. Auerbach, "Interacting Electrons and Quantum Magnetism" (1994), §2.
H. Bethe, Z. Physik 71, 205 (1931).

fetch(::Heisenberg1D, ::GroundStateEnergyDensity; J=1.0) -> Float64

Exact ground-state energy per site of the spin-1/2 antiferromagnetic Heisenberg chain in the thermodynamic limit (N → ∞, PBC):

e₀ = J (1/4 − ln 2) ≈ −0.4431 J

This is one of the earliest and most celebrated results of the Bethe ansatz. The derivation proceeds by solving the Bethe equations for the ground state of

H = J Σᵢ Sᵢ · Sᵢ₊₁

in the limit N → ∞, yielding a linear integral equation for the rapidity distribution whose solution gives the energy via integration.

Finite-size corrections

For a PBC chain of N sites, the ground-state energy density approaches e₀ with corrections of order 1/N² (logarithmic corrections also present):

E₀(N)/N = e₀ + O(1/N²)

See test/verification/test_universality_cross_check.jl for a finite-size extrapolation verification using ED at N = 4, 6, 8.

Arguments

• J::Real: Heisenberg coupling constant (default 1.0; J > 0 AFM)

References

H. Bethe, "Zur Theorie der Metalle. I. Eigenwerte und Eigenfunktionen
  der linearen Atomkette", Z. Physik 71, 205–226 (1931) — original
  Bethe ansatz solution.
L. Hulthén, "Über das Austauschproblem eines Kristalles",
  Ark. Mat. Astron. Fys. 26A, No. 11, 1–106 (1938) — first
  evaluation of e₀ = 1/4 − ln 2 from the Bethe equations.

fetch(::Ising2D, ::CriticalExponents) -> NamedTuple

Backward-compatible alias for fetch(Universality(:Ising), CriticalExponents(); d=2).

fetch(::IsingSquare, ::CriticalTemperature; J=1.0) -> Float64

Exact critical temperature of the 2D Ising model on the square lattice:

T_c = 2J / ln(1 + √2) ≈ 2.269 J

Equivalently, the critical reduced coupling is Kc = J/Tc = ln(1+√2)/2, or sinh(2K_c) = 1.

References

L. Onsager, "Crystal Statistics. I.", Phys. Rev. 65, 117 (1944).

fetch(::IsingSquare, ::PartitionFunction; Lx, Ly, β, J=1.0) -> Real

Exact partition function Z = Tr(T^Lx) for the classical 2D Ising model on an Lx × Ly square lattice with periodic boundary conditions in both directions.

The transfer matrix T acts along the Lx direction (row-to-row transfer), with each row containing Ly spins and PBC along the Ly direction.

Bond-counting convention

The transfer-matrix sum double-counts bonds along any dimension of length 2 (PBC wraparound of a length-2 ring produces the factor σ_1 σ_2 + σ_2 σ_1 = 2 σ_1 σ_2). The brute-force enumeration in test/util/classical_partition.jl is built under the same PBC convention so that Z_transfer-matrix == Z_bruteforce exactly for every (Lx, Ly, β, J). For Lx ≥ 3 and Ly ≥ 3 each physical bond is enumerated exactly once and the result coincides with the standard physical Z of a PBC lattice.

Special values

• β = 0 (any Lx, Ly, J): Z = 2^(Lx·Ly) — all configurations equally weighted
• J = 0 (any β, Lx, Ly): Z = 2^(Lx·Ly) — no interactions, same as β = 0

Automatic differentiation

fetch is generic in β and J so that ForwardDiff.Dual numbers propagate through it. This allows macroscopic thermodynamic quantities to be recovered from Z by differentiation — e.g.

⟨E⟩ = -∂(log Z)/∂β
C_v = β² · ∂²(log Z)/∂β²

See test/verification/test_ising_ad_thermodynamics.jl for a cross-check against direct ensemble averages.

Arguments

• Lx::Int: number of rows (transfer direction)
• Ly::Int: number of columns (row length, PBC)
• β::Real: inverse temperature (β = 1/(k_B T))
• J::Real: Ising coupling constant (default 1.0; J > 0 ferromagnetic)

References

L. Onsager, Phys. Rev. 65, 117 (1944).
B. M. McCoy and T. T. Wu, "The Two-Dimensional Ising Model" (1973).

fetch(::IsingSquare, ::SpontaneousMagnetization; β, J=1.0) -> Float64

Exact spontaneous magnetization of the 2D Ising model on the infinite square lattice:

M(T) = (1 − sinh⁻⁴(2βJ))^{1/8}    for T < T_c  (i.e. sinh(2βJ) > 1)
M(T) = 0                            for T ≥ T_c

The critical exponent β = 1/8 is visible in the approach M → 0 as T → T_c⁻.

Special values:

• T = 0 (β → ∞): M = 1 (fully ordered)
• T = T_c: M = 0 (onset of disorder)

Arguments

• β::Real: inverse temperature (β = 1/(k_B T))
• J::Real: Ising coupling constant (default 1.0; J > 0 ferromagnetic)

References

C. N. Yang, "The spontaneous magnetization of a two-dimensional Ising
model", Phys. Rev. 85, 808 (1952).

fetch(model::KitaevHoneycomb, ::Energy, ::Infinite; rtol=1e-8) -> Float64

Ground state energy per site in the thermodynamic limit:

ε_gs = − 1/(8π²) ∫₀^{2π} dθ₁ ∫₀^{2π} dθ₂ |f(θ₁, θ₂)|

Computed as a nested Gauss-Kronrod quadrature; rtol sets the outer-integral tolerance and 10× rtol the inner.

At the isotropic point Kx = Ky = Kz = 1 this returns ε_gs ≈ −0.78729862... per site (Baskaran–Mandal–Shankar 2007, Eq. 9); for Kitaev's original |K_γ| ≤ 1/2 convention the same call at Kx = Ky = Kz = 0.5 gives ε_gs ≈ −0.39364931..., half of the K = 1 value (H is linear in the couplings). Finite-size PBC sums converge to this TL value within ~10⁻³ at Lx = Ly = 8 and ~10⁻⁶ at Lx = Ly = 64 — see test/models/test_KitaevHoneycomb.jl.

fetch(model::KitaevHoneycomb, ::Energy, ::OBC; Lx, Ly) -> Float64

Per-site ground state energy on a Lx × Ly honeycomb strip with open boundaries in both lattice directions. Uses the flux-free-sector ansatz (u_{ij} = +1) and diagonalises the bipartite hopping matrix M returned by _obc_hopping_matrix. Singular values σ_k of M are the positive Majorana single-particle energies; ground state energy = −Σₖ σₖ; per site = that divided by 2·Lx·Ly.

bc.N is ignored; pass Lx, Ly explicitly.

fetch(model::KitaevHoneycomb, ::Energy, bc::PBC; Lx, Ly) -> Float64

Per-site ground state energy on a Lx × Ly unit-cell torus (PBC in both lattice directions) — enumerates all four topological flux sectors and returns the minimum. Each sector corresponds to a choice of fermion boundary conditions (W_x, W_y ∈ {±1}); Lieb's theorem fixes plaquette fluxes at +1, so the spin-Hamiltonian ground state is one of these four.

Bond connectivity matches Lattice2D.build_lattice(Honeycomb, Lx, Ly; boundary=PeriodicAxis()). bc.N is ignored; pass Lx, Ly as kwargs.

For large L the four sectors converge to the same energy and individual Bloch-sum terms dominate; for small L sector choice is essential (e.g. Lx = Ly = 2 gives distinct sector energies differing by ~10%).

fetch(model::KitaevHoneycomb, ::MassGap, ::Infinite) -> Float64

Single-Majorana gap in the thermodynamic limit.

Δ = 2 · min_k |f(k)|.

In the A/B/C gapless phase (each |Kᵧ| ≤ sum of the other two), f(k) has two linear (Dirac) zeros and Δ = 0. In the gapped Az / Ax / Ay phases (|Kᵧ| exceeds the sum of the other two), |f| is bounded away from zero and `Δ = 2·( |Kγmax| − |Kγother1| − |Kγ_other2| )`.

fetch(::MeanField, ::CriticalExponents) -> NamedTuple

Mean-field critical exponents (exact, Rational{Int}): α=0, β=1/2, γ=1, δ=3, ν=1/2, η=0.

Satisfies Rushbrooke, Widom, and Fisher scaling relations exactly.

fetch(model, quantity, bc; kwargs...)

Return the stored / computed value of quantity for model under boundary condition bc. The canonical signature takes a concrete model struct + concrete quantity struct + BC; a legacy fetch(::Symbol, ::Symbol, bc; kwargs...) shim is also provided in src/deprecate/legacy_fetch.jl for backward compatibility.

Each (model, quantity, bc) triple must be implemented as a separate method; this top-level definition throws an informative error for un-implemented triples.

fetch(::Graphene, ::TightBindingSpectrum; Lx, Ly, t=1.0) -> Vector{Float64}

Sorted single-particle spectrum of the nearest-neighbor tight-binding Hamiltonian on an Lx × Ly honeycomb lattice with periodic boundary conditions in both directions.

The spectrum is obtained in closed form by diagonalizing the 2×2 Bloch Hamiltonian at the Lx·Ly allowed momenta

k_{mn} = (m/Lx) b₁ + (n/Ly) b₂,   m ∈ {0,…,Lx−1}, n ∈ {0,…,Ly−1}

where (b₁, b₂) is the reciprocal basis dual to the real-space basis (a₁, a₂) used by Lattice2D's Honeycomb topology.

Explicitly, using k·a₁ = 2π m/Lx and k·a₂ = 2π n/Ly,

E_{mn,±} = ± t · √(3 + 2cos(2π m/Lx) + 2cos(2π n/Ly) + 2cos(2π(n/Ly − m/Lx)))

The chiral (sublattice) symmetry of the bipartite honeycomb lattice guarantees that the spectrum is symmetric about zero. A small negative floor is applied inside the square root to absorb rounding error at Dirac points where the argument vanishes exactly.

Arguments

• Lx::Int: number of unit cells in the first lattice direction
• Ly::Int: number of unit cells in the second lattice direction
• t::Real: nearest-neighbor hopping amplitude (default 1.0)

Return

A sorted Vector{Float64} of length 2·Lx·Ly.

References

P. R. Wallace, Phys. Rev. 71, 622 (1947).
A. H. Castro Neto et al., Rev. Mod. Phys. 81, 109 (2009).

fetch(::Kagome, ::TightBindingSpectrum; Lx, Ly, t=1.0) -> Vector{Float64}

Sorted single-particle spectrum of the nearest-neighbor tight-binding Hamiltonian on an Lx × Ly kagome lattice with periodic boundary conditions in both directions.

The 3×3 Bloch Hamiltonian is diagonalized numerically at each of the Lx·Ly allowed momenta

k_{mn} = (m/Lx) b₁ + (n/Ly) b₂,   m ∈ {0,…,Lx−1}, n ∈ {0,…,Ly−1}

with b₁, b₂ dual to Lattice2D's kagome primitive basis. In the θ parametrization (θ₁ = k·a₁ = 2π m/Lx, θ₂ = k·a₂ = 2π n/Ly),

H(k) = -2t · [ 0              cos(θ₁/2)           cos(θ₂/2)        ]
             [ cos(θ₁/2)      0                   cos((θ₂−θ₁)/2)   ]
             [ cos(θ₂/2)      cos((θ₂−θ₁)/2)      0                ]

Flat band

One band is exactly flat at +2t, so the returned spectrum contains at least Lx·Ly eigenvalues equal to 2t. At the Γ-point (k = 0) the upper dispersive band touches the flat band, contributing one extra +2t eigenvalue, giving a Lx·Ly + 1-fold degeneracy there. The remaining 2·Lx·Ly − 1 eigenvalues lie in [−4t, +2t).

Arguments

• Lx::Int: number of unit cells in the first lattice direction
• Ly::Int: number of unit cells in the second lattice direction
• t::Real: nearest-neighbor hopping amplitude (default 1.0)

Return

A sorted Vector{Float64} of length 3·Lx·Ly.

References

I. Syôzi, Prog. Theor. Phys. 6, 306 (1951).
D. L. Bergman et al., Phys. Rev. B 78, 125104 (2008).

fetch(::Lieb, ::TightBindingSpectrum; Lx, Ly, t=1.0) -> Vector{Float64}

Sorted single-particle spectrum of the nearest-neighbor tight-binding Hamiltonian on an Lx × Ly Lieb lattice with periodic boundary conditions in both directions.

The closed-form eigenvalues of the 3×3 Bloch Hamiltonian are

E_{mn,±} = ± 2t · √(cos²(π m/Lx) + cos²(π n/Ly)),
E_{mn,0} = 0

for m ∈ {0,…,Lx−1}, n ∈ {0,…,Ly−1}. The E = 0 contribution appears for every (m, n) and forms the Lieb flat band.

Flat band

For a generic (Lx, Ly) the spectrum contains exactly Lx·Ly zero eigenvalues from the flat band. When both Lx and Ly are even the M-point (θ₁, θ₂) = (π, π) becomes an allowed momentum; there the two dispersive bands also collapse to zero, contributing two extra zero modes (a three-fold band touching). The total number of E = 0 eigenvalues is therefore

Lx·Ly + (2 if both Lx and Ly are even else 0).

Arguments

• Lx::Int: number of unit cells in the first lattice direction
• Ly::Int: number of unit cells in the second lattice direction
• t::Real: nearest-neighbor hopping amplitude (default 1.0)

Return

A sorted Vector{Float64} of length 3·Lx·Ly.

References

E. H. Lieb, Phys. Rev. Lett. 62, 1201 (1989).

fetch(::Triangular, ::TightBindingSpectrum; Lx, Ly, t=1.0) -> Vector{Float64}

Sorted single-particle spectrum of the nearest-neighbor tight-binding Hamiltonian on an Lx × Ly triangular lattice with periodic boundary conditions in both directions.

Because the triangular lattice has a single sublattice per unit cell, each allowed momentum contributes exactly one eigenvalue

E_{mn} = -2t · [ cos(2π m/Lx) + cos(2π n/Ly) + cos(2π(n/Ly − m/Lx)) ]

for m ∈ {0,…,Lx−1}, n ∈ {0,…,Ly−1}. Note the lower cut-off is -6t (at the Γ-point) but the upper cut-off is only +3t (at the K-points (1/3, 2/3) and (2/3, 1/3) — allowed when Lx and Ly are divisible by 3).

Arguments

• Lx::Int: number of unit cells in the first lattice direction
• Ly::Int: number of unit cells in the second lattice direction
• t::Real: nearest-neighbor hopping amplitude (default 1.0; positive convention)

Return

A sorted Vector{Float64} of length Lx·Ly.

References

G. H. Wannier, Phys. Rev. 79, 357 (1950).

fetch(m::Symbol, q::Symbol, bc::BoundaryCondition = Infinite(); kwargs...)

Deprecated in v0.13. Prefer the concrete-struct form:

fetch(TFIM(; J=1.0, h=1.0), Energy(), OBC(N=24); beta=5.0)

This symbol-based signature is kept for drop-in backward compatibility — the call is routed through Model(m; kwargs...) + Quantity(q) and forwarded to the concrete fetch(::AbstractQAtlasModel, ::AbstractQuantity, ::BoundaryCondition) method registered by each model. A one-shot informational log is emitted the first time the shim is reached from a given site.

fetch(model::TFIM, ::CentralCharge, ::Infinite) -> Float64

Central charge of the TFIM critical point (h = J): c = 1/2 (Ising CFT). Returns NaN if not at the critical point (|h/J - 1| > 1e-6).

fetch(model::TFIM, ::Energy, ::Infinite; beta, betas) -> Float64 or Vector{Float64}

Energy per site ⟨H⟩/N in the thermodynamic limit (PBC, N → ∞).

ε(β) = -(1/π) ∫₀^π dk  Λ(k)/2 · tanh(β Λ(k) / 2)

where the PBC dispersion is Λ(k) = 2√(J² + h² - 2Jh cos k).

• beta::Float64: return scalar ε(β)
• betas::AbstractVector{Float64}: return vector
• no keyword: return ground-state energy per site (β → ∞)

Uses adaptive Gauss-Kronrod quadrature (QuadGK).

fetch(model::TFIM, ::Energy, bc::OBC; beta, betas) -> Float64 or Vector{Float64}

Total energy ⟨H⟩(β) for the OBC TFIM with N sites.

• N is read from bc.N (OBC(N) / OBC(; N)) or from kwargs[:N] as a legacy fallback.
• beta::Float64: return scalar ⟨H⟩(β)
• betas::AbstractVector{Float64}: return vector, reusing spectrum (O(N³) once)
• no keyword: return ground-state energy E₀ = -Σₙ Λₙ/2 (β → ∞)

Uses the exact BdG formula: ⟨H⟩ = -Σₙ (Λₙ/2) tanh(β Λₙ / 2)

fetch(model::TFIM, ::EnergyLocal, bc::OBC; beta::Float64, kwargs...)
    -> Vector{Float64}

Site-local energy density ε_i of the OBC TFIM at inverse temperature beta, defined so that Σᵢ ε_i = ⟨H⟩_β. Each bond is split symmetrically between its two endpoints:

ε_i = -(J/2) (⟨σᶻ_{i-1} σᶻ_i⟩_β + ⟨σᶻ_i σᶻ_{i+1}⟩_β) - h ⟨σˣ_i⟩_β

with the missing bonds at the i = 1 and i = N boundaries taken to be zero. Bond expectations are read off as Σ(β)[2i, 2i+1] from the Majorana thermal covariance (exact, O(N) after the single 2N×2N diagonalisation).

fetch(model::TFIM, ::FreeEnergy, ::Infinite; beta::Float64, kwargs...)

Per-site free_energy of the TFIM in the thermodynamic limit at inverse temperature beta. Uses adaptive Gauss-Kronrod quadrature over the BdG dispersion Λ(k) = 2√(J² + h² − 2Jh cos k).

fetch(model::TFIM, ::FreeEnergy, bc::OBC; beta::Float64, kwargs...)

Per-site free_energy of the OBC TFIM with N = bc.N sites at inverse temperature beta. Computed exactly via the BdG diagonalisation.

fetch(model::TFIM, ::MagnetizationX, ::Infinite; beta::Float64, kwargs...)

Per-site transverse_magnetization of the TFIM in the thermodynamic limit at inverse temperature beta. Uses adaptive Gauss-Kronrod quadrature over the BdG dispersion Λ(k) = 2√(J² + h² − 2Jh cos k).

fetch(model::TFIM, ::MagnetizationX, bc::OBC; beta::Float64, kwargs...)

Per-site transverse_magnetization of the OBC TFIM with N = bc.N sites at inverse temperature beta. Computed exactly via the BdG diagonalisation.

fetch(model::TFIM, ::MagnetizationXLocal, bc::OBC; beta::Float64, kwargs...)
    -> Vector{Float64}

Site-resolved transverse magnetisation [⟨σˣ_i⟩_β for i = 1:N] of the OBC TFIM at inverse temperature beta, read off from the Majorana thermal covariance as Σ[2i-1, 2i]. N is taken from bc.N.

beta = Inf falls back to the ground state.

fetch(model::TFIM, ::MagnetizationZLocal, bc::OBC; beta::Float64, kwargs...)
    -> Vector{Float64}

Site-resolved longitudinal magnetisation [⟨σᶻ_i⟩_β for i = 1:N]. Identically zero in the OBC TFIM by the Z₂ symmetry σᶻ → −σᶻ of the Hamiltonian (Gaussian state, odd product of Majoranas). Returned as an explicit zero vector so consumers can use it as an exact baseline against finite random-sample estimates that fluctuate around zero.

fetch(model::TFIM, ::MassGap, ::Infinite) -> Float64

Mass gap of the infinite-chain TFIM: the lowest single-quasiparticle excitation energy

Δ = min_k Λ(k),     Λ(k) = 2 √( J² + h² − 2 J h cos k ).

Closed form:

Δ = 2 |h − J|.

Canonical values:

• ordered (h < J): Δ = 2(J − h)
• disordered (h > J): Δ = 2(h − J)
• critical (h = J): Δ = 0 (Ising CFT, Δ ~ π v_F / N on finite chains)

fetch(model::TFIM, ::MassGap, bc::OBC) -> Float64

Single-quasiparticle gap of the N-site OBC TFIM read off the BdG spectrum as Λ_min, the smallest positive eigenvalue of the 2N×2N Bogoliubov-de Gennes Hamiltonian.

This is the one-particle excitation energy. Away from the critical point (|h − J| > O(1/N)) it converges to 2|h − J| exponentially in N. At the critical point h = J the OBC gap scales as Δ(N) ~ π J / N (Ising CFT).

Size is taken from bc.N (or kwargs[:N] as a legacy fallback).

fetch(model::TFIM, ::SpecificHeat, ::Infinite; beta::Float64, kwargs...)

Per-site specific_heat of the TFIM in the thermodynamic limit at inverse temperature beta. Uses adaptive Gauss-Kronrod quadrature over the BdG dispersion Λ(k) = 2√(J² + h² − 2Jh cos k).

fetch(model::TFIM, ::SpecificHeat, bc::OBC; beta::Float64, kwargs...)

Per-site specific_heat of the OBC TFIM with N = bc.N sites at inverse temperature beta. Computed exactly via the BdG diagonalisation.

fetch(model::TFIM, ::SusceptibilityXX, ::Infinite; beta::Float64, kwargs...)

Per-site transverse_susceptibility of the TFIM in the thermodynamic limit at inverse temperature beta. Uses adaptive Gauss-Kronrod quadrature over the BdG dispersion Λ(k) = 2√(J² + h² − 2Jh cos k).

fetch(model::TFIM, ::SusceptibilityXX, bc::OBC; beta::Float64, kwargs...)

Per-site transverse_susceptibility of the OBC TFIM with N = bc.N sites at inverse temperature beta. Computed exactly via the BdG diagonalisation.

fetch(model::TFIM, ::SusceptibilityZZ, bc::OBC;
      beta::Float64) -> Float64

Static uniform longitudinal (q = 0) susceptibility per site,

χ_zz(β) = (β/N) Σ_{i,j} ⟨σᶻ_i σᶻ_j⟩_β.

fetch(model::TFIM, ::ThermalEntropy, ::Infinite; beta::Float64, kwargs...)

Per-site entropy of the TFIM in the thermodynamic limit at inverse temperature beta. Uses adaptive Gauss-Kronrod quadrature over the BdG dispersion Λ(k) = 2√(J² + h² − 2Jh cos k).

fetch(model::TFIM, ::ThermalEntropy, bc::OBC; beta::Float64, kwargs...)

Per-site entropy of the OBC TFIM with N = bc.N sites at inverse temperature beta. Computed exactly via the BdG diagonalisation.

fetch(model::TFIM, ::VonNeumannEntropy, bc::OBC;
      ℓ::Int, beta::Float64 = Inf, kwargs...) -> Float64

Von Neumann entanglement entropy of the first ℓ spins of the N-site OBC TFIM in the thermal state at inverse temperature beta (or the ground state when beta = Inf), computed by Peschel's correlation- matrix method — see equation (1) in the file header.

• N = _bc_size(bc, kwargs) (read from OBC(N) or legacy kwargs[:N]).
• ℓ must satisfy 1 ≤ ℓ ≤ N - 1.
• Cost is O(ℓ³); for typical N = 200, ℓ = 100 this runs in a few milliseconds, whereas the full-ED SVD baseline scales as O(4^N).

The result matches the full-ED reference at every small N (verified to 1e-10 in test/models/test_TFIM_entanglement.jl).

See full derivation in docs/src/calc/jw-tfim-bdg.md.

fetch(model::TFIM, ::XXCorrelation{:dynamic}, bc::OBC;
      i::Int, j::Int, t::Float64, beta::Float64 = Inf) -> ComplexF64

Exact ⟨σˣ_i(t) σˣ_j(0)⟩_β for the OBC TFIM.

fetch(model::TFIM, ::ZZCorrelation{:dynamic}, bc::OBC;
      i::Int, j::Int, t::Float64, beta::Float64 = Inf) -> ComplexF64

Exact ⟨σᶻ_i(t) σᶻ_j(0)⟩_β for the OBC TFIM. beta = Inf (the default) gives the ground-state result. N comes from bc.N.

fetch(model::TFIM, ::ZZCorrelation{:lightcone}, bc::OBC;
      center::Int, times::AbstractVector{<:Real}, beta::Float64 = Inf) -> Matrix{ComplexF64}

Exact spreading correlation C[it, ix] = ⟨σᶻ_ix(t_it) σᶻ_center(0)⟩_β for all sites ix ∈ 1:N and t_it ∈ times.

fetch(model::TFIM, ::ZZCorrelation{:static}, bc::OBC;
      beta::Float64, [i::Int, j::Int]) -> Matrix{Float64} or Float64

Static (equal-time) thermal correlator ⟨σᶻ_i σᶻ_j⟩_β for the OBC TFIM. With both i and j given returns a scalar; otherwise returns the full N×N matrix.

fetch(model::TFIM, ::ZZStructureFactor, bc::OBC;
      beta::Float64, q::Real) -> Float64

Static structure factor S_zz(q, β) for the OBC TFIM at wave vector q.

fetch(::Universality{:Heisenberg}, ::CriticalExponents; d) -> NamedTuple

Critical exponents of the Heisenberg (O(3)) universality class.

• d ≤ 2: Mermin-Wagner theorem — no spontaneous order at T > 0.
• d = 3: Conformal bootstrap / Monte Carlo.
• d ≥ 4: Mean-field.

fetch(::Universality{:Ising}, ::CriticalExponents; d) -> NamedTuple

Critical exponents of the Ising universality class (Z₂ symmetry).

• d = 2: Exact rational values (CFT minimal model M(3,4), c = 1/2).
• d = 3: High-precision numerical estimates from the conformal bootstrap (Kos et al. 2016). Fields α_err, β_err, … give the uncertainty in the last digits.
• d ≥ 4: Mean-field (Landau) exponents (upper critical dimension).

fetch(::Universality{:KPZ}, ::GrowthExponents; d) -> NamedTuple

Scaling exponents of the KPZ universality class.

• d = 1 (1+1D): All three exponents are exact (Rational{Int}). Scaling relations: α + z = 2 (Galilean invariance), β = α / z.

Higher dimensions (d ≥ 2) are not exactly known and are not yet included.

fetch(::Universality{:Percolation}, ::CriticalExponents; d) -> NamedTuple

Critical exponents of the ordinary (site/bond) percolation universality class.

• d = 2: Exact rational values (Coulomb gas mapping).
• d = 3: Numerical (Monte Carlo).
• d ≥ 6: Mean-field.

fetch(::Universality{:Potts3}, ::CriticalExponents; d=2) -> NamedTuple

Exact critical exponents of the 3-state Potts model in d=2 (S₃ symmetry).

fetch(::Universality{:Potts4}, ::CriticalExponents; d=2) -> NamedTuple

Exact critical exponents of the 4-state Potts model (Ashkin–Teller point) in d=2 (S₄ symmetry). This is the marginal case: the transition is second-order with logarithmic corrections.

fetch(::Universality{:XY}, ::CriticalExponents; d) -> NamedTuple

Critical exponents of the XY (O(2)) universality class.

• d = 2: BKT transition — no standard power-law exponents. Returns η(T_c) = 1/4 only.
• d = 3: Conformal bootstrap (Chester et al. 2020).
• d ≥ 4: Mean-field.

fetch(model::XXZ1D, ::CentralCharge, ::Infinite) -> Float64

Central charge of the XXZ chain:

• -1 < Δ < 1 → c = 1 (Luttinger liquid)
• otherwise → NaN (non-critical)

fetch(model::XXZ1D, ::Energy, ::Infinite) -> Float64

Ground-state energy per site of the infinite XXZ chain in units of the Hamiltonian J. Currently exposes the three canonical values:

• Δ = -1 → -J/4 (isotropic FM, saturated)
• Δ = 0 → -J/π (XX, free fermion)
• Δ = 1 → J (1/4 - ln 2) (AF Heisenberg, Hulthén 1938)

For every other Δ a warning is emitted and NaN is returned — the general-Δ Bethe-ansatz integral is tracked as a v0.13 follow-up.

fetch(model::XXZ1D, ::Energy, ::OBC; beta) -> Float64

Per-site thermal energy ⟨H⟩_β / N for the spin-½ OBC chain at finite size, computed by dense ED. Works for any Δ and any N ≤ 12. Intended as a reference for MPS thermal methods (TPQMPS / Purification / METTS).

    ⟨H⟩_β / N = Tr(H exp(-βH)) / (N · Tr(exp(-βH)))

fetch(model::XXZ1D, ::GroundStateEnergyDensity, ::Infinite) -> Float64

Alias for fetch(::XXZ1D, ::Energy, ::Infinite) kept so that the GroundStateEnergyDensity quantity — already exported by Heisenberg.jl — works uniformly across 1D Bethe-ansatz chains.

fetch(model::XXZ1D, ::LuttingerParameter, ::Infinite) -> Float64

Luttinger-liquid parameter K = π / (2(π − γ)), with γ = arccos(Δ), valid for -1 < Δ ≤ 1.

Canonical values:

• Δ = 0 (XX free fermion) → K = 1
• Δ = 1 (AF Heisenberg) → K = 1/2
• Δ → -1 (FM boundary) → K → ∞

fetch(model::XXZ1D, ::LuttingerVelocity, ::Infinite) -> Float64
fetch(model::XXZ1D, ::SpinWaveVelocity,   ::Infinite) -> Float64

Sound velocity of the low-energy Luttinger-liquid mode,

u(Δ) = J · (π/2) · sin(γ)/γ,   γ = arccos(Δ).

Canonical values:

• Δ = 0 (XX) → u = J (= free-fermion v_F)
• Δ = 1 (AF) → u = (π/2) J (des Cloizeaux-Pearson)

SpinWaveVelocity dispatches here via the const SpinWaveVelocity = LuttingerVelocity type alias (both are the same physical quantity for 1D critical spin chains).

E8_MASS_RATIOS

Analytical expressions for the 8 particle masses in the E8 spectrum of the Ising field theory. Normalized by the lightest mass m₁.

References

• Theory: Zamolodchikov, A. B. (1989). Integrals of motion and S-matrix of the (scaled) T= Tc Ising model with magnetic field. International Journal of Modern Physics A, 4(16), 4235-4248.
• Review: Delfino, G. (2004). Integrable field theory and critical phenomena: the Ising model in a magnetic field. Journal of Physics A: Mathematical and General, 37(14), R45-R78.
• Experiment: Coldea, R., Tennant, D. A., Wheeler, E. M., Wawrzynska, E., Prabhakaran, D., Telling, M., ... & Kiefer, K. (2010). Quantum criticality in an Ising chain: experimental evidence for emergent E8 symmetry. Science, 327(5962), 177-180.

pfaffian(A::AbstractMatrix) -> eltype(A)

Compute the Pfaffian of an antisymmetric matrix A.

A must be square of even dimension. Antisymmetry (A[i,j] == -A[j,i]) is assumed; only the strictly upper (or lower) triangle is semantically used — diagonal entries are ignored.

For an odd-dimensional antisymmetric matrix the Pfaffian is zero by definition, and that is returned.

Works for both real and complex numeric element types. The algorithm is numerically stable with partial pivoting and runs in O(n³).