Tutorial: Writing a Custom Solver

This tutorial shows how to implement a custom solver by subclassing the Solver ABC and returning a SolverResult.

The Solver Interface

Every solver must implement one method:

from qvartools.solvers.solver import Solver, SolverResult
from qvartools.hamiltonians.hamiltonian import Hamiltonian

class MySolver(Solver):
    def solve(self, hamiltonian: Hamiltonian, mol_info: dict) -> SolverResult:
        ...

The mol_info dictionary contains at minimum:

  • "name": molecule name (str)

  • "n_qubits": number of qubits (int)

  • "basis": Gaussian basis set (str)

  • "geometry": atomic coordinates (list)

  • "charge": molecular charge (int)

  • "spin": 2S value (int)

Example: Random Sampling Solver

Here’s a complete example of a solver that samples random configurations, builds a projected Hamiltonian, and diagonalizes:

import time
import numpy as np
import torch
from qvartools.solvers.solver import Solver, SolverResult
from qvartools.hamiltonians.hamiltonian import Hamiltonian
from qvartools.diag.eigen.projected_hamiltonian import ProjectedHamiltonianBuilder
from qvartools.diag.eigen.eigensolver import solve_generalized_eigenvalue

class RandomSamplingSolver(Solver):
    """Solver that samples random valid configurations.

    Parameters
    ----------
    n_samples : int
        Number of random configurations to generate.
    seed : int, optional
        Random seed for reproducibility.
    """

    def __init__(self, n_samples: int = 1000, seed: int | None = None):
        self._n_samples = n_samples
        self._seed = seed

    def solve(self, hamiltonian: Hamiltonian, mol_info: dict) -> SolverResult:
        t0 = time.perf_counter()

        if self._seed is not None:
            torch.manual_seed(self._seed)

        n_sites = hamiltonian.num_sites

        # Generate random binary configurations
        configs = torch.randint(0, 2, (self._n_samples, n_sites))
        configs = torch.unique(configs, dim=0)

        # Build projected Hamiltonian
        builder = ProjectedHamiltonianBuilder(hamiltonian)
        H_proj = builder.build(configs)

        # Diagonalize
        eigenvalues, _ = solve_generalized_eigenvalue(
            H_proj, np.eye(len(H_proj)), k=1
        )

        wall_time = time.perf_counter() - t0
        return SolverResult(
            energy=float(eigenvalues[0]),
            diag_dim=len(configs),
            wall_time=wall_time,
            method="RandomSampling",
            converged=True,
            metadata={"n_samples_requested": self._n_samples},
        )

Using the Custom Solver

from qvartools.molecules import get_molecule

hamiltonian, mol_info = get_molecule("H2")

solver = RandomSamplingSolver(n_samples=500, seed=42)
result = solver.solve(hamiltonian, mol_info)
print(result)
# SolverResult(method='RandomSampling', energy=..., ...)

Design Guidelines

When writing custom solvers:

  1. Return ``SolverResult`` with all required fields populated.

  2. Use ``metadata`` for method-specific information (iteration counts, per-step energies, convergence history).

  3. Set ``converged`` to False when the solver reaches its iteration limit without meeting the tolerance.

  4. Measure ``wall_time`` from the start of solve() to the end.

  5. Handle edge cases: empty basis, singular matrices, missing PySCF, etc.

See Extending Solvers for more advanced patterns.