As mentioned the number of CI-states grows exponentially, therefore the number of CI-states is potentially much larger than the number of molecular/spin-orbitals, or even the number of atom coordinates. Essentially the CI algorithm amounts to the diagonalization of this exponentially large many-body chemical Hamiltonian. There are two ways to remedy this: exploit sparsity/symmetry and reduce the wavefunction's search-space (Hilbert Space, the active-space). In practice both should be used. Because it is sparse, due to the sparsity of the SO integrals (defined above) and the Slater-Condon rules, we may save these entries in a convenient format: as a sparse-matrix (for off-diagonal parts) and a dense-vector (for diagonal parts). The algorithm above only requires B Iσp (which scales exponentially due to I CI-state index), and the ladder operators {a, a † } to compute the H CI matrix data entries, and hence are only transiently needed. Ultimately, the workflow begins with selecting the theory (active-space and reference(s) states) and the Hartree-Fock solution. Then using these to compute intermediaries as needed to yield the sparse-matrix {I K , J K , d K } (K goes over all SC values, not a CI-state index) and dense-vector