API reference for NuQuLib

Module contents

Nuclear Quantum computing Library (NuQuLib).

NuQuLib provides a set of tools for quantum simulations of nuclear physics problems. This library enables quantum computing research in nuclear physics by providing:

  • Nuclear Hamiltonian encodings/mappings

  • Quantum circuit implementations for some selected problems and algorithms

  • A small example of a Variational Quantum Eigensolver (VQE) using PennyLane/Qiskit

The library is designed to work in conjunction with existing nuclear physics codes and supports various interaction formats used in nuclear structure calculations.

Submodules

angular_momentum_projection submodule

class nuqulib.angular_momentum_projection.AngularMomentumProjector(hamiltonian, simulator, num_ancillas=1, angle_scale=2.0)[source]

Bases: object

Apply ancilla-assisted Jz, Jx, and iterative J=0 projections.

The projector caches the orbital blocks and the second-quantized KJx gate associated with a NuQuLib Hamiltonian. The current statevector implementation assumes that the ancillas are the least-significant qubits and that the target register immediately follows them.

Parameters:
  • hamiltonian (Hamiltonian) – NuQuLib Hamiltonian defining the ordered M-scheme orbitals.

  • simulator (qiskit.providers.Backend) – Statevector-capable simulator.

  • num_ancillas (int, optional) – Number of least-significant ancilla qubits.

  • angle_scale (float, optional) – Conversion factor between one-particle Givens angles and the parameter accepted by G_gate.

project_j_zero(circuit, ancilla_qubit, target_qubits, angles, num_steps=1, postselect_state=0, verbose=False, history_granularity='stage')[source]

Iteratively apply Jz and Jx projections to approach J=0.

The returned ProjectionResult contains a history entry after every Jz and Jx stage by default. Set history_granularity to "angle" to apply and record each filter angle separately. The angle-resolved mode implements the same ordered filter product while exposing the evolution within each Jz and Jx stage.

project_jx(circuit, ancilla_qubit, target_qubits, angles, postselect_state=0, verbose=False)[source]

Apply Jx projection using KJx_dagger, Jz filters, and KJx.

project_jz(circuit, ancilla_qubit, target_qubits, angles, postselect_state=0, verbose=False)[source]

Apply a sequence of post-selected Jz filters.

class nuqulib.angular_momentum_projection.ProjectionResult(state: Statevector, probability: float, circuit: QuantumCircuit, history: tuple[ProjectionRecord, ...] = ())[source]

Bases: object

Result of a post-selected angular-momentum projection.

circuit: QuantumCircuit
history: tuple[ProjectionRecord, ...] = ()
probability: float
state: Statevector
nuqulib.angular_momentum_projection.print_state_amplitudes(statevector, Nq, Na, filter_0=True, tolerance=0.0001)[source]

Print significant amplitudes for one selected ancilla sector.

Parameters:
  • statevector (array_like) – Statevector of the ancilla-plus-target system.

  • Nq (int) – Number of target qubits.

  • Na (int) – Number of least-significant ancilla qubits.

  • filter_0 (bool, optional) – Select ancilla state 0 when True. When False, select the all-ones ancilla state 2**Na - 1.

  • tolerance (float, optional) – Suppress amplitudes with magnitude below this threshold.

class nuqulib.angular_momentum_projection.ProjectionRecord(step: int, axis: str, state: Statevector, probability: float, cumulative_probability: float, circuit: QuantumCircuit, angle_index: int | None = None, angle: float | None = None)[source]

Bases: object

State and success probability recorded after one projection stage.

step

One-based iteration index of the alternating Jz/Jx projection.

Type:

int

axis

Projected angular-momentum component, either "Jz" or "Jx".

Type:

str

state

Normalized target state after this projection stage.

Type:

qiskit.quantum_info.Statevector

probability

Conditional post-selection probability for this stage.

Type:

float

cumulative_probability

Product of all conditional probabilities up to this stage.

Type:

float

circuit

Circuit preparing state on the target register.

Type:

qiskit.QuantumCircuit

angle_index

One-based index of the filter angle when angle-resolved history is requested. None for stage-level history.

Type:

int or None

angle

Filter angle used for this record. None for stage-level history.

Type:

float or None

ansatz submodule

Quantum ansatz circuits for nuclear simulations.

This module provides various ansatz implementations for nuclear quantum simulations, including Hartree-Fock states (lowest-filling more precisely), Givens rotation-based circuits, and pair Unitary Coupled-Cluster Doubles (pUCCD).

nuqulib.ansatz.circuit_XXYY(qc_ansatz: QuantumCircuit, adopted: str, Nq: int, methods_XXYY: str = 'Google', backend: str | None = None, opt_level=3, verbose: bool = False)[source]

Generates quantum circuits to measure X_iX_j+Y_iY_j terms in the Hamiltonian.

Parameters:
  • qc_ansatz (QuantumCircuit) – Ansatz circuit to copy before adding basis rotations and measurements.

  • adopted (str) – String that indicates the simulator/device type. Acceptable values are “simNISQ”/”simFTQC” and “Real”.

  • Nq (int) – Number of qubits in the ansatz circuit.

  • methods_XXYY (str, optional) – Circuit construction method. Currently only "Google" is supported.

  • backend (optional) – Backend used for transpilation when targeting a real device or noisy simulation.

  • opt_level (int, optional) – Qiskit transpiler optimization level.

  • verbose (bool, optional) – If True, also return intermediate circuits before measurement for inspection.

Returns:

A list of transpiled and/or decomposed quantum circuits with measurements added. If verbose is True, a tuple is returned where the first element is the list of final circuits and the second element is the list of intermediate circuits prepared for drawing.

Return type:

list

Raises:

ValueError – If the provided methods_XXYY is not “Google”, which means that this function is designed for basis rotation (Givens rotation) to diagonalize XX+YY term in the computational basis.

nuqulib.ansatz.lowest_filling_ansatz(Hamiltonian_obj: Hamiltonian, proton_number: int, neutron_number: int, Mtot: int | None = None, return_idxs: bool = False)[source]

Constructs a lowest-filling ansatz circuit for protons and neutrons.

This function creates a quantum circuit by first preferring zero-seniority (jz, -jz) pair fillings and then occupying additional orbitals to satisfy the requested total magnetic quantum number Mtot. The ordering of single-particle states (msps_p and msps_n) can vary with a specific choice of the interaction.

Parameters:
  • Hamiltonian_obj (Hamiltonian) – The Hamiltonian object containing system information.

  • proton_number (int) – Number of protons to be placed in the lowest energy states.

  • neutron_number (int) – Number of neutrons to be placed in the lowest energy states.

  • Mtot – int | None: Total magnetic quantum number (optional).

  • return_idxs (bool) – If True, also return the selected proton and neutron orbital indices.

Returns:

A quantum circuit representing the lowest-filling ansatz. If return_idxs is True, returns (circuit, selected_p, selected_n).

Return type:

QuantumCircuit or tuple

nuqulib.ansatz.naive_filling_ansatz(proton_qubits: Iterable[int], neutron_qubits: Iterable[int], proton_number: int, neutron_number: int, mapping_method='JordanWigner', Hamildict_opform: dict = {}, filepath: str | PathLike = './')[source]

Construct a product-state ansatz by filling the last proton/neutron orbitals.

Parameters:
  • proton_qubits (Iterable[int]) – Qubit indices assigned to proton modes.

  • neutron_qubits (Iterable[int]) – Qubit indices assigned to neutron modes.

  • proton_number (int) – Number of proton orbitals to occupy.

  • neutron_number (int) – Number of neutron orbitals to occupy.

  • mapping_method (str, optional) – Fermion-to-qubit mapping used to create the occupation-flip Pauli strings. Defaults to "JordanWigner".

  • Hamildict_opform (dict, optional) – Operator-form Hamiltonian dictionary containing "SPE" entries for proton and neutron sectors.

  • filepath (str | os.PathLike, optional) – Base path used by mappings that need to load or save mapper data.

Returns:

Circuit preparing the naive filling configuration.

Return type:

QuantumCircuit

nuqulib.ansatz.nucl_ansatz(Hamildict_opform: dict, n_qubit: int, proton_qubits: Iterable[int], neutron_qubits: Iterable[int], proton_number: int, neutron_number: int, params: Iterable[float], method: str = 'HF', mapping_method: str = 'JordanWigner', filepath: str | PathLike = './', return_Gdict: bool = False)[source]

Construct nuclear ansatz circuit for proton-neutron systems.

Creates quantum circuits for nuclear many-body states with separate proton and neutron sectors. Supports Hartree-Fock initial states and Givens rotation-based variational ansätze.

Parameters:
  • Hamildict_opform (dict) – Operator-form Hamiltonian dictionary used by the filling-state preparation.

  • n_qubit (int) – Total number of qubits.

  • proton_qubits (Iterable[int]) – Indices of proton qubits.

  • neutron_qubits (Iterable[int]) – Indices of neutron qubits.

  • proton_number (int) – Number of protons.

  • neutron_number (int) – Number of neutrons.

  • params (Iterable[float]) – Variational parameters.

  • method (str, optional) – Ansatz method. Options: “HF”, “HF+Givens”. Defaults to “HF”.

  • mapping_method (str, optional) – Fermion-to-qubit mapping used for the initial filling state. Defaults to "JordanWigner".

  • filepath (str | os.PathLike, optional) – Base path used by mappings that need to load or save mapper data.

  • return_Gdict (bool, optional) – Whether to return gate type dictionary. Defaults to False.

Returns:

Quantum circuit representing the ansatz.

If return_Gdict=True, returns (circuit, gate_dict).

Return type:

QuantumCircuit or tuple

nuqulib.ansatz.pair_ansatz_pennylane(Hamil_pl, params, Nq: int, Nocc: int, type_of_ansatz: str = 'HF', observable: str = 'Hamil', return_Gdict: str = False, random_init_circ: str = False, idxs_hole_in=[], combination_h_v=[])[source]

Construct pairing model ansatz circuit using PennyLane.

This function is a counterpart of pair_ansatz_qiskit but uses PennyLane and only supports pUCCD ansatze.

nuqulib.ansatz.pair_ansatz_qiskit(params: Iterable[float], Norb: int, Nocc: int, method: str = 'HF', return_Gdict: bool = False, decent_order: bool = True, rotation_XXYY=[], idxs_hole_in=[])[source]

Construct pairing model ansatz circuit using Qiskit.

Creates quantum circuits for pairing Hamiltonian simulations using Hartree-Fock states with optional Givens rotation enhancements. Designed for pairing Hamiltonian or pair-wise form of shell model, such as Hard-core boson formulation focused on seniority-zero states.

Parameters:
  • params (Iterable[float]) – Variational parameters for gates, mainly Givens rotations.

  • Norb (int) – Number of orbitals (qubits).

  • Nocc (int) – Number of occupied orbitals.

  • method (str, optional) – Ansatz method. Options: “HF”, “HF+Givens”. Defaults to “HF”.

  • return_Gdict (bool, optional) – Whether to return gate type dictionary. Defaults to False.

  • decent_order (bool, optional) – Whether to use descending qubit ordering. Defaults to True.

  • rotation_XXYY (list, optional) – List of XX+YY rotation parameters. Defaults to []. This is used to diagonalize XX+YY term in the computational basis.

  • idxs_hole_in (list, optional) – Indices of hole states. Defaults to [].

Returns:

Pairing ansatz circuit. If return_Gdict=True,

returns (circuit, gate_dict).

Return type:

QuantumCircuit or tuple

Note

The circuit starts with Hartree-Fock state preparation (X gates on occupied orbitals) followed by layer of Givens rotations, which can be ragarded as excitations from the reference state.

circuits submodule

Quantum circuit implementations for nuclear simulations.

This module provides various quantum circuit implementations including Givens rotations, controlled gates, and utility functions for Pauli string measurements.

nuqulib.circuits.G2_gate(theta: float, method: str = 'Xanadu')[source]

Create a two-pair Givens rotation gate.

Parameters:
  • theta (float) – Rotation angle parameter.

  • method (str, optional) – Implementation method. Currently only “Xanadu” is supported. Defaults to “Xanadu”.

Returns:

Two-pair Givens rotation gate acting on 4 qubits.

Return type:

Gate

Raises:

ValueError – If invalid method is provided.

nuqulib.circuits.G_gate(theta: float, method: str = 'magic')[source]

Create a Givens rotation gate using specified method.

Parameters:
  • theta (float) – Rotation angle parameter.

  • method (str, optional) – Implementation method. Options: “magic”, “Xanadu”, “iSWAP_Rz”. Defaults to “magic”.

Returns:

Givens rotation gate.

Return type:

Gate

Raises:

ValueError – If invalid method is provided.

nuqulib.circuits.Givens_2_Xanadu(circ: QuantumCircuit, i: int, j: int, k: int, l: int, theta: float)[source]

Two-qubit pair Givens rotation using Xanadu method.

Implements a four-qubit Givens rotation that acts on two qubit pairs simultaneously. This is useful for efficient fermionic simulations.

Parameters:
  • circ (QuantumCircuit) – Quantum circuit to modify.

  • i (int) – First qubit of first pair.

  • j (int) – Second qubit of first pair.

  • k (int) – First qubit of second pair.

  • l (int) – Second qubit of second pair.

  • theta (float) – Rotation angle parameter.

nuqulib.circuits.Givens_Xanadu(circ: QuantumCircuit, i: int, j: int, theta: float)[source]

Givens rotation implementation from Xanadu.

Implementation based on the approach described in Quantum 6, 742 (2022) - doi:10.22331/q-2022-06-20-742.

Parameters:
  • circ (QuantumCircuit) – Quantum circuit to modify.

  • i (int) – First qubit index.

  • j (int) – Second qubit index.

  • theta (float) – Rotation angle parameter.

nuqulib.circuits.Givens_iSWAP_Rz(circ: QuantumCircuit, i: int, j: int, theta: float)[source]

Givens rotation implementation using iSWAP and Rz gates.

This implementation is based on Google’s approach from Science 369, 1084 (2020) - DOI: 10.1126/science.abb9811.

Parameters:
  • circ (QuantumCircuit) – Quantum circuit to modify.

  • i (int) – First qubit index.

  • j (int) – Second qubit index.

  • theta (float) – Rotation angle parameter.

nuqulib.circuits.Givens_magic(circ: QuantumCircuit, i: int, j: int, theta: float)[source]

Optimal Givens rotation using magic basis decomposition.

This implementation is based on optimal quantum circuits for general two-qubit gates from Phys. Rev. A 69, 032315 (2004). It minimizes the number of CNOT gates required.

Parameters:
  • circ (QuantumCircuit) – Quantum circuit to modify.

  • i (int) – First qubit index.

  • j (int) – Second qubit index.

  • theta (float) – Rotation angle parameter.

Reference:

Used in npj Quantum Information (2023) 9:60

nuqulib.circuits.additional_qc(qc_in, pauli_str, register_target, Qiskit_order=True)[source]

Add basis rotation gates for Pauli string measurement.

Applies the necessary basis rotation gates (H for X measurement, SdgH for Y measurement) to measure a Pauli string.

Parameters:
  • qc_in (QuantumCircuit) – Quantum circuit to modify.

  • pauli_str (str) – Pauli string to measure.

  • register_target (list) – List of target qubit indices.

  • Qiskit_order (bool, optional) – Whether to use Qiskit qubit ordering. Defaults to True.

Raises:

ValueError – If invalid Pauli string character is encountered.

nuqulib.circuits.apply_sqrt_iSWAP(qc: QuantumCircuit, i: int, j: int)[source]

Apply square root of iSWAP gate to quantum circuit.

Parameters:
  • qc (QuantumCircuit) – Quantum circuit to modify.

  • i (int) – First qubit index.

  • j (int) – Second qubit index.

nuqulib.circuits.cG1(circ: QuantumCircuit, c_qubit: int, i: int, j: int, theta: float)[source]

Apply controlled Givens rotation to a quantum circuit.

Implements a controlled version of the Givens rotation gate using elementary quantum gates. This is used for e.g. state preparation.

Parameters:
  • circ (QuantumCircuit) – Quantum circuit to modify.

  • c_qubit (int) – Control qubit index.

  • i (int) – First target qubit index.

  • j (int) – Second target qubit index.

  • theta (float) – Rotation angle parameter.

nuqulib.circuits.cG1_gate(theta: float, method: str = 'magic')[source]

Create a controlled Givens rotation gate.

Parameters:
  • theta (float) – Rotation angle parameter.

  • method (str, optional) – Implementation method. Options: “magic”, “Xanadu”, “iSWAP_Rz”. Defaults to “magic”.

Returns:

Controlled Givens rotation gate.

Return type:

Gate

Raises:

ValueError – If invalid method is provided.

nuqulib.circuits.expec_Zstring(res: dict, idx_relevant: list[int], Qiskit_ordering: bool = True, target_qubits: list[int] = [], ancilla_qubit: int | None = None)[source]

Calculate expectation value of Z-string measurement from results.

Computes the expectation value of a Z-string Pauli operator from measurement results, with optional ancilla qubit for some algorithms such as QKrylov.

Parameters:
  • res (dict) – Dictionary of measurement results {bitstring: count/weights}. values can be either raw counts (int) or normalized weights (float)

  • idx_relevant (list) – List of relevant qubit indices for Z measurements.

  • Qiskit_ordering (bool, optional) – Whether to use Qiskit ordering. Defaults to True.

  • target_qubits (list, optional) – List of target qubit indices. Defaults to [].

  • ancilla_qubit (int, optional) – Ancilla qubit index for post-selection. Defaults to None.

Returns:

Tuple containing:
  • exp_val (float): Overall expectation value.

  • exp_val_p0 (float): Expectation value when ancilla=0.

  • exp_val_p1 (float): Expectation value when ancilla=1.

Return type:

tuple

nuqulib.circuits.get_idx_ancilla_in_string(n_qubit: int, ancilla: int | None, Qiskit_ordering: bool) int[source]

Get ancilla qubit index in measurement string.

Parameters:
  • n_qubit (int) – Total number of qubits.

  • ancilla (int or None) – Ancilla qubit index.

  • Qiskit_ordering (bool) – Whether to use Qiskit qubit ordering.

Returns:

Index of ancilla qubit in the measurement string,

or None if no ancilla qubit.

Return type:

int or None

nuqulib.circuits.get_idx_to_measure(pauli_str: str, Qiskit_order: bool = True) list[int][source]

Get qubit indices that need to be measured for a Pauli string.

Returns indices of qubits that have non-identity Pauli operators and therefore need to be measured.

Parameters:
  • pauli_str (str) – Pauli string.

  • Qiskit_order (bool, optional) – Whether to use Qiskit qubit ordering. Defaults to True.

Returns:

List of qubit indices to measure.

Return type:

list

diagonalization submodule

Functions to diagonalize nuclear Hamiltonians

nuqulib.diagonalization.Diagonalize_Hamiltonian(Hamil_mapped: SparsePauliOp, hamil: Hamiltonian, Z: int, N: int, target_parity: int, Zc: int = 0, Nc: int = 0, calc_J2: bool = True, verbose: bool = False, use_basis: str = 'NPM') dict[source]

Diagonalize a mapped Hamiltonian in a particle-number projected basis.

Parameters:
  • Hamil_mapped (SparsePauliOp) – Qubit Hamiltonian to diagonalize.

  • hamil (Hamiltonian) – Nuclear Hamiltonian object defining the model space and proton/neutron split.

  • Z (int) – Proton number of the target nucleus.

  • N (int) – Neutron number of the target nucleus.

  • target_parity (int) – Target parity, either +1 or -1.

  • Zc (int, optional) – Number of inert core protons.

  • Nc (int, optional) – Number of inert core neutrons.

  • calc_J2 (bool, optional) – If True, also compute J^2 expectation values in the eigenvectors.

  • verbose (bool, optional) – If True, print basis and diagonalization details.

  • use_basis (str, optional) – Basis projection mode. "NPM" fixes particle number, parity, and M; "NP" fixes particle number and parity; any other value uses the full computational basis.

Returns:

Diagonalization results, including basis, projected Hamiltonian, eigenvalues, eigenvectors, and optionally J^2 data.

Return type:

dict

nuqulib.diagonalization.Mprojected_hamiltonian(hamiltonian: SparsePauliOp, n_qubits_p: int, n_qubits_n: int, basis: list[int], index: dict[int, int]) ndarray[source]

Build the M-projected Hamiltonian for a Hamiltonian whose Pauli strings may act on protons, neutrons, or both. Treat proton and neutron parts independently by padding shorter Pauli labels appropriately.

nuqulib.diagonalization.Mprojected_hamiltonian_single(hamiltonian, n_qubits, basis, index)[source]

Project a single-species Pauli Hamiltonian onto a fixed basis.

Parameters:
  • hamiltonian (SparsePauliOp) – Pauli Hamiltonian acting on one species.

  • n_qubits (int) – Number of qubits in the single-species register.

  • basis (list[int]) – Basis states retained in the projected subspace.

  • index (dict[int, int]) – Mapping from basis state to matrix index.

Returns:

Dense Hermitian Hamiltonian matrix in the projected basis.

Return type:

np.ndarray

nuqulib.diagonalization.angular_momentum_from_J2(J2_values: ndarray, tol: float = 1e-10) ndarray[source]

Convert J(J+1) values to the corresponding effective J.

nuqulib.diagonalization.apply_pauli_string_pn(pauli: str, state: int, n_qubits_p: int, n_qubits_n: int, debug=False)[source]

Apply a proton-neutron Pauli string to an integer-encoded basis state.

Parameters:
  • pauli (str) – Pauli label with neutron part followed by proton part.

  • state (int) – Integer-encoded combined proton-neutron state.

  • n_qubits_p (int) – Number of proton qubits.

  • n_qubits_n (int) – Number of neutron qubits.

  • debug (bool, optional) – If True, print the state transition and phase.

Returns:

New combined state and accumulated phase.

Return type:

tuple[int, complex]

nuqulib.diagonalization.apply_pauli_string_single(pauli: str, state: int, n_qubits: int) tuple[int, complex][source]

Apply one Pauli string to an integer-encoded basis state.

Parameters:
  • pauli (str) – Pauli label ordered as a Qiskit label.

  • state (int) – Integer-encoded computational basis state.

  • n_qubits (int) – Number of qubits represented by state.

Returns:

New basis state and accumulated phase.

Return type:

tuple[int, complex]

nuqulib.diagonalization.expectation_values(operator_matrix: ndarray, states: ndarray) ndarray[source]

Return <psi|O|psi> for one state vector or column-wise state vectors.

nuqulib.diagonalization.fixed_N_P_M_basis(n_qubits_p: int, n_qubits_n: int, Z: int, N: int, msps: list, parity: int, M_tot: int, verbose: bool = False) tuple[list[int], dict[int, int]][source]

Generate basis states with fixed proton number Z, neutron number N, parity, and total M. Note that M is assumed to be doubled.

nuqulib.diagonalization.fixed_N_P_M_basis_neutron(n_qubits, n_particles, msps, parity, M_tot)[source]

Generate neutron-only basis states with fixed N, parity, and doubled M.

Parameters:
  • n_qubits (int) – Number of neutron qubits.

  • n_particles (int) – Required number of occupied neutron orbitals.

  • msps (list) – M-scheme single-particle states.

  • parity (int) – Target parity, either +1 or -1.

  • M_tot (int) – Target total M in doubled convention.

Returns:

Basis states and state-to-row index.

Return type:

tuple[list[int], dict[int, int]]

nuqulib.diagonalization.fixed_N_P_basis(n_qubits_p: int, n_qubits_n: int, Z: int, N: int, target_parity: int, msps: list) tuple[list[int], dict[int, int]][source]

Generate basis states with fixed proton number Z and neutron number N, without parity or M projection.

nuqulib.diagonalization.selected_ci_sequential(basis_configs, Hmat, sampled_counts, max_iter=100, add_per_iter=15, threshold=None, initial_size=None)[source]

Selected CI: sequentially grow a subspace by adding configurations with the largest coupling to the current subspace.

Parameters:
  • basis_configs (list[int]) – M-scheme basis (bitstring as int).

  • Hmat (np.ndarray) – Full Hamiltonian matrix in basis_configs ordering.

  • sampled_counts (dict[str,int] | list[str] | list[int]) – Sampled configurations. If dict, keys are bitstrings (“0101…”) and counts are their weights. If list, it can be bitstrings or indices.

  • max_iter (int) – Number of selection iterations.

  • add_per_iter (int | None) – How many configs to add per iteration. None -> add all candidates.

  • threshold (float | None) – Threshold on max |H_ij|. If set, add configs with coupling >= threshold.

  • initial_size (int | None) – If sampled_counts is dict, use only top-N by count as initial pool.

Returns:

  • selected_indices (list[int]) – Selected subspace indices in basis_configs order (sequentially grown).

  • history (list[dict]) – Iteration log with added counts and selected size.

nuqulib.diagonalization.total_angular_momentum_squared_matrix(n_qubits_p: int, n_qubits_n: int, basis: list[int], index: dict[int, int], msps: list, M_tot: int | None = None) ndarray[source]

Build the J^2 matrix in the fixed-M M-scheme basis.

The single-particle angular momenta in msps are assumed to use the doubled convention: j=2J and jz=2M. The returned matrix is in physical units, so an eigenstate with angular momentum J has eigenvalue J(J+1).

encoding submodule

Quantum encoding utilities for nuclear Hamiltonians.

This module provides functions for mapping nuclear Hamiltonians to qubit operators using e.g., Jordan-Wigner, and utilities for working with Pauli operators in different quantum computing frameworks.

nuqulib.encoding.check_XXYYterm(hamiltonian_op_XXYY)[source]

Check that XX and YY terms have identical coefficients.

In pairing or pair-wise Hamiltonians, the XX and YY terms should have the same coefficient. This function validates this constraint for debugging and verification.

Parameters:

hamiltonian_op_XXYY (SparsePauliOp) – Hamiltonian containing only XX and YY terms.

Returns:

True if all XX/YY coefficient pairs match.

Return type:

bool

Raises:

AssertionError – If XX and YY terms have different coefficients.

nuqulib.encoding.mapping_of_pn_hamiltonians(op_pn: dict[tuple[str, str], float], n_qubits_p: int, n_qubits_n: int, method: str, Hamildict_specified_p: dict, Hamildict_specified_n: dict, filepath_p: str | PathLike, filepath_n: str | PathLike)[source]

Map proton-neutron coupled Hamiltonians to Pauli operators.

This function handles the mapping of nuclear Hamiltonians that include both proton and neutron sectors with their respective interactions. We here assume that proton indices are lower than neutron indices, and those are to be coupled like neutron part followed by proton part so that one can use them in Qiskit.

Parameters:
  • op_pn (dict) – Dictionary with (proton_str, neutron_str) keys and coefficient values.

  • n_qubits_p (int) – Number of qubits for proton sector.

  • n_qubits_n (int) – Number of qubits for neutron sector.

  • method (str) – Encoding method (e.g., “Jordan-Wigner”, “Bravyi-Kitaev”).

  • Hamildict_specified_p (dict) – Dictionary of Hamiltonian terms in fermionic form for protons.

  • Hamildict_specified_n (dict) – Dictionary of Hamiltonian terms in fermionic form for neutrons.

  • filepath_p (str|os.PathLike) – File path for saving/loading proton mapper.

  • filepath_n (str|os.PathLike) – File path for saving/loading neutron mapper.

Returns:

Combined Pauli operator representing the full Hamiltonian.

Return type:

SparsePauliOp

nuqulib.encoding.mapping_to_Pauli_string(Fermionic_op: FermionicOp, n_qubits: int, init_qubit: int = 0, method: str = 'JordanWigner', Hamildict_specified: dict = {}, filepath: str | PathLike = './', verbose: bool = False)[source]

Map fermionic operators to Pauli strings using specified encoding.

Parameters:
  • Fermionic_op (FermionicOp) – Fermionic operator to be mapped.

  • n_qubits (int) – Number of qubits for the mapping.

  • init_qubit (int) – Starting index for qubits in the mapping.

  • method (str) – Encoding method. Options: “JordanWigner”/”JW”/”Jordan-Wigner” or “BravyiKitaev”/”BK”/”Bravyi-Kitaev”.

  • Hamildict_specified (dict) – Dictionary of Hamiltonian terms in fermionic form. used in the special case, HATTMapper.

  • filepath (str|os.PathLike) – File path for saving/loading mapper in the case of HATTMapper.

  • verbose (bool) – If True, print intermediate mapping results for debugging.

Returns:

Mapped Pauli operator.

Return type:

SparsePauliOp

Raises:

ValueError – If invalid encoding method is provided.

Note

BravyiKitaev mapping has never been tested in this module.

nuqulib.encoding.qpo_from_sparsepauliop(sp_op: SparsePauliOp) QubitPauliOperator[source]

Convert Qiskit SparsePauliOp to PyTKET QubitPauliOperator.

This function converts Qiskit’s SparsePauliOp representation to PyTKET’s QubitPauliOperator format, with automatic reversal of qubit ordering to match PyTKET conventions.

Parameters:

sp_op (SparsePauliOp) – Qiskit SparsePauliOp to convert.

Returns:

PyTKET QubitPauliOperator with reversed qubit order.

Return type:

QubitPauliOperator

Note

This code is based on PyTKET documentation but modified to handle qubit ordering differences between Qiskit and PyTKET. Original reference

nuqulib.encoding.qps_from_sparsepauliop(paulis)[source]

Convert SparsePauliOp Pauli strings to PyTKET QubitPauliString.

Parameters:

paulis (str) – Pauli string (e.g., “IXYZ”).

Returns:

PyTKET QubitPauliString representation.

Return type:

QubitPauliString

Note

Identity operators are automatically filtered out as they don’t contribute to the QubitPauliString representation.

nuqulib.encoding.separate_Hamil_terms(hamiltonian_op: SparsePauliOp)[source]

Separate Hamiltonian, pairing or pair-wise ones, into diagonal and XX+YY terms.

This function partitions a Hamiltonian operator into terms that are diagonal in the computational basis (I and Z terms) and off-diagonal XX+YY terms. This separation is useful for quantum algorithms that treat different types of terms separately.

Parameters:

hamiltonian_op (SparsePauliOp) – Full Hamiltonian operator.

Returns:

Tuple containing:
  • hamiltonian_op_diag (SparsePauliOp): Diagonal terms (I, Z).

  • hamiltonian_op_XXYY (SparsePauliOp): XX and YY terms.

Return type:

tuple

Note

This function validates that XX and YY terms have matching coefficients and assumes real coefficients.

nuqulib.encoding.task_pn_mapping(op_pn_key, op_pn, n_qubits_p, n_qubits_n, method, Hamildict_specified_p, Hamildict_specified_n, filepath_p, filepath_n)[source]

Map one proton-neutron fermionic term to Pauli-label terms.

Parameters:
  • op_pn_key (tuple[str, str]) – (proton_op, neutron_op) fermionic operator strings.

  • op_pn (dict) – Dictionary containing the coefficient for op_pn_key.

  • n_qubits_p (int) – Number of proton qubits.

  • n_qubits_n (int) – Number of neutron qubits.

  • method (str) – Fermion-to-qubit mapping method.

  • Hamildict_specified_p (dict) – Proton-sector Hamiltonian terms used by mappings such as HATTMapper.

  • Hamildict_specified_n (dict) – Neutron-sector Hamiltonian terms used by mappings such as HATTMapper.

  • filepath_p (str | os.PathLike) – Proton-sector mapper cache path.

  • filepath_n (str | os.PathLike) – Neutron-sector mapper cache path.

Returns:

Pauli labels and coefficients generated for this term.

Return type:

list[list]

hatt_mapper submodule

class nuqulib.hatt_mapper.HATTMapper(loader: FermionicOp | list[str], nqubits: int | None = None)[source]

Bases: ModeBasedMapper, FermionicMapper

Heuristic anticommutation tree mapper for fermionic operators.

The mapper can be built directly from a FermionicOp or from a saved list of Pauli-table strings. The resulting table is then reused when mapping additional operators with the same register size.

static load(path: str)[source]

Load a HATT mapper from a saved raw Pauli-table file.

Parameters:

path (str) – Input filename created by save().

Returns:

Mapper initialized with the stored table.

Return type:

HATTMapper

map(second_q_ops: FermionicOp, *, _: int | None = None) SparsePauliOp[source]

Map a fermionic operator using the precomputed HATT Pauli table.

Parameters:
  • second_q_ops (FermionicOp) – Fermionic operator to map.

  • _ (int | None, optional) – Ignored compatibility argument.

Returns:

Qubit operator in Qiskit’s Pauli representation.

Return type:

SparsePauliOp

pauli_table(register_length: int)[source]

Return the Pauli table expected by Qiskit Nature mappers.

Parameters:

register_length (int) – Register length requested by the mapper API. The stored HATT table determines the actual Pauli strings.

Returns:

Majorana-pair Pauli table.

Return type:

list[tuple[Pauli, Pauli]]

save(path: str)[source]

Write the raw Pauli table to a text file.

Parameters:

path (str) – Output filename.

majoranaop submodule

The Majorana-particle Operator.

class nuqulib.majoranaop.MajoranaOp(data: Mapping[str, complex | ParameterExpression], num_modes: int | None = None, *, copy: bool = True, validate: bool = True)[source]

Bases: SparseLabelOp

N-mode Majorana operator.

A MajoranaOp represents a weighted sum of Majorana fermion operator terms. These terms are encoded as sparse labels, which are strings consisting of a space-separated list of expressions. Each expression must look like _<index>, where the <index> is a non-negative integer representing the index of the mode on which the Majorana operator is applied. The maximum value of index is bound by num_modes. Note that, when converting from a FermionicOp there are two modes per spin orbital, i.e. num_modes is 2 * FermionicOp.num_spin_orbitals - 1

Initialization

A MajoranaOp is initialized with a dictionary, mapping terms to their respective coefficients:

from qiskit_nature.second_q.operators import MajoranaOp

op = MajoranaOp(
    {
        "_0 _1": .25j,
        "_1 _0": -.25j,
        "_2 _3": -.25j,
        "_3 _2": .25j,
    },
    num_modes=4,
)

By default, this way of initializing will create a full copy of the dictionary of coefficients. If you have very restricted memory resources available, or would like to avoid the additional copy, the dictionary will be stored by reference if you disable copy like so:

some_big_data = {
    "_0 _1": .25j,
    "_1 _0": -.25j,
    # ...
}

op = MajoranaOp(
    some_big_data,
    num_modes=4,
    copy=False,
)

Note

It is the users’ responsibility, that in the above scenario, some_big_data is not changed after initialization of the MajoranaOp, since the operator contents are not guaranteed to remain unaffected by such changes.

Construction from Fermionic operator

As an alternative to the manual construction above, a more convenient way of initializing a MajoranaOp is, to construct it from an existing FermionicOp:

from qiskit_nature.second_q.operators import FermionicOp, MajoranaOp
f_op = FermionicOp({"+_0 -_1": 1}, num_spin_orbitals=2)
m_op = MajoranaOp.from_fermionic_op(f_op)

Note that each FerminonicOp-term consisting of \(n\) expressions will result in a MajoranaOp-term consisting of \(2^n\) expressions. The conversion uses the convention that

\[a_i = \frac{1}{2}(\gamma_{2i} + i \gamma_{2i+1}), \quad a_i^\dagger = \frac{1}{2}(\gamma_{2i} - i \gamma_{2i+1}) \,,\]

where \(a_i\) and \(a_i^\dagger\) are the Fermionic annihilation and creation operators and \(\gamma_i\) the Majorana operators.

Construction from a ``PolynomialTensor``

Using the from_polynomial_tensor() constructor method, a MajoranaOp can be constructed from a PolynomialTensor. In this case, the underscore character _ is the only allowed character in the keys of the PolynomialTensor. For example,

p_t = PolynomialTensor(
    {
        "_": np.arange(1, 3),
        "__": np.arange(1, 5).reshape((2, 2)),
    }
)
op = MajoranaOp.from_polynomial_tensor(p_t)

# op is then
MajoranaOp({'_0': 1, '_1': 2, '_0 _0': 1, '_0 _1': 2, '_1 _0': 3, '_1 _1': 4}, num_modes=2)

Algebra

This class supports the following basic arithmetic operations: addition, subtraction, scalar multiplication, operator multiplication, and adjoint. For example,

Addition

MajoranaOp({"_1": 1}, num_modes=2) + MajoranaOp({"_0": 1}, num_modes=2)

Sum

sum(MajoranaOp({label: 1}, num_modes=4) for label in ["_0", "_1", "_2 _3"])

Scalar multiplication

0.5 * MajoranaOp({"_1": 1}, num_modes=2)

Operator multiplication

op1 = MajoranaOp({"_0 _1": 1}, num_modes=3)
op2 = MajoranaOp({"_0 _1 _2": 1}, num_modes=3)
print(op1 @ op2)

Tensor multiplication

op = MajoranaOp({"_0 _1": 1}, num_modes=2)
print(op ^ op)

Adjoint

MajoranaOp({"_0 _1": 1j}, num_modes=2).adjoint()

Note

Since Majorana operators are self-adjoined, the adjoint of a MajoranaOp is the original operator with all strings reversed, e.g. "_0 _1" becomes "_1 _0" in the example above, and coefficients become complex conjugated.

Iteration

Instances of MajoranaOp are iterable. Iterating a MajoranaOp yields (term, coefficient) pairs describing the terms contained in the operator.

num_modes

the number of modes on which this operator acts. This is considered a lower bound, which means that mathematical operations acting on two or more operators will result in a new operator with the maximum number of modes of any of the involved operators. When converting from a FermionicOp, this is twice the number of spin orbitals.

Type:

int | None

Note

MajoranaOp can contain qiskit.circuit.ParameterExpression objects as coefficients. However, a MajoranaOp containing parameters does not support the following methods:

  • is_hermitian

compose(other: MajoranaOp, qargs=None, front: bool = False) MajoranaOp[source]

Compose this operator with another Majorana operator.

Parameters:
  • other (MajoranaOp) – Operator to compose with this one.

  • qargs – Included for compatibility with the Qiskit operator API.

  • front (bool, optional) – If True, place other in front of this operator in the composition.

Returns:

Composed operator.

Return type:

MajoranaOp

expand(other: MajoranaOp) MajoranaOp[source]

Return the reverse-order tensor product with other.

classmethod from_fermionic_op(op: FermionicOp, *, simplify: bool = True) MajoranaOp[source]

Constructs the operator from a FermionicOp.

Parameters:
  • op – the FermionicOp to convert.

  • simplify – whether to index order and simplify the resulting operator.

Returns:

The converted MajoranaOp.

classmethod from_polynomial_tensor(tensor: PolynomialTensor) MajoranaOp[source]

Construct a Majorana operator from a polynomial tensor.

Parameters:

tensor (PolynomialTensor) – Tensor whose keys contain only underscores and can be expanded into Majorana labels.

Returns:

Operator containing the tensor entries as sparse labels.

Return type:

MajoranaOp

classmethod from_terms(terms: Sequence[tuple[list[tuple[str, int]], complex | ParameterExpression]]) MajoranaOp[source]

Construct an operator from split label terms.

Parameters:

terms (Sequence) – Sequence of (label_terms, coefficient) pairs, where each label term contains a mode index.

Returns:

Operator assembled from the provided terms.

Return type:

MajoranaOp

index_order() MajoranaOp[source]

Convert to the equivalent operator with the terms of each label ordered by index.

Returns a new operator (the original operator is not modified).

Note

You can use this method to achieve the most aggressive simplification. simplify() does not reorder the terms. For instance, using only simplify() will reduce _2 _0 _1 _0 _0 to _2 _0 _1 but cannot deduce this label to be identical to _0 _1 _2. Calling this method will reorder the former label to _0 _0 _0 _1 _2, after which simplify() will be able to correctly collapse these two labels into one.

Returns:

The index ordered operator.

is_hermitian(atol: float | None = None) bool[source]

Checks whether the operator is hermitian.

Parameters:

atol – Absolute numerical tolerance. The default behavior is to use self.atol.

Returns:

True if the operator is hermitian up to numerical tolerance, False otherwise.

Raises:

ValueError – Operator contains parameters.

property register_length: int

Return the number of Majorana modes required by this operator.

simplify(atol: float | None = None) MajoranaOp[source]

Collapse duplicate Majorana labels and remove near-zero coefficients.

Parameters:

atol (float | None, optional) – Absolute tolerance for removing numerically zero coefficients. Defaults to self.atol.

Returns:

Simplified operator.

Return type:

MajoranaOp

tensor(other: MajoranaOp) MajoranaOp[source]

Return the tensor product of this operator with other.

terms() Iterator[tuple[list[tuple[str, int]], complex | ParameterExpression]][source]

Provides an iterator analogous to items() but with the labels already split into pairs of operation characters and indices.

Yields:

A tuple with two items; the first one being a list of pairs of the form (‘’, int) where the empty string is for compatibility with other SparseLabelOp and the integer corresponds to the mode index on which the operator gets applied; the second item of the returned tuple is the coefficient of this term.

transpose() MajoranaOp[source]

Return a copy with every Majorana label string reversed.

mode_based_mapper submodule

Mode Based Mapper.

class nuqulib.mode_based_mapper.ModeBasedMapper[source]

Bases: QubitMapper

Mapper from SparseLabelOp to a qubit operator using a Pauli table.

mode_based_mapping(second_q_op: SparseLabelOp, register_length: int | None = None) SparsePauliOp[source]

Utility method to map a SparseLabelOp to a qubit operator using a pauli table.

Parameters:
  • second_q_op – the SparseLabelOp to be mapped.

  • register_length – when provided, this will be used to overwrite the register_length attribute of the operator being mapped. This is possible because the register_length is considered a lower bound.

Returns:

The qubit operator corresponding to the problem-Hamiltonian in the qubit space.

Raises:

QiskitNatureError – If number length of pauli table does not match the number of operator modes, or if the operator has unexpected label content

abstractmethod pauli_table(register_length: int) list[tuple[PauliList | SparsePauliOp | Pauli | list | str, PauliList | SparsePauliOp | Pauli | list | str]][source]

Generates a Pauli-lookup table mapping from modes to Pauli operators or pairs of Pauli operators.

This table is a list of tuples \((P, Q)\) of two Pauli operators, corresponding to the real part \(P\) and imaginary part \(Q\) for the respective mode index. These Pauli operators are used to construct the creation and annihilation operators \((P \pm i Q)/2\).

The generated table is processed by sparse_pauli_operators().

Parameters:

register_length – the register length for which to generate the table.

Returns:

A list of tuples of two Pauli string operators.

sparse_pauli_operators(register_length: int) tuple[list[SparsePauliOp], list[SparsePauliOp]][source]

Generates the SparsePauliOp terms.

This uses pauli_table() to construct a list of operators used to translate the second-quantization symbols into qubit operators.

Parameters:

register_length – the register length for which to generate the operators.

Returns:

Two lists stored in a tuple, consisting of the creation and annihilation operators, applied on the individual modes.

myutils submodule

Utility functions for NuQuLib.

This module provides various utility functions for NuQuLib, including format conversions between different quantum computing frameworks, LaTeX formatting for nuclear notation, etc.

class nuqulib.myutils.Orbit_nlj(n: int, l: int, j: int)[source]

Bases: object

Single particle state in the model space with n, l, j quanta.

class nuqulib.myutils.Orbit_nljjztz(n: int, l: int, j: int, jz: int, tz: int)[source]

Bases: object

Single particle state in the model space with n, l, j, jz, tz quanta.

class nuqulib.myutils.Orbit_nljtz(n: int, l: int, j: int, tz: int)[source]

Bases: object

Single particle state in the model space with n, l, j, tz quanta.

nuqulib.myutils.count_msps(emax: int, vemin: int = 0, vemax: int = 100)[source]

Count the number of single particle states in the model space up to emax.

nuqulib.myutils.get_operator_from_QiskitStr(pauli_str: str) Operator[source]

Convert a Qiskit Pauli string to a PennyLane operator.

Parameters:

pauli_str (str) – Pauli string in Qiskit format (e.g., “IXYZ”).

Returns:

PennyLane operator corresponding to the Pauli string.

Return type:

qml.operation.Operator

Raises:

ValueError – If the Pauli string contains invalid characters.

Note

The function reverses the Pauli string to match PennyLane’s qubit ordering.

nuqulib.myutils.get_spsidx_from_nljtz(single_particle_states: list, n, l, j, tz)[source]

Find a single-particle-state index from quantum numbers.

Parameters:
  • single_particle_states (list) – Sequence of states with n, l, j, and tz attributes.

  • n (int) – Radial quantum number.

  • l (int) – Orbital angular momentum.

  • j (int) – Twice the total angular momentum.

  • tz (int) – Isospin projection.

Returns:

Index of the matching single-particle state.

Return type:

int

Raises:

ValueError – If no matching state exists.

nuqulib.myutils.latex_nuc(nuc: str) str[source]

Convert nuclear notation to LaTeX format.

Parameters:

nuc (str) – Nuclear notation string (e.g., “6He”, “12C”).

Returns:

LaTeX formatted string with mass number as superscript.

Return type:

str

Example

>>> latex_nuc("6He")
'${}^{6}$He'
nuqulib.myutils.read_QiskitPauli(ops_qiskit, coeffs_qiskit)[source]

Convert Qiskit Pauli operators to PennyLane format.

Parameters:
  • ops_qiskit (list) – List of Qiskit Pauli operator strings.

  • coeffs_qiskit (list) – List of coefficients corresponding to operators.

Returns:

Tuple containing:
  • coeffs (list[float]): List of float coefficients.

  • obs (list): List of PennyLane operators.

Return type:

tuple

nuqulib.myutils.t_count(eps_tol=1e-09, model='ross-selinger', precision=None)[source]

Estimate T-count needed for single-qubit rotation synthesis.

Parameters:
  • eps_tol (float, optional) – Target synthesis error tolerance.

  • model (str, optional) – Estimate model. Supported values are "ross-selinger" and "bocharov".

  • precision – Deprecated compatibility argument; currently unused.

Returns:

Estimated T-count.

Return type:

int

Raises:

ValueError – If eps_tol is non-positive or the model is unknown.

nuqulib.myutils.transform_pytket_counts_to_qiskit(counts_pytket)[source]

Convert PyTKET measurement counts to Qiskit format.

Parameters:

counts_pytket (dict) – Dictionary of measurement outcomes from PyTKET.

Returns:

Dictionary with bit strings reversed to match Qiskit ordering.

Return type:

dict

Note

PyTKET and Qiskit use different qubit ordering conventions. This function reverses the bit string keys to convert between them.

nuqulib.myutils.transform_qiskitOps_to_pennylane(Qiskit_Ops: SparsePauliOp)[source]

Transform Qiskit SparsePauliOp operator(s) to PennyLane format.

Parameters:

Qiskit_Ops (SparsePauliOp) – Qiskit sparse Pauli operator(s).

Returns:

Tuple containing:
  • coeffs (list[float]): List of float coefficients.

  • obs (list): List of PennyLane operators.

Return type:

tuple

nuqulib.myutils.write_out_pytketCircuit_from_Qiskit(qc: QuantumCircuit, fname: str, circ_name: str = '')[source]

Save a single Qiskit circuit as PyTKET circuit to a numpy file.

Parameters:
  • qc (QuantumCircuit) – Qiskit quantum circuit to convert and save.

  • fname (str) – Output filename (will be saved as .npy file).

  • circ_name (str, optional) – Name for the circuit. Defaults to empty string.

Note

The circuit is saved as a dictionary with ‘circuit’ and ‘name’ keys using numpy’s save function with pickle enabled.

nuqulib.myutils.write_out_pytketCircuits_from_Qiskit(qc: list[QuantumCircuit], fname: str, circ_names: list[str] = None)[source]

Save a list of Qiskit circuits as PyTKET circuits to a numpy file.

Parameters:
  • qc (list[QuantumCircuit]) – List of Qiskit quantum circuits.

  • fname (str) – Output filename (will be saved as .npy file).

  • circ_names (list[str], optional) – Names for the circuits. Defaults to empty strings if not provided.

Note

The circuits are saved as a dictionary with ‘circuit’ and ‘name’ keys using numpy’s save function with pickle enabled.

nft_opt_method submodule

Nakanishi-Fujii-Todo (NFT) optimization method for VQE.

This module implements the Nakanishi-Fujii-Todo optimization algorithm for Variational Quantum Eigensolver (VQE) problems. The NFT method is a derivative-free and sequential optimization algorithm that uses discrete Fourier transform to optimize the gate angles by exploiting the periodic structure on the parameter space.

Ref: KM.Nakanishi, K.Fujii, T.Todo, Phys. Rev. Research 2, 043158 (2020),

nuqulib.nft_opt_method.DiscreteFT(it: int, Es: Iterable[float], spot: Iterable[float], theta_scale_for_DFT: float, verbose: bool, method: str = 'scipy', eps: float = 1e-08)[source]

Perform discrete Fourier transform optimization step.

Uses DFT to fit energy measurements and find optimal parameter value by minimizing the reconstructed energy landscape.

Parameters:
  • it (int) – Current iteration number.

  • Es (array) – Energy measurements.

  • spot (array) – Parameter values where energies were measured.

  • theta_scale_for_DFT (float) – Scaling factor for DFT.

  • verbose (bool) – Whether to print debug information.

  • method (str, optional) – Optimization method. Options: “scipy”, “naive”. Defaults to “scipy”.

  • eps (float, optional) – Numerical precision parameter. Defaults to 1.0e-8.

Returns:

Tuple containing:
  • new_theta (float): Optimal parameter value.

  • y_min (float): Minimum energy value.

Return type:

tuple

nuqulib.nft_opt_method.NFTmethod(it, Ecurrent, hamiltonian_op, params, Nq, Nocc, method_ansatz, method_measure, where_is_G_or_cG1, which_Gate, verbose, pairinghamiltonian=True, proton_number=0, neutron_number=0, proton_qubits=[], neutron_qubits=[])[source]

Perform one step of NFT method optimization

Optimizes a single parameter using the NFT method by measuring energies at periodic intervals and using DFT to find the optimal parameter value.

Parameters:
  • it (int) – Current iteration number.

  • Ecurrent (float) – Current energy value.

  • hamiltonian_op – Hamiltonian operator.

  • params (array) – Current parameter values.

  • Nq (int) – Number of qubits.

  • Nocc (int) – Number of occupied orbitals.

  • method_ansatz (str) – Ansatz method.

  • method_measure (str) – Measurement method.

  • where_is_G_or_cG1 (list) – List indicating gate types (“G” or “cG1”).

  • which_Gate (int) – Index of gate to optimize.

  • verbose (bool) – Whether to print debug information.

  • pairinghamiltonian (bool, optional) – Whether using pairing Hamiltonian. Defaults to True.

  • proton_number (int, optional) – Number of protons. Defaults to 0.

  • neutron_number (int, optional) – Number of neutrons. Defaults to 0.

  • proton_qubits (list, optional) – Proton qubit indices. Defaults to [].

  • neutron_qubits (list, optional) – Neutron qubit indices. Defaults to [].

Returns:

Updated parameter values with optimized parameter.

Return type:

array

nuqulib.nft_opt_method.cost_func(hamiltonian_op, params: Iterable[float], Nq: int, Nocc: int, method_ansatz: str, method_measure: str, pairinghamiltonian: bool = True, proton_number: int = 0, neutron_number: int = 0, proton_qubits=[], neutron_qubits=[])[source]

Evaluate cost function for VQE optimization with NFT method.

Computes the expectation value of the Hamiltonian for given parameters using the specified ansatz and measurement method.

Parameters:
  • hamiltonian_op – Hamiltonian operator to evaluate.

  • params (array) – Variational parameters.

  • Nq (int) – Number of qubits.

  • Nocc (int) – Number of occupied orbitals.

  • method_ansatz (str) – Ansatz method to use.

  • method_measure (str) – Measurement method (“statevector” or “Aer”).

  • pairinghamiltonian (bool, optional) – Whether using pairing Hamiltonian. Defaults to True.

  • proton_number (int, optional) – Number of protons. Defaults to 0.

  • neutron_number (int, optional) – Number of neutrons. Defaults to 0.

  • proton_qubits (list, optional) – Proton qubit indices. Defaults to [].

  • neutron_qubits (list, optional) – Neutron qubit indices. Defaults to [].

Returns:

Energy expectation value.

Return type:

float

Raises:

ValueError – If invalid measurement method is provided.

nuqulib.nft_opt_method.draw_DFT_curve(xplot, DFTcurve, spot, Es, x_min, y_min)[source]

Draw DFT curve with measured points and minimum candidate.

Parameters:
  • xplot (array) – X-axis values for plotting.

  • DFTcurve (array) – DFT curve values.

  • spot (array) – Measurement points.

  • Es (array) – Energy values at measurement points.

  • x_min (float) – X-coordinate of minimum candidate.

  • y_min (float) – Y-coordinate of minimum candidate.

nuqulib.nft_opt_method.eval_DFT_coeff(Es, spot, theta_scale_for_DFT: float = 1.0)[source]

Evaluate discrete Fourier transform coefficients.

Computes DFT coefficients from energy measurements at specific points to reconstruct the periodic energy landscape.

Parameters:
  • Es (array) – Energy measurements.

  • spot (array) – Parameter values where energies were measured.

  • theta_scale_for_DFT (float) – Optional scaling factor for theta in DFT.

Returns:

DFT coefficients.

Return type:

array

nuqulib.nft_opt_method.func_DFT_coeff(DFT_coef: Iterable[float], theta_scale_for_DFT=1.0)[source]

Create function from DFT coefficients for energy landscape reconstruction.

Parameters:
  • DFT_coef – DFT coefficients.

  • theta_scale_for_DFT – Scaling factor for theta. Defaults to 1.0.

Returns:

Function that evaluates DFT curve at given parameter values.

Return type:

function

nuqulib.nft_opt_method.optimize_params_with_NFT(it_max: int, hamiltonian_op, params: Iterable[float], Nq: int, Nocc: int, ngate: int, where_is_G_or_cG1: dict, method_ansatz: str, method_measure: str, verbose: bool = False, pairinghamiltonian: bool = True, proton_number: int = 0, neutron_number: int = 0, proton_qubits=[], neutron_qubits=[])[source]

Optimize ansatz parameters with sequential NFT updates.

Parameters:
  • it_max (int) – Number of outer optimization sweeps.

  • hamiltonian_op – Hamiltonian operator used by the cost function.

  • params (Iterable[float]) – Initial variational parameters.

  • Nq (int) – Number of qubits.

  • Nocc (int) – Number of occupied orbitals for pairing-style ansätze.

  • ngate (int) – Number of variational gates/parameters to update.

  • where_is_G_or_cG1 (dict) – Mapping from parameter index to gate type.

  • method_ansatz (str) – Ansatz construction method passed to cost_func.

  • method_measure (str) – Measurement method passed to cost_func.

  • verbose (bool, optional) – If True, print energy and parameter updates.

  • pairinghamiltonian (bool, optional) – Whether to use pairing-Hamiltonian cost-function conventions.

  • proton_number (int, optional) – Proton number for nuclear ansätze.

  • neutron_number (int, optional) – Neutron number for nuclear ansätze.

  • proton_qubits (list, optional) – Proton qubit indices.

  • neutron_qubits (list, optional) – Neutron qubit indices.

Returns:

Optimized parameter vector.

Return type:

Iterable[float]

pairing_hamiltonian submodule

Pairing Hamiltonian implementation

This module implements the pairing Hamiltonian, a simplified model used in nuclear physics, condensed matter physics, and quantum computing. The pairing Hamiltonian describes the interactions between pairs of particles in a many-body system, making it a valuable tool for studying phenomena like superconductivity and nuclear pairing effects.

This implementation is not optimized for performance but serves as a pedagogical example of how to construct and manipulate a pairing Hamiltonian in Python. More detailed and efficient implementations can be found as a Julia package PairingHamiltonians.jl by the same author.

class nuqulib.pairing_hamiltonian.PairingHamiltonian(Norb, Nocc, gval, delta_eps=1.0)[source]

Bases: object

Pairing Hamiltonian model for nuclear many-body systems.

This class implements a pairing Hamiltonian with single-particle energies with a fixed spacing and a constant pairing interaction strength.

Parameters:
  • Norb (int) – Number of single-particle orbitals.

  • Nocc (int) – Number of occupied orbitals.

  • gval (float) – Pairing interaction strength.

  • delta_eps (float, optional) – Single-particle energy spacing. Defaults to 1.0.

Norb

Number of orbitals.

Type:

int

Nocc

Number of occupied states.

Type:

int

gval

Pairing strength parameter.

Type:

float

delta_eps

Energy level spacing.

Type:

float

basis

List of basis states as orbital occupation tuples.

Type:

list

epsilon

Single-particle energies.

Type:

list

Hmat

Hamiltonian matrix.

Type:

numpy.ndarray

diagonalize()[source]

Diagonalize the Hamiltonian matrix.

Performs exact diagonalization to find eigenvalues and eigenvectors of the Hamiltonian matrix.

Returns:

Tuple containing:
  • evals (numpy.ndarray): Eigenvalues in ascending order.

  • evecs (numpy.ndarray): Corresponding eigenvectors as columns.

Return type:

tuple

encoding()[source]

Encode the Hamiltonian as Pauli operators for quantum computing.

Maps the pairing Hamiltonian to a sum of Pauli operators suitable for quantum algorithms. The encoding includes:

  • Identity terms for constant energy shifts

  • Z terms for single-particle energies

  • XX and YY terms for pairing interactions

Returns:

Qiskit SparsePauliOp representing the Hamiltonian.

Return type:

SparsePauliOp

Note

The qubit ordering is reversed to match Qiskit conventions.

eval_Hmat()[source]

Construct the Hamiltonian matrix.

Builds the Hamiltonian matrix in the many-body basis including:

  • Single-particle energy terms (diagonal)

  • Pairing interaction terms (off-diagonal for single excitations)

Returns:

Hamiltonian matrix of size (dim, dim) where dim is the basis size.

Return type:

numpy.ndarray

eval_epsilon()[source]

Calculate single-particle energies.

Sets up a linear spacing of single-particle energies with epsilon_i = 2 * i * delta_eps for orbital i.

Returns:

List of single-particle energies.

Return type:

list

make_basis()[source]

Generate basis states for the pairing model.

Creates all possible configurations of Nocc particles in Norb orbitals using combinations.

Returns:

List of tuples representing occupied orbital configurations.

Return type:

list

nuqulib.pairing_hamiltonian.tuple_to_bitstring(tup, Norb, rev=True)[source]

Convert an orbital occupation tuple to bitstring notation.

Parameters:
  • tup (tuple) – Tuple of occupied orbital indices.

  • Norb (int) – Total number of orbitals.

  • rev (bool, optional) – Whether to reverse the bitstring. Defaults to True.

Returns:

Bitstring representation in ket notation (e.g.,

"|0110>").

Return type:

str

Example

>>> tuple_to_bitstring((0, 2), 4, rev=True)
'|0101>'

operators submodule

Functions to build operators in the M-scheme

pairwise submodule

Pairwise Hamiltonian utilities

This module provides utilities for constructing and manipulating pairwise Hamiltonians, which are often used in quantum computing and many-body physics. It includes functions for reading Hamiltonian data from files, generating configurations, and evaluating Hamiltonian matrices in many-body bases.

These utilities are designed to work with pairing Hamiltonians or pairwise representations of shell model (e.g., time-reversal pairs).

nuqulib.pairwise.add_Vmonopole(Hamil: dict, config: str, n_qubits: int, Nocc: int, Dict_qubits_to_sps: dict, neutron=True, verbose=False)[source]

Add monopole term to the Hamiltonian matrix.

This function calculates the monopole term of the Hamiltonian based on the given configuration and the Hamiltonian data. One should consider the monopole term only if the number of occupied states (Nocc) under a shell-model type Hamiltonian, is greater than 1, as it involves pairs of particles.

Parameters:
  • Hamil – Dictionary containing Hamiltonian matrix elements.

  • config – Configuration string representing the occupation of orbitals.

  • n_qubits – Number of qubits (or orbitals).

  • Nocc – Number of occupied states.

  • Dict_qubits_to_sps – Dictionary mapping qubits to single-particle states.

  • neutron (optional) – If True, consider neutron interactions; otherwise, proton interactions.

  • verbose (optional) – If True, print detailed information.

nuqulib.pairwise.eval_Hflat_eigen(n_qubits: int, Nocc: int, Hamil: dict, h1b: ndarray, h2b: ndarray, Dict_qubits_to_sps: dict, verbose=False)[source]

Evaluate the Hamiltonian matrix in the many-body basis.

Constructs the Hamiltonian matrix in the many-body basis explicitly, and then computes its eigenvalues and eigenvectors.

Parameters:
  • n_qubits – Number of qubits.

  • Nocc – Number of occupied states.

  • Hamil – Dictionary containing Hamiltonian matrix elements.

  • h1b – One-body Hamiltonian matrix.

  • h2b – Two-body Hamiltonian matrix.

  • Dict_qubits_to_sps – Dictionary mapping qubits to single-particle states.

  • verbose – If True, print detailed information.

Returns:

Tuple containing:
  • evals: Eigenvalues of the Hamiltonian matrix.

  • evecs: Eigenvectors of the Hamiltonian matrix.

Return type:

tuple

nuqulib.pairwise.generate_config_bitstr_list(Nq: int, Nocc: int, rev=False)[source]

Generate a list of strings representing configurations of occupied qubits.

nuqulib.pairwise.get_pairwise_Hamil(Hamil: dict, adopted_sps: list, neutron=True)[source]

Generate the pairwise Hamiltonian matrix elements.

This function constructs the one-body and two-body Hamiltonian matrices based on the provided Hamiltonian data and pair-wise representation of single-particle states.

nuqulib.pairwise.get_possible_time_reversal_pairs(sps: dict, verbose=False)[source]

Generate all time-reversal pairs of single-particle states

Parameters:
  • sps (dict) – Dictionary of single-particle states, where keys are state identifiers and values are tuples containing (n, l, j, mz, SPE).

  • verbose (bool, optional) – If True, print each detected pair.

Note

sps should be a dictionary for a single species (either protons or neutrons).

nuqulib.pairwise.get_v0_pairs(msps_p: Orbit_nljjztz, msps_n: Orbit_nljjztz)[source]

Get the list of time-reversal pairs with Jz=0

nuqulib.pairwise.make_pw_hamil_qiskit(Hamil: dict, h1b: ndarray, Nq: int, Nocc: int, Dict_qubits_to_sps: dict)[source]

Construct the pairwise Hamiltonian in Qiskit format.

This function constructs the pair-wise Hamiltonian in the Qiskit format, which includes single-particle-like terms, pair terms, and monopole terms. Since we are working within the so-called zero-seniority approximation a.k.a. Hard-core boson representation, the Jordan-Wigner transformation leads to only the following terms:

  • Single-particle-like terms: (e_i + Vii) * (I-Z)/2

  • Pair terms: (XX + YY)

  • Monopole terms: (I-Zi)(I-Zj)/4

Parameters:
  • Hamil – Dictionary containing Hamiltonian matrix elements.

  • h1b – One-body Hamiltonian matrix elements.

  • Nq – Number of qubits (or orbitals).

  • Nocc – Number of occupied states.

  • Dict_qubits_to_sps – Dictionary mapping qubits to single-particle states.

nuqulib.pairwise.read_msnt(fn, A: int, Acore: int, pow_A: float = -0.3, massdep: int = 0, degenerate_SPE: bool = False, verbose: bool = False, neutron: bool = True, pn_system: bool = False, only_monopole: bool = False)[source]

Read a file in the msnt format, snt format with explicit jz(m), and return the dictionary of the Hamiltonian, single-particle states, and dictionaries to convert between single-particle states and qubits.

Note

The current implementation assumes that the file is formatted correctly and outputs will be used for systems consisting of single species nucleons (either protons or neutrons).

postselection submodule

Implementation of post-selection methods for quantum simulations.

This module provides post-selection and measurement analysis tools for quantum nuclear physics simulations, including particle number conservation constraints and expectation value calculations for different types of Hamiltonians.

nuqulib.postselection.eval_Ediag(adopted: str, Nq: int, Nocc: int, hamiltonian_op_diag: SparsePauliOp, qc_ansatz: QuantumCircuit, backend, sampler, nshot: int, postselection_diag=True)[source]

Evaluate diagonal energy terms using quantum circuit measurements.

Runs a quantum circuit to measure diagonal Hamiltonian terms and computes their expectation values with optional post-selection.

Parameters:
  • adopted – Simulation mode (“simNISQ”, “Real”, etc.).

  • Nq – Number of qubits.

  • Nocc – Number of occupied orbitals.

  • hamiltonian_op_diag (SparsePauliOp) – Diagonal Hamiltonian terms.

  • qc_ansatz (QuantumCircuit) – Ansatz circuit to evaluate.

  • backend – Quantum backend for transpilation/execution.

  • sampler – Quantum sampler for measurements.

  • nshot – Number of measurement shots.

  • postselection_diag (bool, optional) – Whether to apply particle number post-selection. Defaults to True.

Returns:

Expectation value of diagonal Hamiltonian terms.

Return type:

float

nuqulib.postselection.eval_Energy_using_GoogleCircuit(Nq: int, Nocc: int, hamiltonian_op_XXYY: SparsePauliOp, qc_list_XXYY: list, sampler, nshot: int, num_experiment: int, using_noisy_simulation: bool, postselection_XXYY: bool = True, debug_mode: bool = False, verbose: bool = False)[source]

Evaluate XX+YY energy terms using Google’s circuit approach.

This function evaluates the expectation values of XX and YY Hamiltonian terms using the circuit design from Google’s quantum simulation approach. It handles the specific qubit ordering and measurement strategies for optimal performance.

Parameters:
  • Nq (int) – Number of qubits.

  • Nocc (int) – Number of occupied orbitals.

  • hamiltonian_op_XXYY (SparsePauliOp) – Hamiltonian with XX+YY terms.

  • qc_list_XXYY (list) – List of quantum circuits for measuring XX+YY terms.

  • sampler – Quantum sampler for measurements.

  • nshot (int) – Number of measurement shots per circuit.

  • num_experiment (int) – Number of experimental repetitions.

  • using_noisy_simulation (bool) – Whether using noisy quantum simulation.

  • postselection_XXYY (bool, optional) – Whether to apply post-selection. Defaults to True.

  • debug_mode (bool, optional) – Enable debug output. Defaults to False.

  • verbose (bool, optional) – Enable verbose output. Defaults to False.

Returns:

Energy values from multiple experimental runs.

Return type:

list

Note

  • Even circuits measure XX+YY for i=0,2,4,…

  • Odd circuits measure XX+YY for i=1,3,5,…

  • For odd Nq systems, special handling of X_0X_1 term is required

  • Qubit indices are flipped due to Qiskit ordering conventions

nuqulib.postselection.eval_XXYY_from_Hcounts(counts, Nq, hamiltonian_op_XXYY, return_counts_as_well)[source]

Evaluate XX+YY contribution from counts measured in the H basis.

Parameters:
  • counts (dict) – Bitstring counts or probabilities.

  • Nq (int) – Number of qubits.

  • hamiltonian_op_XXYY (SparsePauliOp) – Hamiltonian terms containing XX/YY pair operators.

  • return_counts_as_well (bool) – If True, return (energy, counts).

Returns:

Estimated energy contribution, optionally with counts.

Return type:

float or tuple

nuqulib.postselection.eval_XXYY_w_Hgates(adopted, Nq: int, Nocc: int, hamiltonian_op_XXYY: SparsePauliOp, qc_ansatz: QuantumCircuit, sampler, nshot: int, using_noisy_simulation=False, backend=None, return_counts_as_well=False)[source]

Evaluate XX+YY terms by inserting Hadamard gates before measurement.

Parameters:
  • adopted (str) – Execution target label, such as "Real".

  • Nq (int) – Number of qubits.

  • Nocc (int) – Target occupation number used by count reweighting.

  • hamiltonian_op_XXYY (SparsePauliOp) – Hamiltonian terms to evaluate.

  • qc_ansatz (QuantumCircuit) – Ansatz circuit to measure.

  • sampler – Qiskit sampler-like object.

  • nshot (int) – Number of measurement shots.

  • using_noisy_simulation (bool, optional) – If True, read counts from the noisy-simulator result layout.

  • backend (optional) – Backend for transpilation.

  • return_counts_as_well (bool, optional) – If True, return the reweighted counts together with the energy.

Returns:

Estimated XX+YY contribution, optionally with counts.

Return type:

float or tuple

nuqulib.postselection.eval_expec_diag_postselection(hamiltonian_op_diag: SparsePauliOp, counts: dict, Nq: int, Nocc: int)[source]

Evaluate expectation value of diagonal Hamiltonian terms with post-selection.

Computes the expectation value of Hamiltonian terms that are diagonal in the computational basis (I and Z terms) using post-selected measurement results.

Parameters:
  • hamiltonian_op_diag (SparsePauliOp) – Hamiltonian with only diagonal terms.

  • counts (dict) – Dictionary of measurement results.

  • Nq (int) – Number of qubits.

  • Nocc (int) – Number of occupied orbitals for post-selection.

Returns:

Expectation value of the diagonal Hamiltonian terms.

Return type:

float

nuqulib.postselection.even_swap(labeling)[source]

Swap elements at even indices with their next neighbors.

Parameters:

labeling (list) – Input list to perform swaps on.

Returns:

New list with even-indexed elements swapped.

Return type:

list

Example

>>> even_swap([0, 1, 2, 3, 4, 5])
[1, 0, 3, 2, 5, 4]
nuqulib.postselection.odd_swap(labeling)[source]

Swap elements at odd indices with their next neighbors.

Parameters:

labeling (list) – Input list to perform swaps on.

Returns:

New list with odd-indexed elements swapped.

Return type:

list

Example

>>> odd_swap([0, 1, 2, 3, 4, 5])
[0, 2, 1, 4, 3, 5]
nuqulib.postselection.reweight_results(sampler_results: dict, Norb: int, Nocc: int, normalize: bool = True, verbose: bool = False, make_number_conservation: bool = True)[source]

Reweight quantum measurement results to enforce particle number conservation.

This function filters and reweights measurement results to include only configurations with the correct number of particles, which is essential for nuclear quantum simulations.

Parameters:
  • sampler_results – Raw measurement results from quantum simulation. The keys are bitstrings or hexadecimal strings, and the values are counts or probabilities.

  • Norb – Total number of orbitals (qubits).

  • Nocc – Number of occupied orbitals (target particle number).

  • normalize (optional) – Whether to normalize results to unity. Defaults to True.

  • verbose (optional) – Whether to print valid counts ratio. Defaults to False.

  • make_number_conservation (optional) – Whether to enforce particle number conservation. Defaults to True.

Returns:

Reweighted measurement results with proper particle number.

Return type:

dict

Note

Automatically handles both bitstring and hexadecimal input formats.

nuqulib.postselection.single_eval_XXYY_Google(Nq: int, Nocc: int, qc_list_XXYY: list, list_counts_G: list, hamiltonian_op_XXYY: SparsePauliOp, postselection_XXYY: bool)[source]

Evaluate XX+YY energy from Google-style basis-rotation counts.

Parameters:
  • Nq (int) – Number of qubits.

  • Nocc (int) – Target occupation number used for postselection.

  • qc_list_XXYY (list) – Circuits corresponding to the XX+YY measurement schedule.

  • list_counts_G (list) – Counts or probability dictionaries for each circuit in qc_list_XXYY.

  • hamiltonian_op_XXYY (SparsePauliOp) – Hamiltonian terms containing XX/YY pair operators.

  • postselection_XXYY (bool) – If True, enforce number conservation when reweighting counts.

Returns:

Estimated XX+YY energy contribution.

Return type:

float

quantum_algorithms submodule

class nuqulib.quantum_algorithms.ODMD(Uprep: QuantumCircuit, HamiltonianOps: SparsePauliOp, delta_t: float, max_iterations: int, trotter_rank: int, trotter_steps: int, sampler, ancilla_qubits, target_qubits, num_shot: int = 1, using_statevector: bool = True, dim_Hankel: int = 8, tol_SVD: float = 1e-08, verbose: bool = False, plot_lambda: bool = False, tol_lambda: float = 0.01)[source]

Bases: object

Online dynamic mode decomposition solver for quantum time snapshots.

construct_A_from_XY(X, Y)[source]

Fit the reduced linear propagator from Hankel matrices.

Parameters:
  • X (np.ndarray) – Source Hankel matrix.

  • Y (np.ndarray) – Shifted target Hankel matrix.

Returns:

Least-squares propagator A = Y X^+ after SVD truncation.

Return type:

np.ndarray

construct_X_and_Y(snapshots, dim_Hankel)[source]

Construct shifted Hankel matrices from a snapshot sequence.

Parameters:
  • snapshots (np.ndarray) – Complex time-snapshot sequence.

  • dim_Hankel (int) – Number of Hankel rows.

Returns:

Shifted Hankel matrices X and Y.

Return type:

tuple[np.ndarray, np.ndarray]

estimate_resource(tol=1e-10, model='ross-selinger', verbose=False)[source]

Estimate T-count for the ODMD overlap measurements.

Parameters:
  • tol (float, optional) – Rotation synthesis tolerance.

  • model (str, optional) – Rotation synthesis estimate model.

  • verbose (bool, optional) – If True, print intermediate estimates.

Returns:

Estimated T-count.

Return type:

float

get_results()[source]

Return cached ODMD results from the most recent run() call.

run()[source]

Run ODMD snapshot generation, fitting, and energy extraction.

Returns:

Extracted energies and selected unit-circle eigenvalues.

Return type:

tuple[list[float], np.ndarray]

nuqulib.quantum_algorithms.QuantumKrylov(Uprep: QuantumCircuit, hamiltonian_op: SparsePauliOp, sampler, ancilla_qubits, target_qubits, delta_t=0.01, max_iterations=10, trotter_rank=2, trotter_steps=1, num_shot=10000, using_statevector=False, do_simulation=True, Bosonic=True, verbose=False, tol_eig=1e-06)[source]

Backward-compatible wrapper around QuantumKrylovSolver.

class nuqulib.quantum_algorithms.QuantumKrylovSolver(Uprep: QuantumCircuit, hamiltonian_op: SparsePauliOp, sampler, ancilla_qubits, target_qubits, delta_t: float = 0.01, max_iterations: int = 10, trotter_rank: int = 2, trotter_steps: int = 1, num_shot: int = 10000, using_statevector: bool = False, do_simulation: bool = True, Bosonic: bool = False, verbose: bool = False, tol_eig: float = 1e-06)[source]

Bases: object

Design:
  • prepare_* methods: circuit construction only

  • evaluate_* methods: execute circuits and compute matrix elements

  • run: iterative Krylov loop

  • estimate_resource: coarse resource estimate without simulation

estimate_resource(tol=1e-10, model='ross-selinger', reduction_factor=1, verbose=False)[source]

Coarse resource estimate without running simulation.

evaluate_diag_element(qcs, idxs_circuit)[source]

Evaluate diagonal matrix element H[it, it].

evaluate_offdiag_element(gate_cUi, gate_cUj, term_types, idxs_circuit)[source]

Evaluate off-diagonal matrix element H[it, j].

prepare_iteration_circuits(it)[source]

Build all circuits needed at iteration it.

run()[source]

Execute Quantum Krylov iterations and return (Hsub, Nsub, ws).

solve_projected_eigenproblem(Nsub, Hsub)[source]

Solve projected generalized eigenvalue problem for current Krylov subspace.

nuqulib.quantum_algorithms.T_formula_QFT(N_ancilla, eps_rotation: float = 1e-10)[source]

Estimate T-count for the inverse QFT rotation gates.

Parameters:
  • N_ancilla (int) – Number of QPE ancilla qubits.

  • eps_rotation (float, optional) – Rotation synthesis tolerance.

Returns:

Estimated T-count.

Return type:

float

nuqulib.quantum_algorithms.circuit_HadamardTest(Norb, Uprep, Hamiltonian_op, t, trotter_steps, using_statevector=True)[source]

Build a Hadamard-test circuit for a time-evolution operator.

Parameters:
  • Norb (int) – Number of target orbitals/qubits.

  • Uprep (QuantumCircuit) – State-preparation circuit on the target register.

  • Hamiltonian_op (SparsePauliOp) – Hamiltonian used for Pauli evolution.

  • t (float) – Evolution time.

  • trotter_steps (int) – Number of Suzuki-Trotter repetitions.

  • using_statevector (bool, optional) – If True, omit measurement and return a statevector-ready circuit.

Returns:

Decomposed Hadamard-test circuit.

Return type:

QuantumCircuit

nuqulib.quantum_algorithms.circuit_my_QPE(n_ancilla: int, Norb: int, Hamiltonian_op: SparsePauliOp, Uprep: QuantumCircuit, time: float, measure=False, trotter_order: int = 2, trotter_steps: int = 1, repeat: bool = False)[source]

Construct a basic quantum phase-estimation circuit.

Parameters:
  • n_ancilla (int) – Number of counting qubits.

  • Norb (int) – Number of target qubits.

  • Hamiltonian_op (SparsePauliOp) – Hamiltonian used for time evolution.

  • Uprep (QuantumCircuit) – State-preparation circuit.

  • time (float) – Base evolution time.

  • measure (bool, optional) – If True, measure the ancilla register.

  • trotter_order (int, optional) – Suzuki-Trotter order.

  • trotter_steps (int, optional) – Number of Suzuki-Trotter repetitions.

  • repeat (bool, optional) – If True, build powers by repeated composition instead of scaling the evolution time.

Returns:

QPE circuit.

Return type:

QuantumCircuit

nuqulib.quantum_algorithms.get_idx_circuit(op_string, term_types)[source]

Find the prepared measurement circuit matching a Pauli string type.

Parameters:
  • op_string (str) – Pauli label to classify.

  • term_types (list[str]) – Available measurement basis descriptors.

Returns:

Index of the matching circuit type.

Return type:

int

nuqulib.quantum_algorithms.lambda_plot(lam, Ens)[source]

Save a complex-plane plot of selected ODMD eigenvalues.

Parameters:
  • lam (Iterable[complex]) – Eigenvalues to plot.

  • Ens (Iterable[float]) – Energies used in legend labels.

nuqulib.quantum_algorithms.make_Circ_forNondiagH(term_types, Ntar, ancilla_qubits, target_qubits, gate_cUi, gate_cUj, qcs_re, qcs_im, using_statevector)[source]

Append measurement circuits for non-diagonal Hamiltonian term types.

Parameters:
  • term_types (list[str]) – Measurement basis descriptors.

  • Ntar (int) – Number of target qubits.

  • ancilla_qubits (list) – Ancilla qubits used by the controlled gates.

  • target_qubits (list) – Target qubit indices.

  • gate_cUi – Controlled unitary for the first state.

  • gate_cUj – Controlled unitary for the second state.

  • qcs_re (list) – Output list receiving real-part circuits.

  • qcs_im (list) – Output list receiving imaginary-part circuits.

  • using_statevector (bool) – If False, add measurements.

Returns:

The output lists are modified in place.

Return type:

None

nuqulib.quantum_algorithms.make_cU(Uprep, Ui, Ntar)[source]

Create a controlled gate for Uprep followed by Ui.

nuqulib.quantum_algorithms.make_overlap_qc(Ntar, gate_cUi, gate_cUj, ancilla_qubits, target_qubits, using_statevector)[source]

Build circuits for real and imaginary parts of an overlap test.

Parameters:
  • Ntar (int) – Number of target qubits.

  • gate_cUi – Controlled unitary for the first state.

  • gate_cUj – Controlled unitary for the second state.

  • ancilla_qubits (list) – Ancilla qubits used by the controlled gates.

  • target_qubits (list) – Target qubits acted on by the controlled gates.

  • using_statevector (bool) – If False, add measurements.

Returns:

Circuits for real and imaginary overlap components.

Return type:

tuple[QuantumCircuit, QuantumCircuit]

nuqulib.quantum_algorithms.measure_overlap(num_shot, Ntar, gate_cUi, gate_cUj, ancilla_qubits, target_qubits, sampler, using_statevector, do_simulation=True)[source]

Estimate an overlap from Hadamard-test circuits.

Parameters:
  • num_shot (int) – Number of shots for sampled execution.

  • Ntar (int) – Number of target qubits.

  • gate_cUi – Controlled unitary for the first state.

  • gate_cUj – Controlled unitary for the second state.

  • ancilla_qubits (list) – Ancilla qubits used by the controlled gates.

  • target_qubits (list) – Target qubits acted on by the controlled gates.

  • sampler – Sampler or backend-like object used for transpilation/execution.

  • using_statevector (bool) – If True, use statevector simulation.

  • do_simulation (bool, optional) – If False, print resource information and skip execution.

Returns:

Estimated overlap, or None for resource-only mode.

Return type:

complex | None

class nuqulib.quantum_algorithms.myTextBookQPE(n_ancilla, Norb, Hamiltonian_op, Uprep, time, trotter_order=2, trotter_steps=1)[source]

Bases: object

Small wrapper around textbook quantum phase estimation.

construct_circuit()[source]

Construct the measured textbook QPE circuit.

estimate_resource(tol=1e-10, verbose=False)[source]

Estimate the T-count for the stored QPE circuit parameters.

Parameters:
  • tol (float, optional) – Rotation synthesis tolerance.

  • verbose (bool, optional) – If True, print intermediate estimates.

Returns:

The estimate is printed to stdout.

Return type:

None

nuqulib.quantum_algorithms.prepare_qc_for_QKrylov(Hamiltonian_op, Uprep, Ui, Ntar, Bosonic, using_statevector=False, verbose=True)[source]

Prepare quantum circuits for evaluating Hamiltonian terms in QKrylov method. This could be also used for VQE-type algorithms where one needs to evaluate the expectation value of Hamiltonian terms.

nuqulib.quantum_algorithms.reorder_based_on_layout(res, qlayout)[source]

Transpilers sometimes change the order of qubits, so we need to reorder the bitstrings according to the mapping given by qlayout.

nuqulib.quantum_algorithms.trans_XYloc_str(Xloc, Yloc, Bosonic)[source]

Translate X/Y support lists into a measurement-basis descriptor string.

Parameters:
  • Xloc (list[int]) – Qubit locations with X operators.

  • Yloc (list[int]) – Qubit locations with Y operators.

  • Bosonic (bool) – If True, collapse all-X/all-Y terms to "XX" or "YY" descriptors.

Returns:

Descriptor consumed by QKrylov measurement-circuit builders.

Return type:

str

nuclear_hamiltonian submodule

Nuclear Hamiltonian construction and manipulation module.

This module provides functionality for constructing and manipulating nuclear Hamiltonians for quantum many-body calculations. It supports both valence-space and no-core shell model (NCSM) interactions in the standard snt format used by nuclear physics codes such as KSHELL, NuHamil, and NuclearToolkit.jl.

The module handles: - Two-body and three-body nuclear force interactions - M-scheme and JT-coupled representations - Mapping to fermionic operators and Pauli strings for quantum computing - Matrix element calculations and transformations - Quantum resource estimation for nuclear Hamiltonians

Key Classes: - Hamiltonian: Main class for nuclear Hamiltonian construction and manipulation

class nuqulib.nuclear_hamiltonian.Hamiltonian(fn_NN: str | PathLike, Z: int, N: int, jz: int | None = None, fn_3NF: str | PathLike | None = None, ncsm: bool = False, emax_truncate: int = 5, Qiskit_order: bool = True, verbose: bool = False, e3max: int | None = None, mapping_method: str = 'Jordan-Wigner', single_spiecies: int = 0)[source]

Bases: object

Nuclear Hamiltonian constructor and manager for quantum many-body calculations.

This class provides comprehensive functionality for constructing nuclear Hamiltonians from standard interaction files, supporting both valence-space and no-core shell model (NCSM) calculations. It handles two-body and three-body nuclear forces, performs basis transformations, and maps to quantum operator representations.

The class reads interaction files in the snt format used by nuclear physics codes (KSHELL, NuHamil, NuclearToolkit.jl) and provides methods for: - Single-particle basis management - Matrix element calculations and transformations - M-scheme and JT-coupled representations - Mapping to fermionic operators and Pauli strings - Three-body force handling and optimization

fn_NN

Path to two-body interaction file.

Type:

str

Z

Number of protons.

Type:

int

N

Number of neutrons.

Type:

int

Anum

Mass number (Z + N).

Type:

int

nucleus

Nuclear symbol (e.g., “16O”).

Type:

str

emax

Maximum excitation energy for truncation.

Type:

int

ncsm

Whether using no-core shell model.

Type:

bool

hw

Harmonic oscillator frequency (for NCSM).

Type:

float

n_qubits

Total number of qubits needed.

Type:

int

n_qubits_p

Number of proton qubits.

Type:

int

n_qubits_n

Number of neutron qubits.

Type:

int

msps

M-scheme single-particle states.

Type:

list

v1b

One-body matrix elements.

Type:

list

v2b

Two-body matrix elements.

Type:

list

v3b_pn

Three-body matrix elements in p-n representation.

Type:

dict

v3b_Mscheme

Three-body matrix elements in M-scheme.

Type:

dict

Example

>>> # For valence-space calculation
>>> hamil = Hamiltonian("interaction.snt", Z=8, N=8)
>>> H_dict = hamil.get_mscheme_H(opform=True)
>>>
>>> # For NCSM with three-body forces
>>> hamil = Hamiltonian("ncsm_int.snt", Z=4, N=4,
...                     fn_3NF="3nf.me3j.gz", ncsm=True, emax_truncate=2)

References

    1. Miyagi et al., Eur. Phys. J. A 59, 150 (2023)

  • Standard snt format used by KSHELL and NuclearToolkit.jl for nuclear interactions

delta_jj(morb_a, morb_b)[source]

Check if two orbitals are identical in terms of their n, l, j, and tz quantum numbers.

Parameters:
Returns:

1 if orbitals are identical, 0 otherwise.

Return type:

int

extract_hw()[source]

Extract harmonic oscillator frequency from NCSM interaction filename.

Parses the harmonic oscillator frequency (hw) parameter from the filename of NCSM interaction files, which typically contain “hw{value}” in the name.

Returns:

Harmonic oscillator frequency in MeV.

Return type:

float

Note

Required for NCSM calculations to properly scale kinetic energy terms.

get_3NF_Mscheme(verbose=False)[source]

Convert stored M-scheme 3NF matrix elements to operator strings.

Parameters:

verbose (bool, optional) – If True, print each generated proton and neutron operator string.

Returns:

Mapping (proton_op, neutron_op) pairs to matrix elements. Returns an empty dictionary when no M-scheme 3NF data is available.

Return type:

dict

get_V3_p_n()[source]

Split pure-proton and pure-neutron three-body terms.

Returns:

Proton-only and neutron-only 3NF operator dictionaries extracted from self.Hamildict["V3"].

Return type:

tuple[dict, dict]

get_cG(sps_i, sps_j, J)[source]

Calculate Clebsch-Gordan coefficient for angular momentum coupling.

Computes the Clebsch-Gordan coefficient for coupling two single-particle states to a given total angular momentum J.

Parameters:
  • sps_i (Orbit_nljjztz) – First single-particle state in M-scheme.

  • sps_j (Orbit_nljjztz) – Second single-particle state in M-scheme.

  • J (int) – Total angular momentum quantum number.

Returns:

Clebsch-Gordan coefficient <j1 m1 j2 m2 | J M>.

Return type:

float

get_midx_from_nljjztz(n, l, j, jz, tz)[source]

Find M-scheme single-particle state index from quantum numbers.

Searches for the index of a M-scheme single-particle state with specified quantum numbers in the model space.

Parameters:
  • n (int) – Principal quantum number (radial).

  • l (int) – Orbital angular momentum quantum number.

  • j (int) – Total angular momentum quantum number (twice the actual value).

  • jz (int) – Magnetic quantum number (twice the actual value).

  • tz (int) – Isospin projection (-1 for protons, +1 for neutrons).

Raises:

ValueError – If the requested M-scheme state is not found.

get_mscheme_H(opform=False, verbose=0)[source]

Construct Hamiltonian matrix elements in M-scheme representation.

Transforms the JT-coupled nuclear interaction matrix elements into M-scheme (magnetic quantum number resolved) representation, which is required for quantum computing applications. The method handles both one-body and two-body terms with proper angular momentum recoupling.

Parameters:
  • opform (bool, optional) – If True, returns operators in fermionic

  • False (string form. If)

  • False. (returns matrix element lists. Defaults to)

  • verbose (int | bool, optional) – If truthy, print additional progress and basis information while constructing the M-scheme Hamiltonian.

Returns:

Dictionary containing Hamiltonian components:
  • ”SPE”: Single-particle energy terms (1-body)

  • ”Vpp”: Proton-proton interaction terms (2-body)

  • ”Vnn”: Neutron-neutron interaction terms (2-body)

  • ”Vpn”: Proton-neutron interaction terms (2-body)

Return type:

dict

Note

For opform=True, terms are returned as fermionic operator strings (e.g., “+_0 -_1” for creation/annihilation operators). For opform=False, terms are returned as matrix element lists.

The M-scheme transformation includes proper Clebsch-Gordan coefficients and normalization factors for antisymmetrized matrix elements.

get_mscheme_sps()[source]

Generate M-scheme single-particle states from proton-neutron J-coupled states.

Constructs the M-scheme single-particle states from the JT-coupled orbital basis. This includes expanding each j-orbital into its 2j+1 magnetic substates and organizing them by particle type.

Returns:

A tuple containing:
  • proton_register (range): Qubit indices for proton states

  • neutron_register (range): Qubit indices for neutron states

  • mps (list): List of Orbit_nljjztz objects representing M-scheme states

  • dict_sps2msps (dict): Mapping from JT states to M-scheme state lists

Return type:

tuple

Note

The proton states are placed first in the qubit register, followed by neutron states. This convention is used throughout the quantum mappings.

guess_emax_from_fn(fn_3NF)[source]

Extract the model space information from three-body force filename.

Attempts to parse the e1max/e2max/e3max for one-body, two-body, and three-body sectors from the three-body force filename pattern.

Parameters:

fn_3NF (str) – Path to three-body force file.

Returns:

A tuple containing (e1max_file, e2max_file, e3max_file).

If parsing fails, returns default values based on emax_truncate.

Return type:

tuple

Note

Expects filename pattern with “_ms{e1max}_{e2max}_{e3max}” substring. Falls back to emax_truncate * {1, 2, 3} if parsing fails.

mapping_3NF_Mscheme(method='Jordan-Wigner', filepath: str | PathLike = './tmp', verbose=False)[source]

Map the 3NF matrix elements to the specified scheme.

mapping_opform(mapping_method: str = 'Jordan-Wigner', filepath: str | PathLike = './')[source]

Map nuclear Hamiltonian to qubit operators using specified fermion-to-qubit mapping.

Transforms the fermionic nuclear Hamiltonian into qubit operators using standard fermion-to-qubit mappings such as Jordan-Wigner or Bravyi-Kitaev. Handles the separate proton and neutron subsystems with appropriate mappings.

Parameters:
  • mapping_method (str) – Fermion-to-qubit mapping method (e.g., “Jordan-Wigner”).

  • filepath (str|os.PathLike) – Path to the file where the mapping results will be saved.

  • "HATTMapper". (This is currently used when the mapping_method is)

Returns:

A tuple containing mapped Pauli operators:
  • H_1b_p (SparsePauliOp): One-body terms for protons

  • H_1b_n (SparsePauliOp): One-body terms for neutrons

  • H_jz_p (SparsePauliOp): Jz terms for protons

  • H_jz_n (SparsePauliOp): Jz terms for neutrons

  • H_pp (SparsePauliOp): Proton-proton terms

  • H_nn (SparsePauliOp): Neutron-neutron terms

  • H_pn (SparsePauliOp): Proton-neutron terms

Return type:

tuple

Note

The proton-neutron terms require special handling due to the tensor product structure of the proton-neutron Hilbert space.

neutron_number_constraint(N0, mapping_method, filepath='./tmp')[source]

Build a quadratic neutron-number penalty operator.

Parameters:
  • N0 (int | float) – Target neutron number.

  • mapping_method (str) – Fermion-to-qubit mapping method.

  • filepath (str | os.PathLike, optional) – Base path used by mappings that load or save mapper data.

Returns:

Qubit operator for (N_n - N0)^2.

Return type:

SparsePauliOp

op_dict_T1_permutations(op_dict, aa, bb, cc, dd, J, v)[source]

Add two-body operator terms with proper permutation symmetry.

Adds fermionic operator terms to the operator dictionary, including the necessary permutations for antisymmetrized matrix elements.

Parameters:
  • op_dict (dict) – Dictionary to store operator strings and coefficients.

  • aa (int) – First bra-state index.

  • bb (int) – Second bra-state index.

  • cc (int) – First ket-state index.

  • dd (int) – Second ket-state index.

  • J (int) – Total angular momentum quantum number.

  • v (float) – Matrix element value.

Note

We use the convention of the order of creation and annihilation operators like “+_aa +_bb -_dd -_cc” for the operator string. For identical bra and ket configurations, only one term is added. For different configurations, both direct and exchanged terms are included to ensure proper antisymmetrization of fermionic operators.

perms(i1, i2, i3)[source]

Return signed permutations of three single-particle indices.

Parameters:
  • i1 (int) – First index.

  • i2 (int) – Second index.

  • i3 (int) – Third index.

Returns:

Permutations mapped to their antisymmetry signs.

Return type:

dict[tuple[int, int, int], int]

proton_number_constraint(N0, mapping_method, filepath='./tmp')[source]

Build a quadratic proton-number penalty operator.

Parameters:
  • N0 (int | float) – Target proton number.

  • mapping_method (str) – Fermion-to-qubit mapping method.

  • filepath (str | os.PathLike, optional) – Base path used by mappings that load or save mapper data.

Returns:

Qubit operator for (N_p - N0)^2.

Return type:

SparsePauliOp

read_3NF_readable(verbose=False)[source]

Read the NuHamil 3NF file in readable.text fmt and store JT-coupled 3NF.

In readable.text, V_{3N} is given as a function of a set of quanta, {a,b,c, Jab, Tab, d,e,f, Jde, Tde, Jabc, Tabc, Jabc, Tabc} where a~f are the single particle states having {n,l,j} quanta. Those matrix elements correspond to Eq.(41) in NuHamil paper, T.Miyagi, EPJA (2023)59:150. Then, proton-neutron matrix elements are obtained through the Clebsch-Gordan coefficients as detailed in Eq.(36) of NuHamil paper.

Returns:

A dictionary containing the 3NF matrix elements in the pn-coupled form.

Return type:

dict

read_snt_file(fn_NN_in)[source]

Read nuclear interaction file in snt format.

Parses nuclear interaction files containing single-particle energies and two-body matrix elements in the standard snt format used by nuclear physics codes like KSHELL and NuclearToolkit.jl.

Parameters:

fn_NN_in (str or PathLike) – Path to the snt format interaction file.

Returns:

A tuple containing:
  • nsp_p (int): Number of proton single-particle states

  • nsp_n (int): Number of neutron single-particle states

  • core_p (int): Number of core protons

  • core_n (int): Number of core neutrons

  • single_particle_states (list): List of Orbit_nljtz objects

  • v1b (list): One-body matrix elements as (a, b, value) tuples

  • v2b (list): Two-body matrix elements as (a, b, c, d, J, value) tuples

Return type:

tuple

Note

For NCSM calculations, one-body matrix elements represent kinetic terms \(T_{n,n'}\) with \(\hbar\omega (A-1)/A\) scaling factor.

The snt format supports both valence-space and NCSM interactions. Mass-dependence factors \((A/A_{ref})^p\) are handled for empirical shell model interactions.

Raises:

FileNotFoundError – If the interaction file cannot be found.

separate_proton_and_neutron(im_a: int, im_b: int, im_c: int, im_d: int, im_e: int, im_f: int)[source]

Separate the proton and neutron operators for a given set of indices.

set_mscheme_3NF(verbose=False)[source]

Calculate the 3NF matrix elements in the mscheme.

Returns:

A dictionary containing the 3NF matrix elements in the mscheme.

Return type:

dict

Note

Q. Why do we need to devide the 3NF matrix elements by 36? A. All the permutations of the 3NF matrix elements are explicitly counted below. It would be more efficient to avoid redundancy in the first place. It is left for the future.

Q. Why do we need to multiply the 3NF matrix elements by 9 for pnn and npp? A. At first, proton and neutron operators commute. We are now generating all the permutations not taking into account of this. While transforming the pn JT-coupled form to J-coupled form, we generate something like <pnn|V|npn>, <pnn|V|nnp>, etc. The factor of 9 is to take account of this.

task_3NF_pn(p_str, n_str, coeff_overall, method: str, filepath: str | PathLike)[source]

Map one p-n 3NF term into Pauli-label tuples.

This path is split into two stages: 1) map each sector (with cache), 2) build output tuples from cartesian products.

write_v3b_Mscheme(fn_out: str | PathLike = 'v3b_Mscheme.txt')[source]

Write the 3NF matrix elements in M-scheme to a file.

Parameters:

fn_out (str|os.PathLike) – Output filename to write the 3NF matrix elements.

class nuqulib.nuclear_hamiltonian.JTcoupledOrbitals(emax)[source]

Bases: object

Class to handle JT-coupled orbitals in the model space.

This class manages orbital states in the JT (angular momentum - isospin) coupled basis. It provides methods to add orbitals and convert between different single-particle state representations.

emax

Maximum excitation energy for orbitals.

Type:

int

orbitals

Dictionary mapping orbital indices to Orbit_nlj objects.

Type:

dict

dict_sps2JTorbitals

Mapping from single-particle states to JT orbitals.

Type:

dict

add_orbital(idx, n, l, j)[source]

Add an orbital to the JT-coupled basis.

Parameters:
  • idx (int) – Index identifier for the orbital.

  • n (int) – Principal quantum number (radial).

  • l (int) – Orbital angular momentum quantum number.

  • j (int) – Total angular momentum quantum number (twice the actual value).

get_dict_sps2JTsps(sps_defined_in_NN)[source]

Convert single particle states to JT-coupled states.

Maps the single particle states defined in the nuclear interaction files to the JT-coupled orbital basis for consistent matrix element handling.

Parameters:

sps_defined_in_NN (list) – Single particle states from NN interaction file.

Returns:

Mapping from single-particle state indices to JT orbital indices.

Return type:

dict

Note

The single-particle states are defined when the NN interaction file is loaded. This method ensures compatibility between different representations.

class nuqulib.nuclear_hamiltonian.ReadThBME_me3jgz(single_particle_states: dict, JT_orbitals: dict, Z: int, N: int, filename, e1max: int, e1max_file: int, e2max: int, e3max: int, e3max_file: int, CG_dict_fromHamiltonian: dict, single_species: bool = False, verbose: bool = False)[source]

Bases: object

Read three-body matrix elements from a me3j.gz file.

This class handles the reading and processing of three-body matrix elements stored in a compressed me3j.gz file format.

count_nreads(sps_3b: sps_3Blab, mode='File')[source]

Count the number of elements of floating-point numbers to be read from the interaction file.

In some interaction files, odd indices correspond to proton and even indices correspond to neutron. One should take this into account when processing the data when working in JT-coupled schemes.

In the nested loops, the number of reads is counted by considering possible (allowed) <a, b, c [Jab]|V| d, e, f [Jde]> in isospin space, i.e. a~e are indices for Orbit_nlj. The factors ((twoJCMax - twoJCMin) // 2 + 1) and 5 represents the number of allowed J_3 and T_3. Note that T_3 = 1 case has 4 components (Tab, Tde) = (0, 0), (0, 1), (1, 0), (1, 1), while T_3 = 0 case has only 1 component. The Tindex will be 2 * Tab + Tde + ((T2 - 1) // 2), i.e. 0: (0, 0, 1), 1: (0, 1, 0), 2: (1, 0, 0), 3: (1, 1, 1), 4: (1, 1, 3).

Nreads:

is a list of the head indices relevant to the first index in bra/ket for 3NF matrix elements. This would be helpful if one wants to efficiently access the relevant indices via multi-threading.

get_modelspace(e1max, e1max_file, e2max, e3max, e3max_file)[source]

Construct model-space labels for file and truncated 3NF spaces.

Parameters:
  • e1max (int) – One-body truncation used in the active model space.

  • e1max_file (int) – One-body truncation represented in the input file.

  • e2max (int) – Two-body truncation used in the active model space.

  • e3max (int) – Three-body truncation used in the active model space.

  • e3max_file (int) – Three-body truncation represented in the input file.

Returns:

Model-space descriptor used by the me3j reader.

Return type:

sps_3Blab

loop_over_JT(sps_3b: sps_3Blab, what_you_need='', ME_is_double=True, array_in: ndarray = array([0.]))[source]

Method to loop over JT-coupled three-body matrix elements.

This method iterates over the relevant quantum numbers and constructs the necessary data structures for the three-body matrix elements.

read_me3jgz(filename, count_ME_file)[source]

Read three-body matrix elements from a me3j.gz file.

The ordering of the matrix elements are uniquely determined by the e*max values defining the adopted model space.

nuqulib.nuclear_hamiltonian.RecouplingCG(idx_abc, ja2, jb2, jc2, Jab_in, Jab, J2, dWS) float[source]

Return the angular-momentum recoupling coefficient for a 3NF permutation.

Parameters:
  • idx_abc (int) – Permutation index for the three incoming orbitals.

  • ja2 (int) – Twice the angular momentum of orbital a.

  • jb2 (int) – Twice the angular momentum of orbital b.

  • jc2 (int) – Twice the angular momentum of orbital c.

  • Jab_in (int) – Input coupled angular momentum of a and b.

  • Jab (int) – Output coupled angular momentum of the permuted pair.

  • J2 (int) – Twice the total three-body angular momentum.

  • dWS (prep_dicts_for_WignerSymbols) – Precomputed Wigner-symbol tables.

Returns:

Recoupling coefficient.

Return type:

float

nuqulib.nuclear_hamiltonian.Trans_JT_into_pn(Z: int, N: int, a: int, b: int, c: int, Jab: int, Tab: int, d: int, e: int, f: int, Jde: int, Tde: int, J3: int, T3: int, ME: float, pnME_3NF: dict, JTorbitals: dict, single_particle_states: list, CGdict: dict)[source]

Transforms the JT-coupled 3NF into the proton-neutron (single-particle) one.

Note

Developers should be aware that pnME_3NF dict is destructively updated in this method.

nuqulib.nuclear_hamiltonian.export_encoded_Hamil(Hamil: SparsePauliOp, fn: str | PathLike)[source]

Export the encoded Hamiltonian to a text file.

Parameters:
  • Hamil (SparsePauliOp) – Encoded Hamiltonian to write.

  • fn (str | os.PathLike) – Output text filename.

nuqulib.nuclear_hamiltonian.get_CGs_from_dict(j1, m1, j2, m2, J, M, CG_dict: dict)[source]

Retrieve or calculate Clebsch-Gordan coefficient.

Gets a Clebsch-Gordan coefficient from cache or calculates and stores it. This provides efficient access to frequently used CG coefficients in nuclear many-body calculations.

Parameters:
  • j1 (int) – First angular momentum (twice the actual value).

  • m1 (int) – First magnetic quantum number (twice the actual value).

  • j2 (int) – Second angular momentum (twice the actual value).

  • m2 (int) – Second magnetic quantum number (twice the actual value).

  • J (int) – Coupled angular momentum (twice the actual value).

  • M (int) – Coupled magnetic quantum number (twice the actual value).

  • CG_dict (dict) – Cache dictionary for storing calculated coefficients.

Returns:

Clebsch-Gordan coefficient <j1 m1 j2 m2 | J M>.

Return type:

float

Note

The coefficient is calculated using sympy if not found in cache. All angular momenta are provided as integers (twice their actual values).

nuqulib.nuclear_hamiltonian.get_Hamiltonian(fn_NN, Z: int, N: int, fn_3NF=None, emax: int = 20, e3max: int = 0, ncsm=False, mapping_method='Jordan-Wigner', single_spiecies=0)[source]

Get nuclear Hamiltonian from interaction files.

This is a convenience wrapper function that constructs a Hamiltonian object and returns the mapped quantum operators ready for use in quantum algorithms. It handles both two-body and three-body nuclear interactions.

Parameters:
  • fn_NN (str) – Path to the snt format interaction file.

  • Z (int) – Number of protons in the nucleus.

  • N (int) – Number of neutrons in the nucleus.

  • fn_3NF (str, optional) – Path to three-body force file. Defaults to “”.

  • emax (int, optional) – Maximum excitation energy for truncation. Defaults to 20.

  • e3max (int, optional) – Maximum excitation energy for three-body forces. Defaults to 0.

  • ncsm (bool, optional) – Whether using no-core shell model. Defaults to False.

  • mapping_method (str, optional) – The mapping method for converting fermionic operators to Pauli operators. Defaults to “Jordan-Wigner”.

  • single_spiecies (int, optional) – Select a single-species calculation mode. 0 keeps both proton and neutron states, 1 keeps protons only, and 2 keeps neutrons only.

Returns:

A tuple containing:
  • hamil (Hamiltonian): The Hamiltonian object

  • H_mapped (SparsePauliOp): Mapped Hamiltonian as Pauli operators

  • proton_qubits (list): Qubit indices corresponding to proton states

  • neutron_qubits (list): Qubit indices corresponding to neutron states

Return type:

tuple

Example

>>> hamil, H_mapped, p_qubits, n_qubits = get_Hamiltonian("interaction.snt", 8, 8)
>>> print(f"Total qubits: {len(p_qubits) + len(n_qubits)}")
nuqulib.nuclear_hamiltonian.get_canonical_order_6j(j1: int, j2: int, j3: int, j4: int, j5: int, j6: int) int[source]

Return the canonical column order code for a 6j-symbol key.

nuqulib.nuclear_hamiltonian.get_key6j_sym(j1: int, j3: int, j5: int, j2: int, j4: int, j6: int) int[source]

Return a symmetry-canonical packed key for a Wigner 6j symbol.

nuqulib.nuclear_hamiltonian.get_nkey6(a, b, c, d, e, f)[source]

Pack six non-negative integer labels into one integer key.

nuqulib.nuclear_hamiltonian.get_nkey6_shift(a, b, c, d, e, f, int_shift=3)[source]

Pack six possibly negative integer labels into one shifted integer key.

Parameters:
  • a (int) – First label.

  • b (int) – Second label.

  • c (int) – Third label.

  • d (int) – Fourth label.

  • e (int) – Fifth label.

  • f (int) – Sixth label.

  • int_shift (int, optional) – Shift added before packing.

Returns:

Packed key.

Return type:

int

nuqulib.nuclear_hamiltonian.hat(a)[source]

Angular momentum hat notation: sqrt(2*j + 1).

Parameters:

a (int or float) – Angular momentum quantum number (can be half-integer).

Returns:

sqrt(2*a + 1), commonly used in nuclear physics calculations.

Return type:

float

Note

This is the standard “hat” notation used in nuclear physics and atomic physics for angular momentum algebra.

nuqulib.nuclear_hamiltonian.j_col_score(j1: int, j2: int) int[source]

Score one 6j-symbol column for canonical ordering.

nuqulib.nuclear_hamiltonian.permutation_parity(lst)[source]

Calculate the parity (even/odd) of a permutation.

Determines whether a permutation is even (0) or odd (1) by counting the number of inversions in the list.

Parameters:

lst (list) – A permutation as a list of elements.

Returns:

0 for even permutation, 1 for odd permutation.

Return type:

int

Note

Used in nuclear physics calculations where antisymmetrization requires tracking permutation signs.

class nuqulib.nuclear_hamiltonian.prep_dicts_for_WignerSymbols(emax)[source]

Bases: object

Precompute Clebsch-Gordan and Wigner 6j lookup tables.

prep_d6j_int(emax, jmax_in)[source]

Precompute integer-grid 6j symbols for 3NF transformations.

prep_d6j_lj(jmax2)[source]

Precompute 6j symbols for l-j recoupling up to jmax2.

prep_dcg_spin()[source]

Precompute spin Clebsch-Gordan coefficients used in 3NF recoupling.

nuqulib.nuclear_hamiltonian.process_op(args)[source]

Process a single proton-neutron operator term for quantum mapping.

This function processes individual terms from the three-body force dictionary, mapping them from fermionic operators to Pauli strings using the specified fermion-to-qubit mapping. Designed for parallel processing.

Parameters:

args (tuple) – A tuple containing: - p_str (str): Proton fermionic operator string - n_str (str): Neutron fermionic operator string - coeff_overall (complex): Overall coefficient - n_qubits_p (int): Number of proton qubits - n_qubits_n (int): Number of neutron qubits - method (str): Fermion-to-qubit mapping method

Returns:

List of (pauli_label, coefficient) tuples for this operator.

Return type:

list

Note

This function is designed to be used with multiprocessing for efficient parallel processing of large three-body force matrices.

nuqulib.nuclear_hamiltonian.removing_redundant_ops(op_list)[source]

Remove redundant Pauli operators by combining terms with same labels.

Consolidates a list of Pauli operators by combining coefficients for operators with identical Pauli labels, removing near-zero terms.

Parameters:

op_list (list) – List of (pauli_label, coefficient) tuples.

Returns:

Consolidated Pauli operator with redundant terms removed.

Return type:

SparsePauliOp

Note

Terms with coefficients smaller than 1e-16 are discarded to avoid numerical precision issues in quantum simulations.

nuqulib.nuclear_hamiltonian.removing_redundant_terms(ops: SparsePauliOp, verbose=False)[source]

Remove redundant terms from existing SparsePauliOp.

Similar to removing_redundant_ops but works on existing SparsePauliOp objects. Combines terms with identical Pauli labels and removes negligible coefficients.

Parameters:
  • ops (SparsePauliOp) – Input Pauli operator.

  • verbose (bool, optional) – If True, print the retained Pauli strings and coefficients after consolidation.

Returns:

Cleaned Pauli operator with redundant terms combined.

Return type:

SparsePauliOp

Note

This is used for final cleanup of mapped Hamiltonians before quantum algorithm applications. Without this step, imaginary terms may persist and affect the results.

nuqulib.nuclear_hamiltonian.set_op_list_from_op_dict_3b(op_dict_3b, n_qubits_p, n_qubits_n, method='Jordan-Wigner')[source]

Map three-body force dictionary to Pauli operator list using parallel processing.

Efficiently converts a dictionary of proton-neutron operator pairs to a list of Pauli operators using the specified fermion-to-qubit mapping. Employs multiprocessing for computational efficiency on large three-body matrices.

Parameters:
  • op_dict_3b (dict) – Dictionary of (proton_op, neutron_op) -> coefficient pairs.

  • n_qubits_p (int) – Number of proton qubits.

  • n_qubits_n (int) – Number of neutron qubits.

  • method (str, optional) – Fermion-to-qubit mapping method. Defaults to “Jordan-Wigner”.

Returns:

List of (pauli_label, coefficient) tuples representing the mapped operators.

Return type:

list

Note

Uses multiprocessing with fork context for efficient parallel processing. The number of processes is set to cpu_count - 2 to avoid system overload.

nuqulib.nuclear_hamiltonian.sort_3_orbits(a_in, b_in, c_in)[source]

Sort three orbital labels into descending canonical order.

Even labels are shifted down by one before sorting, matching the indexing convention used by the three-body matrix-element reader.

Parameters:
  • a_in (int) – First orbital label.

  • b_in (int) – Second orbital label.

  • c_in (int) – Third orbital label.

Returns:

Sorted labels and the permutation index.

Return type:

tuple[int, int, int, int]

class nuqulib.nuclear_hamiltonian.sps_3Blab(e1max: int, e1max_file: int, e2max_file: int, e3max: int, e3max_file: int, norbits_ms: int, norbits_file: int, sps, sps_file)[source]

Bases: object

Single-particle state manager for three-body matrix element calculations.

This class manages single-particle states and model space parameters needed for three-body force calculations. It stores both the truncated model space parameters and the parameters from the original interaction files.

e1max

Maximum single-particle excitation energy for model space.

Type:

int

e1max_file

Maximum single-particle excitation from original file.

Type:

int

e2max_file

Maximum two-body excitation from original file.

Type:

int

e3max

Maximum three-body excitation energy for model space.

Type:

int

e3max_file

Maximum three-body excitation from original file.

Type:

int

norbits_ms

Number of orbitals in model space.

Type:

int

norbits_file

Number of orbitals in original file.

Type:

int

sps

Single-particle states for model space.

Type:

dict

sps_file

Single-particle states from original file.

Type:

dict

nuqulib.nuclear_hamiltonian.sum_over_J(Hamil)[source]

Sum two-body matrix elements over total angular momentum J.

Combines two-body matrix elements with the same orbital indices but different total angular momentum quantum numbers, which is needed for certain nuclear structure calculations.

Parameters:

Hamil (dict) – Dictionary of Hamiltonian components with matrix element lists.

Returns:

Modified Hamiltonian dictionary with J-summed matrix elements.

Return type:

dict

Note

Single-particle energies (SPE) are passed through unchanged. Two-body terms (Vpp, Vnn, Vpn) are combined over J quantum numbers.

nuqulib.nuclear_hamiltonian.tri_check(a, b, c)[source]

Return True if three angular momenta satisfy triangle inequalities.

nuqulib.nuclear_hamiltonian.truncate_v3bme(dim_v3bme_file: int, dim_v3bme_ms: int, sps_3b: sps_3Blab, ThBME, dWS, dict_idxheadThBME_file, dict_idxheadThBME_MS)[source]

Truncate JT-coupled 3NF matrix elements with specified model space truncations.

nuqulib.nuclear_hamiltonian.unhash_key6j(i)[source]

Unpack a key produced by get_nkey6() into six labels.

nuqulib.nuclear_hamiltonian.valid_check(ea, eb, ec, ed, ee, ef, e1max, e2max, e3max)[source]

Check if three-body matrix element satisfies model space truncation limits.

Validates that a three-body matrix element with excitation energies ea-ef is within the specified model space truncation parameters.

Parameters:
  • ea (int) – Excitation energy of the first bra orbital.

  • eb (int) – Excitation energy of the second bra orbital.

  • ec (int) – Excitation energy of the third bra orbital.

  • ed (int) – Excitation energy of the first ket orbital.

  • ee (int) – Excitation energy of the second ket orbital.

  • ef (int) – Excitation energy of the third ket orbital.

  • e1max (int) – Maximum single-particle excitation energy.

  • e2max (int) – Maximum two-body excitation energy.

  • e3max (int) – Maximum three-body excitation energy.

Returns:

True if the matrix element is within model space limits.

Return type:

bool

Note

This enforces the standard nuclear physics truncation scheme where single-particle, two-body, and three-body excitations are separately limited to maintain computational tractability.

vqe_example_pennylane submodule

VQE example implementation using PennyLane.

This module provides a Variational Quantum Eigensolver (VQE) example using PennyLane for nuclear quantum simulations. It includes functions for setting up nuclear Hamiltonians in OpenFermion format and running VQE optimization.

nuqulib.vqe_example_pennylane.define_Hamil_in_OpenFermion(Hamil_dict, proton_number=1, neutron_number=1)[source]

Convert Hamiltonian dictionary to OpenFermion format.

This function takes a nuclear Hamiltonian in dictionary format and converts it to OpenFermion FermionOperator objects for use in quantum algorithms.

Parameters:
  • Hamil_dict (dict) – Dictionary containing Hamiltonian matrix elements with keys: - “SPE”: Single-particle energies - “Vpp”: Proton-proton interactions - “Vnn”: Neutron-neutron interactions - “Vpn”: Proton-neutron interactions

  • proton_number (int, optional) – Number of protons. Defaults to 1.

  • neutron_number (int, optional) – Number of neutrons. Defaults to 1.

Returns:

Tuple containing:
  • H1b (FermionOperator): One-body Hamiltonian terms.

  • Hpp (FermionOperator): Proton-proton interaction terms.

  • Hnn (FermionOperator): Neutron-neutron interaction terms.

  • Hpn (FermionOperator): Proton-neutron interaction terms.

Return type:

tuple

Note

Orbital indices in the input dictionary are assumed to be 1-based and are converted to 0-based for OpenFermion compatibility.

nuqulib.vqe_example_pennylane.vqe_example_pennylane(Hdict, proton_number, neutron_number, n_qubits_p, n_qubits_n, using_chs, mapping_method='JordanWigner', save_plot=False)[source]

Run VQE optimization for nuclear systems using PennyLane.

This function performs a Variational Quantum Eigensolver (VQE) calculation for nuclear systems using PennyLane and OpenFermion. It supports both Jordan-Wigner and Bravyi-Kitaev fermion-to-qubit mappings.

Parameters:
  • Hdict (dict) – Dictionary containing Hamiltonian matrix elements.

  • proton_number (int) – Number of protons in the system.

  • neutron_number (int) – Number of neutrons in the system.

  • n_qubits_p (int) – Number of qubits for proton sector.

  • n_qubits_n (int) – Number of qubits for neutron sector.

  • using_chs (list) – List of interaction channels to include (e.g., [“1b”, “pp”, “nn”, “pn”]).

  • mapping_method (str, optional) – Fermion-to-qubit mapping method. Choices are “JordanWigner” or “BravyiKitaev”. Defaults to “JordanWigner”.

  • save_plot (bool, optional) – If True, save the optimization history plot.

Returns:

Tuple containing:
  • params_opt (numpy.ndarray): Optimized variational parameters.

  • min_energy (float): Minimum energy found during optimization.

Return type:

tuple

Note

The ansatz uses Hartree-Fock initial state followed by Givens rotations for single excitations. The optimization uses Adam optimizer with 150 iterations.