---
title: Semistochastic Heat-Bath CI (SHCI)
url: https://www.emergentmind.com/topics/semistochastic-heat-bath-configuration-interaction-shci-8ccedf34-15cf-4a21-87b2-a36c9146065b
type: topic
---

# Semistochastic Heat-Bath CI (SHCI)

Semistochastic Heat-Bath Configuration Interaction (SHCI) is a selected configuration interaction plus perturbation theory method that combines a heat-bath criterion for building a compact variational determinant space with a semistochastic evaluation of the multireference Epstein–Nesbet second-order correction. It is the semistochastic extension of Heat-Bath Configuration Interaction (HCI), systematically approaches the Full Configuration Interaction (FCI) limit as its thresholds are tightened, and is widely used as a high-accuracy classical reference for strongly correlated electronic structure, often alongside the Density Matrix Renormalization Group (DMRG) [1606.07453][1610.06660][2601.06935].

## 1. Origins and formal placement

SHCI belongs to the family of selected configuration interaction plus perturbation theory (SCI+PT) methods. In this family, one first constructs a compact variational expansion over a selected subset of determinants and then adds a multireference second-order perturbative correction from the enormous external space. HCI introduced the heat-bath selection rule as an efficient deterministic analogue of heat-bath sampling; SHCI retained that variational stage and replaced the fully deterministic perturbative step with a semistochastic algorithm designed to remove the severe memory bottleneck of the original method [1606.07453][1610.06660].

Within the broader SCI+PT landscape, HCI was reported to achieve sub-millihartree accuracy in variational spaces containing up to \(10^{20}\) determinants, while SHCI was described as allowing millihartree accuracy for active spaces of over \(100\) orbitals [2601.06935]. In later applications, SHCI was used for problems with effective FCI spaces as large as \(10^{38}\) determinants, while still delivering FCI-quality energies in the chosen one-electron basis [1803.05941][2511.15570]. These results explain why SHCI is repeatedly treated as a near-exact reference method for multireference molecules, transition-metal systems, heavy atoms, and benchmark Hamiltonian suites [2004.10059][2508.10873].

A central point in current practice is that SHCI is not merely an approximate solver of convenience. In several recent studies it is explicitly positioned as part of the classical “gold standard” for strongly correlated electrons, and as the baseline against which emerging quantum algorithms are judged [2601.06935][2508.10873].

## 2. Variational selected-CI stage and the heat-bath criterion

The variational SHCI wavefunction is a CI expansion over a selected determinant space \(V\),
$$
|\Psi_V\rangle = \sum_{D_i \in V} c_i |D_i\rangle,
$$
with the variational energy obtained by diagonalizing the Hamiltonian in that truncated space [2109.10271]. The practical question is therefore how to enlarge \(V\) efficiently.

HCI and SHCI use a heat-bath selection rule. A determinant \(D_a\) connected to the current variational space is added if there exists at least one \(D_i \in V\) such that
$$
|H_{ai} c_i| \ge \epsilon_1,
$$
where \(H_{ai}=\langle D_a|\hat H|D_i\rangle\) and \(\epsilon_1\) is the variational selection threshold [1606.07453][2109.10271]. In another equivalent formulation, the importance function is
$$
f_{\rm HCI}(D_a)=\max_{D_i\in V}|H_{ai}c_i|,
$$
and the determinant is selected if \(f_{\rm HCI}(D_a)>\epsilon_1\) [1606.07453].

This criterion differs from CIPSI-style selection rules that use an approximate first-order coefficient or second-order energy estimate. The HCI criterion discards the denominator and replaces the full sum by the largest product \(|H_{ai}c_i|\). The data describe this as a major efficiency gain: because most nonzero matrix elements arise from double excitations and depend only on the four orbitals involved, SHCI can precompute and sort these couplings and then generate only those excitations that can pass the cutoff [1606.07453][1708.07544]. In consequence, the algorithm does not spend time generating obviously unimportant determinants.

The determinant search is therefore integral-driven rather than determinant-driven. For each pair of occupied spin-orbitals, SHCI walks through a pre-sorted list of double excitations and stops once the matrix element is too small to satisfy the threshold. Single excitations are treated similarly, though they are less dominant in cost [1606.07453]. As \(\epsilon_1 \to 0\), the selected space approaches the full determinant space and the method reduces to FCI in the chosen basis [2511.15570].

For excited states and state-averaged calculations, the same logic is generalized to multiple roots. One form used in excited-state SHCI is
$$
\max_i \left| H_{ai} \; \max_n \left(|c_i^n|\right)\right| > \epsilon_{\rm V},
$$
so that a common variational space contains determinants important for at least one targeted state [1803.05941]. Relativistic and multistate variants also use state-averaged forms based on norms over several states rather than a single root, which is especially important when fine splittings or near-degeneracies are involved [2301.05794][2110.09946].

## 3. Semistochastic Epstein–Nesbet perturbation theory

After diagonalization in the selected space, SHCI adds a second-order Epstein–Nesbet correction from determinants outside \(V\). With
$$
|\Psi^{(0)}\rangle=\sum_{i\in V} c_i |D_i\rangle
$$
and variational energy \(E^{(0)}\), the standard SHCI correction is
$$
E^{(2)} = \sum_{a \notin V}
\frac{\left(\sum_{i \in V} H_{ai} c_i\right)^2}{E^{(0)} - H_{aa}},
$$
and the total energy is
$$
E_{\rm SHCI}=E^{(0)}+E^{(2)}.
$$
This expression appears throughout the SHCI literature and in later application papers as the canonical perturbative correction for the method [1610.06660][2601.06935][2511.15570].

A first practical simplification is to screen the inner sum with a perturbative threshold \(\epsilon_2\), retaining only terms with \(|H_{ai}c_i|>\epsilon_2\) [1606.07453][2109.10271]. The remaining challenge is that even the screened perturbative space can be too large to store explicitly. The original deterministic HCI PT2 step therefore became memory-limited for realistic \(\epsilon_2\) values, which directly motivated SHCI [1610.06660].

SHCI resolves this by sampling determinants from the variational wavefunction with probability
$$
p_i=\frac{|c_i|}{\sum_j |c_j|},
$$
using the Alias method rather than a Metropolis–Hastings random walk [1610.06660][2511.15570]. The perturbative energy is then evaluated semistochastically: a deterministic contribution is computed with a looser threshold, and a stochastic correction estimates the remaining tail. In its simplest form,
$$
\Delta E_2[\epsilon_2]
=
\left(\Delta E_2^S[\epsilon_2]-\Delta E_2^S[\epsilon_2^{\rm d}]\right)
+\Delta E_2^D[\epsilon_2^{\rm d}],
$$
where \(\epsilon_2^{\rm d}>\epsilon_2\), the superscript \(D\) denotes the deterministic part, and the two stochastic terms are evaluated with the same samples so that their fluctuations largely cancel [1610.06660].

Several practical properties follow directly from this construction. The perturbative calculation is embarrassingly parallel; there is no sign problem; independent Alias-method sampling avoids autocorrelation issues associated with Markov-chain approaches; and memory can be traded against computer time by varying \(\epsilon_2^{\rm d}\), \(\epsilon_2\), and the sampling effort [1610.06660]. For many systems, if a stochastic error of \(0.1\) mHa is acceptable, semistochastic PT2 is faster than the deterministic variant [1610.06660].

Modern large-scale SHCI workflows refine this further with deterministic, pseudo-stochastic, and fully stochastic PT2 stages controlled by \(\epsilon_2^{\rm dtm}\), \(\epsilon_2^{\rm psto}\), and \(\epsilon_2\) [2004.10059][2109.10271]. This suggests a mature implementation pattern: the largest contributions are handled exactly, intermediate contributions are estimated from sampled perturbative batches, and the smallest contributions are sampled over both variational and perturbative spaces.

## 4. Orbital optimization, extrapolation, and implementation practice

At finite \(\epsilon_1\), SHCI is not invariant under unitary orbital rotations. This has a major practical consequence: better orbitals can make the selected expansion dramatically more compact, reduce the PT2 correction, and improve extrapolation to the FCI limit [1708.07544][2104.02587].

Orbital optimization in SHCI resembles CASSCF in that both CI coefficients and orbitals are optimized, but it differs in a fundamental way: in SHCI there is no predefined inactive/active/virtual partition, so effectively all orbitals are active and most rotations are nonredundant [2104.02587]. The paper on orbital optimization in selected CI methods therefore distinguishes uncoupled, fully coupled, and quasi-fully coupled optimization strategies, and finds that taking CI–orbital coupling into account is crucial for fast convergence [2104.02587].

Two quasi-fully coupled methods are specifically recommended for SHCI applications: accelerated diagonal Newton and BFGS [2104.02587]. Starting from natural orbitals, these methods were shown to converge much faster than uncoupled orbital updates while avoiding the cost and nonconvexity problems of fully coupled Newton steps [2104.02587]. In parallel, HCISCF work established that active-space orbitals obtained from HCI-based optimization can be relatively insensitive to the accuracy of the underlying HCI solve, making cheap orbital optimization followed by tighter SHCI energetics a practical workflow [1708.07544].

Extrapolation is another standard part of SHCI methodology. In many applications, one computes \(E^{\rm SHCI}\) at several \(\epsilon_1\) values and fits the total energy against the residual perturbative correction. Weighted quadratic extrapolation in \(-\Delta E^{(2)}\) is used in thermochemistry and benchmark studies [2004.10059], while linear or quadratic extrapolations in related variables are used in excited-state, atomic, and multicomponent settings [1803.05941][2511.15570][2110.09946]. The common idea is that \(\Delta E^{(2)}\) acts as a proxy for the remaining correlation missing from the variational space.

Implementation-wise, SHCI has been used through the Arrow-SHCI code and closely related selected-CI infrastructures. In one HCISCF example, an Fe-porphyrin model complex with an active space of \((44e,44o)\) required 412 seconds per iteration on a single node containing 28 cores, of which 185 seconds were spent in the HCI calculation and 227 seconds mainly in integral transformation [1708.07544]. This suggests that, once the perturbative bottleneck is controlled, integral transformation and sparse Hamiltonian handling can dominate the wall time.

## 5. Benchmarking performance and scientific applications

SHCI has been applied across strongly correlated molecules, excited states, thermochemistry, transition-metal chemistry, heavy atoms, and bioinorganic benchmark systems. The following representative cases illustrate the method’s range.

| Domain | Representative result | Paper |
|---|---|---|
| Large active-space correlation | Better than \(1\) mHa accuracy for F\(_2\) \((14e,108o)\), Mn–Salen \((28e,22o)\), and Cr\(_2\) \((12e,190o)\) with wall times of 55 s, 37 s, and 56 min | [1610.06660] |
| Excited states | Hexatriene in ANO-L-pVDZ, with Hilbert space \(>10^{38}\), gave \(5.58\) and \(5.59\) eV for \(2^1A_g\) and \(1^1B_u\); ozone ring minimum lies \(1.3\) eV above the open-ring minimum and is separated by a \(1.11\) eV barrier | [1803.05941] |
| Thermochemistry | CBS-extrapolated atomization energies for the G2 set gave MADs of \(0.46\) kcal/mol without and \(0.51\) kcal/mol with the basis-set correction | [2004.10059] |
| Transition metals and monoxides | SHCI plus density-based basis-set correction was reported to converge total, ionization, and dissociation energies to the CBS limit within chemical accuracy | [2109.10271] |
| Heavy atoms | First ionization potentials of Cr, Mo, and W were reported as \(6.70 \pm 0.01\), \(7.076 \pm 0.002\), and \(8.124 \pm 0.004\) eV | [2511.15570] |
| Quantum-algorithm benchmarks | For \(N_2\) and \([2Fe\!-\!2S]\), HCI/SHCI serve as the classical reference, with HCI selected spaces reaching tens to hundreds of millions of determinants | [2601.06935] |

In strongly correlated bond breaking and transition-metal clusters, SHCI is used because it recovers both static and dynamic correlation in a systematically improvable way. In the \(N_2\) dissociation benchmark with an all-electron \((14e,28o)\) active space and a CAS size of \(1.40\times10^{12}\) determinants, HCI with \(\epsilon=10^{-6}\) was treated as essentially FCI-quality, while in the \([2Fe_2S_2(SCH_3)_4]^{2-}\) cluster with a \((30e,20o)\) active space and a CAS size of \(2.40\times10^8\) determinants, tighter HCI thresholds brought the energy within \(0.200\) mHa of CASCI [2601.06935].

In thermochemistry, the G2-set study demonstrates the distinction between “exact within basis” and the CBS limit. SHCI makes the former essentially available, and then basis extrapolation or basis-set correction becomes the dominant remaining issue [2004.10059]. In transition-metal atoms, ions, and monoxides, this division of labor is even more explicit: SHCI provides near-FCI energies in each basis, while density-based basis corrections accelerate convergence to the CBS limit [2109.10271].

The method is also noteworthy as a benchmark generator. In the QB ground-state energy benchmark, “optimized SHCI” was the only solver reported to attain a solvability score of \(1.0000\), with 148 tasks solved out of 226 attempted, and SHCI was used directly or indirectly to define FCI-quality references across the dataset [2508.10873].

## 6. Extensions, benchmarking role, and limits

SHCI has been extended well beyond nonrelativistic ground-state molecular CI. A relativistic formulation for arbitrary two-component and four-component Hamiltonians was presented for systems such as \(\mathrm{AuH_2^-}\) and \(\mathrm{NpO_2^{2+}}\), correlating more than \(100\) spinors in both cases [2301.05794]. Earlier one-step SOC treatments had already shown that HCI and its semistochastic extension could treat spin–orbit coupling and correlation on an equal footing in large active spaces; for the Au atom, converged excitation energies were obtained with just over \(10^7\) determinants in an active space containing \((150o,25e)\), whose full determinant count exceeds \(10^{30}\) [1710.00259].

Atomic-structure applications likewise broadened the scope of SHCI. In calculations of first ionization potentials of Cr, Mo, and W, SHCI was combined with orbital optimization, effective core potentials, \(\epsilon_1\)-based extrapolation to the FCI-in-basis limit, and basis-set extrapolation to the CBS limit [2511.15570]. The same paper emphasizes a practical parameterization in which
$$
\epsilon_2 = 10^{-3}\epsilon_1,\qquad
\epsilon_2^{\rm d}=10^{-2}\epsilon_1,
$$
linking perturbative thresholds to the variational threshold [2511.15570].

Methodological generalizations show that the underlying SHCI logic is portable. Multicomponent HCI for protonic excited states introduced deterministic Epstein–Nesbet PT2 in a product space of electronic and protonic determinants and explicitly identified a semistochastic extension as natural but not yet implemented there [2110.09946]. Vibrational heat-bath CI later imported the same semistochastic PT2 machinery into vibrational structure theory, reporting stochastic errors controllable to less than \(1\) cm\(^{-1}\) [2307.13246].

The method’s present limits are also clear. The original SHCI paper explicitly states that, after removing the PT2 memory bottleneck, the variational Hamiltonian becomes the dominant memory object [1610.06660]. Application papers and benchmarks likewise emphasize memory constraints and runtime blow-ups for very large active spaces, very strong correlation, and large perturbative spaces [2508.10873][2511.15570]. In one comparative discussion, CI diagonalization cost is stated to scale roughly as \(O(N^3)\) in the number of determinants \(N\), so reductions in selected-space size are directly consequential [2601.06935].

A further limitation is methodological rather than algorithmic: benchmark universality can be misleading. The QB benchmark reports near-universal solvability for optimized SHCI on its current dataset, but also states that many Hamiltonians originate from prior SHCI-oriented studies, introducing a bias that favors classical selected-CI solvers [2508.10873]. This means that SHCI’s apparent universality on such datasets should not be read as universal dominance across all chemically relevant or physically difficult Hamiltonians. Systems with extremely large active spaces, metallic behavior, long-range entanglement, or especially unfavorable determinant connectivity remain plausible hard cases [2508.10873].

Taken together, these developments define SHCI as a mature, systematically improvable SCI+PT framework: fast determinant selection by the heat-bath criterion, semistochastic multireference Epstein–Nesbet PT2, strong synergy with orbital optimization and extrapolation, and a demonstrated ability to deliver near-FCI reference data across electronic, relativistic, atomic, multicomponent, and vibrational settings [1610.06660][2104.02587][2301.05794].

Source: https://www.emergentmind.com/topics/semistochastic-heat-bath-configuration-interaction-shci-8ccedf34-15cf-4a21-87b2-a36c9146065b