Building a VQE with PennyLane: A Practical Guide

The Variational Quantum Eigensolver (VQE) is one of the most important near-term quantum algorithms. It finds the ground state energy of a molecule or material — a calculation that's exponentially hard for classical computers but tractable on NISQ devices. This guide walks through a complete PennyLane VQE implementation.

What VQE Does

VQE finds the minimum eigenvalue of a Hamiltonian H (usually representing a molecule's energy). It works by:

1. Preparing a parameterized trial state |ψ(θ)⟩ using a quantum circuit
2. Measuring the expectation value ⟨ψ(θ)|H|ψ(θ)⟩ on the QPU
3. Using a classical optimizer to update θ to minimize the energy
4. Repeating until convergence

The variational principle guarantees ⟨ψ(θ)|H|ψ(θ)⟩ ≥ E₀ for any state — so minimizing this quantity gives an upper bound on the true ground state energy E₀.

Setting Up

Install PennyLane with its chemistry plugin:

pip install pennylane pennylane-qchem

For hydrogen molecule (H₂) — the classic VQE benchmark — we need two electrons and four spin-orbitals (4 qubits):

import pennylane as qml
from pennylane import numpy as np
import pennylane.qchem as qchem

# H2 at equilibrium bond length (Angstrom)
symbols = ["H", "H"]
coordinates = np.array([[0.0, 0.0, -0.6614], [0.0, 0.0, 0.6614]])

# Build the qubit Hamiltonian
H, qubits = qchem.molecular_hamiltonian(
    symbols,
    coordinates,
    basis="sto-3g"
)
print(f"Hamiltonian: {len(H.ops)} terms, {qubits} qubits")
# Hamiltonian: 15 terms, 4 qubits

Defining the Ansatz

The ansatz is the parameterized circuit that prepares the trial state. For chemistry problems, the UCCSD (Unitary Coupled-Cluster Singles and Doubles) ansatz is standard:

# Get UCCSD circuit parameters
electrons = 2
singles, doubles = qchem.excitations(electrons, qubits)
s_wires, d_wires = qchem.excitations_to_wires(singles, doubles, wires=range(qubits))

# Initial Hartree-Fock state (reference state)
hf_state = qchem.hf_state(electrons, qubits)

dev = qml.device("default.qubit", wires=qubits)

@qml.qnode(dev)
def circuit(weights, wires, s_wires=[], d_wires=[], hf_state=hf_state):
    # Prepare HF reference state
    qml.BasisState(hf_state, wires=wires)

    # Apply UCCSD excitations
    qml.UCCSD(weights, wires, s_wires=s_wires, d_wires=d_wires, init_state=hf_state)

    return qml.expval(H)

Running the VQE Optimization

With PennyLane's automatic differentiation, we can use gradient-based optimizers directly:

# Initial parameters (all zeros = Hartree-Fock state)
init_params = np.zeros(len(singles) + len(doubles), requires_grad=True)

# Adam optimizer (works well for VQE)
opt = qml.AdamOptimizer(stepsize=0.4)

# Optimization loop
energy_history = []
params = init_params.copy()

for step in range(200):
    params, energy = opt.step_and_cost(
        lambda p: circuit(p, range(qubits), s_wires=s_wires, d_wires=d_wires),
        params
    )
    energy_history.append(energy)

    if step % 20 == 0:
        print(f"Step {step:3d}: E = {energy:.6f} Ha")

print(f"\nVQE ground state energy: {energy:.6f} Ha")
print(f"Reference (exact): -1.136189 Ha")

Typical output:

Step   0: E = -1.117498 Ha
Step  20: E = -1.133254 Ha
Step  40: E = -1.135901 Ha
Step  60: E = -1.136140 Ha
...
VQE ground state energy: -1.136174 Ha
Reference (exact): -1.136189 Ha

VQE reaches within ~0.015 mHa of the exact energy — chemical accuracy for H₂.

Using a Gradient-Free Optimizer

For noisy hardware, gradient-free optimizers like COBYLA or SPSA are often better since hardware gradients are noisy:

from scipy.optimize import minimize

# Objective function (no gradient needed)
def objective(params):
    return float(circuit(params, range(qubits), s_wires=s_wires, d_wires=d_wires))

result = minimize(
    objective,
    x0=init_params,
    method="COBYLA",
    options={"maxiter": 500, "rhobeg": 0.1}
)

print(f"COBYLA energy: {result.fun:.6f} Ha")

Running VQE with HLQuantum

HLQuantum includes a built-in VQE implementation that works across all backends:

import hlquantum as hlq

# Define the Hamiltonian in HLQuantum's format
H = hlq.hamiltonians.h2_molecule(bond_length=1.32)

# Run VQE on any backend
result = hlq.algorithms.vqe(
    hamiltonian=H,
    ansatz="uccsd",
    optimizer="adam",
    max_iterations=200,
    backend="pennylane",   # or "qiskit", "cudaq"
)

print(f"Ground state energy: {result.energy:.6f} Ha")
print(f"Optimal parameters: {result.params}")
print(f"Converged in {result.iterations} iterations")

Tips for Real Hardware

When running VQE on real QPUs (IBM Quantum, IonQ), several additional considerations apply:

Use fewer shots per step. 1000 shots per optimization step is usually enough for the gradient estimate. Don't use 10,000 shots at every step — it wastes QPU time.

Start with shallow circuits. Fewer CNOT gates = less noise. For hardware, consider hardware-efficient ansatz circuits rather than UCCSD.

Enable error mitigation. HLQuantum's error_mitigation="zne" applies Zero Noise Extrapolation, which can significantly improve results on noisy hardware:

result = hlq.run(vqe_circuit, backend="qiskit", device="ibm_sherbrooke",
                 error_mitigation="zne", shots=2048)

Check the full PennyLane guide and the HLQuantum algorithms reference for more details on quantum chemistry simulations.