Crossing the 12,000-atom barrier with heterogeneous quantum-classical supercomputing: quantum chemistry of protein-ligand complexes

Ab initio wavefunction methods provide accurate molecular simulations but their computational scaling restricts applications to small systems. We develop a workflow combining quantum embedding to decompose a molecule into fragments with a heterogeneous quantum-classical (HQC) method to simulate fragments. We sample fragment electronic configurations on two 156-qubit quantum processors (ibm$\_$cleveland, ibm$\_$kobe), using up to 94 qubits, running 9,200 circuits for over 100 hours, collecting $1.3 \cdot 10^9$ measurement outcomes - the most resource-intensive HQC computation for quantum chemistry to date. We compute fragment wavefunctions via optimized subspace diagonalization on two supercomputers (Fugaku, Miyabi-G), achieving 72.5$\%$ parallel efficiency with scalable distributed linear algebra kernels. We simulate two protein-ligand complexes spanning dispersion- and electrostatics-dominated regimes (11,608 and 12,635 atoms), demonstrate $>40\times$ increase in system size and up to $210\times$ improvement in accuracy over the previous state-of-the-art, with HQC matching coupled-cluster (CCSD) accuracy in fragment energies, and establish a scalable pathway for systematically improvable biomolecular simulations.