SBD bindings (sbd)

Direct access to SBD’s two diagonalization methods. TPB (tensor-product basis) diagonalizes in the product of an alpha and a beta determinant list. GDB (general determinant basis) diagonalizes in an explicit list of determinants, and is experimental. For use inside a qiskit-addon-sqd loop, see sbd.sbd_solver instead.

The backend initializes itself on first use, so init() is needed only to choose a device or communicator explicitly.

Diagonalization

tpb_diag(fcidump, adet, bdet, sbd_data, loadname='', savename='', device=None)

Perform TPB diagonalization with data structures.

Parameters:

• fcidump – FCIDump object.

• adet – Alpha determinants.

• bdet – Beta determinants.

• sbd_data – TPB_SBD configuration object.

• loadname – Path prefix to load an initial wavefunction from (optional).

• savename – Path prefix to save the final wavefunction to (optional). See the note below: for reading amplitudes back, sbd_data.dump_matrix_form_wf is usually the one you want.

• device – Backend to run on: ‘cpu’, ‘gpu’ (Thrust), ‘gpu-omp’ (OpenMP offload), or an alias accepted by init(). None uses the default; see get_backend().

Returns:

energy, density, carryover_adet, carryover_bdet, one_p_rdm, two_p_rdm.

Return type:

dict with keys

Note

The returned dict holds no wavefunction amplitudes. There are two ways to get them, and they write different formats – do not mix them up:

sbd_data.dump_matrix_form_wf = "wf.bin" writes the full len(adet) x len(bdet) subspace, gathered onto one rank, as a bare row-major array of float64 with no header – so np.fromfile(path).reshape(len(adet), len(bdet)) reads it. Naming it .dat or .txt instead writes a text table in which every amplitude is labelled with its own alpha and beta bitstring, which is the easiest form to inspect by hand. This is the mechanism this package’s own SQD integration uses.

savename instead writes one file per b_comm position, f"{savename}{rank_b:06d}.bin", each holding only that rank’s block of the subspace: three size_t headers (adet_range, bdet_range, det_length), then that block’s alpha and then beta determinant words, then adet_range x bdet_range float64 amplitudes. It pairs with loadname for restarting a run, so it is the better choice for checkpointing and the worse one for analysis.

Amplitudes are ordered by adet/bdet as the diagonalizer saw them. sort_bitarray removes duplicates as well as sorting, so pass lists already in that form if you need to label the rows and columns yourself.

Both files use SBD’s sign convention, which orders a determinant’s creation operators orbital by orbital (alpha 0, beta 0, alpha 1, …). PySCF, and so qiskit-addon-sqd’s SCIState, puts all alpha operators before all beta ones, so some determinants differ in sign: multiply each amplitude by (-1) ** n, where n counts, over the occupied beta orbitals p, the occupied alpha orbitals above p, before handing it to PySCF. Magnitudes, and so selection by weight, are the same either way. sbd.sbd_solver.solve_sci() converts the amplitudes in its SCIState already.

tpb_diag_from_files(fcidumpfile, adetfile, sbd_data, loadname='', savename='', device=None)

Perform TPB diagonalization from files.

The one determinant file supplies both the alpha and the beta list, so the subspace is the product of that list with itself. For different alpha and beta lists, as an open-shell case needs, load them and call tpb_diag().

Parameters:

• fcidumpfile – Path to FCIDUMP file.

• adetfile – Path to the determinant file, used for both alpha and beta.

• sbd_data – TPB_SBD configuration object.

• loadname – Path prefix to load an initial wavefunction from (optional).

• savename – Path prefix to save the final wavefunction to (optional).

• device – Backend to run on: ‘cpu’, ‘gpu’ (Thrust), ‘gpu-omp’ (OpenMP offload), or an alias accepted by init(). None uses the default; see get_backend().

Returns:

energy, density, carryover_adet, carryover_bdet, one_p_rdm, two_p_rdm.

Return type:

dict with keys

Note

As with tpb_diag(), no amplitudes are returned in the dict; see that function for the two ways of writing them to disk.

gdb_diag(fcidump, det, sbd_data, loadname='', savename='', device=None, determinant_distribution='', determinant_grid_a=0, determinant_grid_b=0)

Perform GDB diagonalization over an explicit list of determinants.

Unlike TPB, which spans the subspace with the Cartesian product of an alpha and a beta determinant list, GDB spans it with the given determinants themselves, so an arbitrary sparse subspace can be diagonalized.

Each determinant is a 2 * norb-bit configuration packed into words of sbd_data.bit_length bits. Bit 2 * i is the occupation of spin-alpha orbital i and bit 2 * i + 1 that of spin-beta orbital i. from_strings() packs a list of bitstrings into the expected (ndets, words) array; a nested list is accepted and converted.

The shard contract. sbd_data.b_comm_size decides what det means:

• 1 — every rank passes the whole basis.

• > 1 — every rank passes its own shard, and the union over b_comm positions is the basis. This is the only way GDB’s memory scales, because b_comm is the only dimension that divides the basis (and with it the excitation lookup); the derived helper dimension divides work but not storage.

The shard index is rank % b_comm_size, and ranks sharing one must pass identical determinants — the in-memory path does not broadcast the list. Shards must be globally sorted and disjoint: shard i strictly below shard i + 1. All of this is checked, and a violation raises rather than silently diagonalizing the wrong subspace. What cannot be checked is completeness — that the union is the basis you meant — so compare the returned global_dim against what you expect.

t_comm_size must not exceed b_comm_size: GDB runs one task per basis-ring station and there are exactly b_comm_size of them. t_comm_size * b_comm_size must divide the rank count exactly. On the Thrust backend the helper dimension must be 1, i.e. every rank goes to t_comm_size * b_comm_size.

Parameters:

• fcidump – FCIDump object.

• det – Determinants spanning the subspace, or this rank’s shard of them, as an (ndets, words) array. Must be distinct.

• sbd_data – GDB_SBD configuration object. method must be 0 or 1 — GDB implements Davidson only, so TPB’s Lanczos methods 2 and 3 do not exist here and are rejected.

• loadname – Path to load an initial wavefunction (optional).

• savename – Path prefix to save the final wavefunction to (optional). SBD writes one file per b_comm position, f"{savename}{rank_b:06d}.bin", each holding only that shard; rank 0’s file is not the whole wavefunction. GDB has no combined matrix-form dump.

• device – Backend to run on: ‘cpu’, ‘gpu’ (Thrust), ‘gpu-omp’ (OpenMP offload), or an alias accepted by init(). None uses the default; see get_backend(). ‘gpu-omp’ is accepted, but SBD has no OpenMP-offload GDB code, so the run executes on the host; use ‘gpu’ for GDB on GPUs.

• determinant_distribution – How to place determinants across b_comm — one of input (keep the shards as given), equal-bra-a (default; equal distinct-alpha count per rank, which is what balances the matvec), count, count-sorted, grid-cyclic or grid-cyclic-balanced. Ignored when b_comm_size is 1.

• determinant_grid_a – Grid rows for the grid-cyclic schemes. With determinant_grid_b, must multiply to b_comm_size; give both or neither, and the default is the factor pair nearest square.

• determinant_grid_b – Grid columns; see determinant_grid_a.

Returns:

dict with keys energy, density, carryover_det, one_p_rdm, two_p_rdm, local_dim, global_dim and determinant_distribution. What is replicated and what is sharded:

• energy — replicated, bit-identical on every rank.

• density, one_p_rdm, two_p_rdm — replicated on every rank.

• carryover_det — this rank’s shard, as an (n, words) array, not the whole list. For carryover_type 1 it is split over b_comm and duplicated across the helper dimension, so gathering means taking one representative per rank % b_comm_size. For types 2 and 3 it is split over the world communicator with no duplication, so a plain allgather is correct; those types return the parents together with the new candidates, i.e. the next subspace rather than only the additions.

• local_dim / global_dim — determinants on this rank, and summed over b_comm.

Configuration

TPB_SBD(device=None)

Create a TPB_SBD configuration object, for tpb_diag().

Its fields are listed in the tables below.

GDB_SBD(device=None)

Create a GDB_SBD configuration object, for gdb_diag().

Its fields are listed in the tables below.

Both return a configuration object whose fields are set as attributes before it is passed to a diagonalization, for example config.eps = 1e-8. The defaults are upstream SBD’s; sbd.sbd_solver.solve_sci() and solve_sci_batch() set their own unless overridden through sbd_config.

Fields shared by TPB_SBD and GDB_SBD:

Field	Meaning
method	0 Davidson, 1 Davidson with the Hamiltonian stored. TPB also has 2 and 3, the Lanczos variants; GDB rejects them.
eps	Davidson stopping tolerance on the residual norm, not on the energy.
max_it	Maximum Davidson iterations. Reaching it before eps returns a partly converged vector without warning.
max_nb	Maximum number of Davidson basis vectors.
max_time	Wall-clock limit, in seconds.
init	Starting vector: 0 puts unit weight on the first determinant of the basis, 1 is random.
seed	Seed for the random starting vector; used only with init = 1.
do_rdm	0 returns the orbital density only, 1 also the 1- and 2-particle RDMs.
carryover_type	How SBD selects determinants for a next iteration. 0 is off. See the per-method tables below for what the other values do.
ratio, threshold	Carryover selection parameters. Which of the two a given carryover_type reads differs, and SBD silently ignores the other.
bit_length	Bits per word when packing determinants. At most 63; 64 is undefined behavior in SBD’s multi-rank redistribution. For GDB it must also be even whenever a determinant spans more than one word: an odd value breaks SBD’s alpha/beta conversion, which can crash or return a wrong energy, in plain diagonalization as well as in the heatbath expansion. The examples/gdb drivers require an even value of at most 62.

Fields of TPB_SBD only:

Field	Meaning
adet_comm_size, bdet_comm_size	MPI ranks spanning the alpha and beta determinant dimensions.
task_comm_size	MPI ranks spanning task-level parallelism. The rank count must be a multiple of task_comm_size * adet_comm_size * bdet_comm_size; the remainder becomes the derived helper dimension.
carryover_type	1 ranks half-determinants by marginal weight; 2 also adds their single excitations; 3 ranks whole determinants by amplitude, then adds single excitations.
do_shuffle	Shuffle the determinants loaded from file before use.
dump_matrix_form_wf	Path to which SBD writes the wavefunction as an |adet| x |bdet| matrix.
use_precalculated_dets	Thrust backend only. Precompute a determinant index for the whole subspace on the GPU: faster, but far more memory.
max_memory_gb_for_determinants	Thrust backend only, and used only when use_precalculated_dets is off: a per-thread buffer cap in GB, or -1 for none.

Fields of GDB_SBD only:

Field	Meaning
b_comm_size	MPI ranks the determinant basis is sharded across. It is the only dimension that divides memory. With more than 1, each rank passes only its own shard to gdb_diag().
t_comm_size	MPI ranks spanning tasks. Must not exceed b_comm_size. The rank count must be a multiple of t_comm_size * b_comm_size; on the Thrust backend it must equal it, since GDB on Thrust has no helper dimension.
carryover_type	1 truncates by weight; 2 and 3 expand by heatbath selection.
heatbath_cutoff	Heatbath selection cutoff: a candidate is added when |c_i H_ij| exceeds it.
heatbath_truncation	Weight below which parents are dropped before expansion. Large values remove the parents the expansion needs.
heatbath_batch_size	Number of candidates processed per batch during expansion.

Inputs

LoadFCIDump(filename, device=None)

Load an FCIDUMP file.

Parameters:

• filename – Path to the FCIDUMP file.

• device – Backend whose FCIDump type to build; None uses the default.

Returns:

The FCIDump object that tpb_diag() and gdb_diag() take.

LoadAlphaDets(filename, bit_length, total_bit_length, device=None)

Load determinants from a file of bitstrings, one per line.

Each line is packed as from_strings() does. Upstream sorts the determinants and removes duplicates, so the result is ready to use as a TPB determinant list and may hold fewer determinants than the file has lines.

The format is chosen by extension: .txt for one bitstring per line, .bin for upstream’s binary format. Any other extension returns an empty list without an error.

Parameters:

• filename – Path to the file.

• bit_length – Bits per packed word. Must match the config’s bit_length.

• total_bit_length – Length of each bitstring, norb for the alpha or beta strings of TPB.

• device – Backend to use; None uses the default.

Returns:

A list of packed determinants, each a list of words, sorted and distinct.

sort_bitarray(dets, device=None)

Sort determinants into canonical order, removing duplicates.

Diagonalization requires its determinant lists to be in this order. gdb_diag() sorts its own input, so this reproduces the order it works in.

from_strings(strings, bit_length, total_bit_length, device=None)

Pack many bitstrings into one (n, words) array.

The bulk form of sbd.from_string, looping in C++ instead of once per determinant across the Python boundary. Prefer it for anything larger than a handful: the per-call form costs tens of thousands of round trips on the vendored inputs alone.

Each row holds total_bit_length bits packed bit_length to a word, so a bitstring longer than bit_length occupies more than one word. The packing runs from the right end of the string: its trailing bit_length characters become word 0, the next bit_length to the left become word 1, and so on, so the leading characters of a multi-word string land in the last word rather than the first:

from_strings(["0011"], 20, 4)           # -> array([[3]], dtype=uint64)
from_strings(["1" + "0" * 69], 62, 70)  # -> array([[  0, 128]], ...)
from_strings(["0" * 69 + "1"], 62, 70)  # -> array([[1, 0]], ...)

Within a word the rightmost character is the least significant bit, which is why "0011" packs to 3 rather than 12. Getting this backwards produces a valid-looking but wrong subspace rather than an error, so check round trips with sbd.makestring when building determinants by hand.

bit_length must match the bit_length of the config passed to the diagonalizers, since it sets how many words each determinant occupies.

sort_bitarray_array(dets, device=None)

Sort packed determinants into canonical order, array in and array out.

The array form of sort_bitarray(), for a determinant list held as an (n, words) array. Deduplicates, like the list form.

Backends and devices

init(device='cpu', comm_backend='mpi')

Initialize SBD with MPI and set the default compute device.

Calling init() explicitly is optional — SBD auto-initializes on first use with device='cpu' and comm_backend='mpi'. Call it explicitly only when you need GPU or want to control startup timing.

The device can be overridden per call via the device parameter of tpb_diag(), tpb_diag_from_files(), gdb_diag() and get_backend().

Parameters:

• device – Default compute device — ‘cpu’, ‘gpu’, ‘gpu-omp’, or ‘auto’. ‘gpu’ is the NVIDIA-only Thrust backend; ‘gpu-omp’ is OpenMP target offload and serves NVIDIA and AMD alike. Aliases for ‘gpu’: ‘gpu-thrust’, ‘gpu-nvidia’, ‘cuda’. Aliases for ‘gpu-omp’: ‘gpu-omp-offload’, ‘gpu-nvhpc-omp’, ‘gpu-nvidia-omp’, ‘gpu-amd-omp’, ‘gpu-rocm-omp’, ‘rocm’.

• comm_backend – Communication backend — ‘mpi’.

Raises:

RuntimeError – If MPI is not available or no backends are compiled.

finalize()

Finalize SBD and reset session state.

Synchronizes GPU (if used) but does NOT call MPI_Finalize — mpi4py handles MPI lifecycle automatically.

After finalize(), init() can be called again.

available_backends()

Devices whose extension is built and structurally loadable.

Each candidate .so is inspected statically – present on disk, no unresolved shared-library dependencies, matching architecture – without importing it. Importing would run MPI_Init, which hangs on an MPI lacking PMIx support when the process was not started under a launcher.

A listed backend is therefore “present and structurally sound”, not “guaranteed to load”: an ABI/symbol mismatch only shows up on a real dlopen, and surfaces from get_backend(). Reasons for anything excluded are in backend_load_errors(); a backend that is merely not built is recorded there but does not warn, since building a subset is normal.

backend_load_errors()

Why each unusable backend is unusable.

Maps device -> reason, distinguishing “not built” from a missing shared library, an architecture mismatch, or a failed import.

get_backend(device=None)

Get the backend module for the given device.

Auto-initializes SBD if needed. Passing device overrides the default, which is how you switch backends within one process. The one combination to avoid is the CPU backend together with the OpenMP-offload one (‘gpu-omp’): they share an OpenMP runtime, and loading both leaves offload regions running on the host. sbd.has_backend_conflict() reports it.

Parameters:

device – ‘cpu’, ‘gpu’ (Thrust), ‘gpu-omp’ (OpenMP offload), ‘auto’, an alias accepted by init(), or None. None means the device set by init(). If nothing has initialized yet it means ‘auto’, which prefers Thrust, then OpenMP offload, then CPU, whereas the first call of tpb_diag() or gdb_diag() initializes with ‘cpu’. Call init() with the device you want to avoid depending on call order.

Returns:

The pybind11 backend module (_core_cpu, _core_gpu_thrust, or _core_gpu_omp_offload).

get_device_id(device=None)

Get the device index this rank will use, or -1 if there is none.

Reports what the backend’s own rank-to-device rule yields rather than recomputing it here, so a caller putting another library on the same card cannot drift out of step with SBD. Selects nothing and creates no context, so it is safe to call before any GPU work.

Uses the communicator rank, i.e. the same gpu_id = rank % num_gpus convention SBD’s own diagonalization applies, documented in examples/tpb/README.md.

MPI

get_rank()

Get MPI rank of current process.

get_world_size()

Get total number of MPI processes.