Atom Computing Nuclear Physics Simulation with 1000-Qubit Neutral Atom System

Atom Computing’s Phoenix processor achieves the 1180-qubit milestone by trapping ytterbium-171 atoms in a 2D optical lattice. Unlike rubidium-based systems, ytterbium-171 has a nuclear spin of one-half, making it a natural qubit with coherence times measured in seconds rather than milliseconds. The qubit is encoded in two nuclear spin states of the electronic ground manifold, isolated from magnetic field noise by the absence of electronic angular momentum. Two-qubit gates are implemented via the Rydberg CZ mechanism: a laser pulse temporarily excites both atoms to a Rydberg state, where the strong van der Waals interaction between Rydberg electrons creates a conditional phase, then the atoms are de-excited. Gate fidelities on Phoenix reach approximately 99.5% for two-qubit CZ gates in the best-performing atom pairs, though maintaining this fidelity uniformly across 1180 sites as the system scales remains the central engineering challenge.

The nuclear shell model describes nucleons (protons and neutrons) occupying discrete energy shells around the nucleus, analogous to electron shells in atomic physics. For Carbon-12 and Oxygen-16, the relevant valence nucleons occupy the p-shell and sd-shell respectively. The many-body Hamiltonian is:

H = sum_i epsilon_i * a_i^dagger a_i + (1/4) * sum_ijkl V_ijkl * a_i^dagger a_j^dagger a_k a_l

where epsilon_i are single-particle energies and V_ijkl are two-body matrix elements from the USD or USDA interaction. Mapping to qubits uses the Jordan-Wigner or Bravyi-Kitaev transformation, converting fermionic creation/annihilation operators to Pauli strings. For Carbon-12 in the p-shell (6 single-particle states for protons, 6 for neutrons), the qubit count is manageable at 12 logical qubits, but realistic sd-shell calculations for Oxygen-16 require 24 qubits with circuits of several hundred gates.

from qiskit import QuantumCircuit, transpile
from qiskit.quantum_info import SparsePauliOp
import numpy as np

# Nuclear shell model Hamiltonian for Carbon-12 (p-shell, Jordan-Wigner)
# Simplified 6-qubit model for demonstration (3 proton + 3 neutron p-shell orbitals)

def p_shell_hamiltonian() -> SparsePauliOp:
    """
    Construct a simplified p-shell Hamiltonian for C-12.
    Single-particle energies (MeV) for p1/2 and p3/2 orbitals.
    Returns SparsePauliOp over 6 qubits (Jordan-Wigner encoded).
    """
    # Single-particle energy differences (Cohen-Kurath interaction)
    epsilon = [-4.2, -4.2, -4.2, -1.0, -1.0, -1.0]  # MeV, p3/2 then p1/2

    terms = []
    coeffs = []

    # Diagonal single-particle terms: epsilon_i * (I - Z_i) / 2
    for i, e in enumerate(epsilon):
        pauli_str = ["I"] * 6
        pauli_str[i] = "Z"
        terms.append("".join(reversed(pauli_str)))
        coeffs.append(-e / 2)
        terms.append("I" * 6)
        coeffs.append(e / 2)

    # Two-body terms (placeholder for Cohen-Kurath matrix elements)
    # In practice, load V_ijkl from nuclear data tables
    for i in range(5):
        for j in range(i + 1, 6):
            pauli_str = ["I"] * 6
            pauli_str[i] = "Z"
            pauli_str[j] = "Z"
            terms.append("".join(reversed(pauli_str)))
            coeffs.append(-0.5)  # placeholder coupling

    return SparsePauliOp(terms, coeffs=coeffs).simplify()

H = p_shell_hamiltonian()
print(f"Hamiltonian has {len(H)} Pauli terms over {H.num_qubits} qubits")

# VQE ansatz: hardware-efficient layered RZ-CNOT circuit
def shell_model_ansatz(n_qubits: int, depth: int) -> QuantumCircuit:
    qc = QuantumCircuit(n_qubits)
    # Hartree-Fock reference: fill lowest orbitals for C-12 (4 protons, 4 neutrons in p-shell)
    for i in range(4):
        qc.x(i)
    # Variational layers
    params = []
    from qiskit.circuit import ParameterVector
    theta = ParameterVector("theta", depth * n_qubits)
    idx = 0
    for _ in range(depth):
        for q in range(n_qubits):
            qc.ry(theta[idx], q)
            idx += 1
        for q in range(n_qubits - 1):
            qc.cx(q, q + 1)  # Rydberg CZ compiled to CNOT
        qc.cx(n_qubits - 1, 0)  # periodic boundary
    return qc

ansatz = shell_model_ansatz(6, depth=3)
print(f"Ansatz depth: {ansatz.depth()}, gate count: {ansatz.size()}")

# Transpile for Atom Computing (CZ native gate, OpenQASM export)
from qiskit.circuit.library import CZGate
transpiled = transpile(ansatz, basis_gates=["rz", "ry", "cz"], optimization_level=3)
print(f"Transpiled depth: {transpiled.depth()}")
qasm_str = transpiled.qasm()
print(qasm_str[:300])

The 1000-qubit milestone matters not just as a headline number but because it opens a new regime of quantum simulation. Nuclear physics calculations that are intractable for exact diagonalization (Hilbert space dimension grows exponentially with valence nucleon count) and expensive for coupled-cluster methods become accessible targets for variational quantum eigensolvers running on Phoenix. Benchmark results on Phoenix showed that entangled states involving all 1180 qubits could be prepared and measured, with coherence times on the order of 10 seconds allowing 100+ sequential gate operations before decoherence dominates. The primary bottleneck at scale is not coherence but two-qubit gate fidelity variation across the array: outer-edge atoms experience slightly different beam intensities, leading to fidelity gradients that error mitigation techniques must compensate. Future work focuses on mid-circuit measurement and feed-forward operations to enable quantum error correction protocols on Phoenix-scale devices.