Dow Chemical Quantum Simulation for Polymer Design

The Problem

Designing specialty polymers for electronics packaging requires accurate prediction of glass transition temperatures (Tg) and mechanical properties such as tensile modulus. These properties depend on subtle quantum mechanical effects: the conformational energy landscape of polymer chains, van der Waals interactions between chains, and electronic polarizability at interfaces. Classical molecular dynamics can simulate large polymer systems, but it relies on force fields that are parameterized from experiment or from lower-level quantum chemistry. For novel monomers with no experimental data, this parameterization is unreliable.

Density functional theory (DFT) improves accuracy, but standard DFT functionals systematically underestimate dispersion interactions, which are critical for predicting Tg in amorphous polymers. The gold-standard coupled cluster method CCSD(T) gives the right answer but scales as O(N^7) in system size, making it impractical for systems larger than a dozen heavy atoms. Dow needed a path to CCSD(T)-quality binding energies for polymer building blocks without the classical computational cost.

The Quantum Approach

Dow joined the IBM Quantum Network and used Qiskit Nature with PySCF as the classical integral backend. The strategy was active space approximation: instead of treating all electrons in a polymer monomer quantum mechanically, identify the frontier orbitals most relevant to chain-chain binding and restrict the quantum simulation to that subspace. For a representative epoxy monomer used in electronics encapsulants, this yielded an active space of 8 electrons in 6 orbitals, mapping to a 12-qubit Hamiltonian via Jordan-Wigner transformation.

VQE with a UCCSD ansatz was used to find the ground state energy of the active space Hamiltonian. The remaining electrons were treated at the Hartree-Fock level, with their contribution added back classically. This hybrid classical-quantum workflow is called CASSCF-VQE (complete active space self-consistent field with VQE as the eigensolver).

from qiskit_nature.second_q.drivers import PySCFDriver
from qiskit_nature.second_q.transformers import ActiveSpaceTransformer
from qiskit_nature.second_q.mappers import JordanWignerMapper
from qiskit_nature.second_q.circuit.library import UCCSD
from qiskit_algorithms.minimum_eigensolvers import VQE
from qiskit_algorithms.optimizers import L_BFGS_B
from qiskit.primitives import Estimator
import numpy as np

# Epoxy monomer fragment (bisphenol-A diglycidyl ether simplified unit)
driver = PySCFDriver(
    atom="""
    C   0.000   0.000   0.000;
    C   1.400   0.000   0.000;
    C   2.100   1.212   0.000;
    C   1.400   2.424   0.000;
    O   2.050   3.600   0.000;
    C   1.750   4.800   0.000;
    O   0.400   4.850   0.000;
    H   2.800   1.212   0.000
    """,
    basis="cc-pVDZ",
    charge=0,
    spin=0,
)

problem = driver.run()

# Active space: 8 electrons in 6 spatial orbitals (12 qubits after JW mapping)
transformer = ActiveSpaceTransformer(
    num_electrons=8,
    num_spatial_orbitals=6,
)
reduced_problem = transformer.transform(problem)

mapper = JordanWignerMapper()
hamiltonian = mapper.map(reduced_problem.second_q_ops()[0])

print(f"Active space Hamiltonian: {hamiltonian.num_qubits} qubits")
print(f"Pauli terms: {len(hamiltonian)}")

# UCCSD ansatz: captures single and double excitations from HF reference
ansatz = UCCSD(
    qubit_mapper=mapper,
    num_particles=reduced_problem.num_particles,
    num_spatial_orbitals=reduced_problem.num_spatial_orbitals,
)

optimizer = L_BFGS_B(maxiter=600)
estimator = Estimator()

vqe = VQE(estimator=estimator, ansatz=ansatz, optimizer=optimizer)
result = vqe.compute_minimum_eigenvalue(hamiltonian)

vqe_energy = result.eigenvalue.real
print(f"VQE ground state energy: {vqe_energy:.6f} Hartree")

Comparison to DFT and the Scaling Challenge

The VQE result for the monomer active space came within 2 kcal/mol of the CCSD(T) reference, a significant improvement over B3LYP-D3 DFT which underestimated the binding energy by 4-7 kcal/mol on the same fragment. For parameterizing force fields, the 2 kcal/mol accuracy is sufficient to distinguish polymers that will have meaningfully different Tg values.

The fundamental scaling challenge is that the active space approximation becomes less reliable as the polymer chain grows. A single monomer unit might require 12 qubits; a dimer requires roughly 24; a representative oligomer with 10 repeat units would need several hundred logical qubits. Crucially, these must be logical qubits with error correction, not the physical qubits available on IBM Eagle. Current estimates place the crossover to fault-tolerant hardware needed for full-chain accuracy in the early 2030s.

Hybrid Classical-Quantum Workflow

Dow’s production workflow combines the quantum eigensolver with classical embedding: the quantum computer solves the active space problem for a small fragment, and those energies are used to reparameterize a classical force field for full-chain molecular dynamics. This means the quantum calculation runs once per monomer type, not once per simulation step, making the approach practical even with current quantum hardware access constraints.

The workflow also incorporated error mitigation via zero-noise extrapolation on the IBM Eagle processor, reducing the impact of two-qubit gate noise on the energy estimates. Even with mitigation, simulation on Qiskit Aer remained the primary workhorse for production runs, with real hardware used for validation of selected configurations.

Learn more: Qiskit Nature Reference | VQE Algorithm Guide