Hamiltonian Evolution for a Water Molecule
This tutorial demonstrates the ability of the Classiq synthesis engine to reduce depth and cx-counts in approximated quantum functions for Hamiltonian evolution, focusing on Suzuki-Trotter (ST) and qDRIFT (qD) product formulas and their controlled operations. In addition, it is compared to the equivalent quantum examples in Qiskit.
The demonstration is for the Hamiltonian of a water molecule, which has 551 terms and a dimension of 12 qubits.
Define a molecule and get the Hamiltonian as a list of Pauli strings and coefficients:
from classiq.applications.chemistry import Molecule, MoleculeProblem
molecule_H2O = Molecule(
atoms=[("O", (0.0, 0.0, 0.0)), ("H", (0, 0.586, 0.757)), ("H", (0, 0.586, -0.757))]
)
gs_problem = MoleculeProblem(
molecule=molecule_H2O,
basis="sto3g",
mapping="jordan_wigner",
z2_symmetries=False,
freeze_core=True,
)
pauli_list = gs_problem.generate_hamiltonian().pauli_list
Transfer the Hamiltonian to the relevant structure in Classiq:
from typing import cast
from classiq import Pauli, PauliTerm
my_list = {"I": Pauli.I, "X": Pauli.X, "Y": Pauli.Y, "Z": Pauli.Z}
def pauli_str_to_enums(pauli):
return [my_list[s] for s in pauli]
def pauli_list_to_hamiltonian(pauli_list):
return [
PauliTerm(
pauli=pauli_str_to_enums(pauli), coefficient=cast(complex, coeff).real
)
for pauli, coeff in pauli_list
]
hamiltonian = pauli_list_to_hamiltonian(pauli_list)
These cases are examined:
ORDERS_for_ST = [1, 2, 4]
REPETITIONS_for_ST = [6, 4, 1]
N_QDS_for_qDRIFT = [1000, 2000]
classiq_depths = []
classiq_cx_counts = []
# transpilation_options = {"classiq": "custom", "qiskit": 3}
transpilation_options = {"classiq": "auto optimize", "qiskit": 1}
1. Approximating with Suzuki-Trotter formulas
from classiq import (
CustomHardwareSettings,
Preferences,
QArray,
QuantumProgram,
allocate,
create_model,
qfunc,
set_preferences,
suzuki_trotter,
synthesize,
)
preferences = Preferences(
custom_hardware_settings=CustomHardwareSettings(basis_gates=["cx", "u"]),
transpilation_option=transpilation_options["classiq"],
)
for k in range(len(ORDERS_for_ST)):
@qfunc
def main() -> None:
qbv = QArray("qbv")
allocate(len(pauli_list[0][0]), qbv)
suzuki_trotter(
pauli_operator=hamiltonian,
evolution_coefficient=1,
order=ORDERS_for_ST[k],
repetitions=REPETITIONS_for_ST[k],
qbv=qbv,
)
qmod = create_model(main)
qmod = set_preferences(qmod, preferences=preferences)
qprog = synthesize(qmod)
circuit = QuantumProgram.from_qprog(qprog)
classiq_depths.append(circuit.transpiled_circuit.depth)
classiq_cx_counts.append(circuit.transpiled_circuit.count_ops["cx"])
2. Implementing Controlled Hamiltonian Dynamics
from classiq import QBit, control
for k in range(2):
@qfunc
def main() -> None:
qbv = QArray("qbv")
allocate(len(pauli_list[0][0]), qbv)
ctrl = QBit("ctrl")
allocate(1, ctrl)
control(
ctrl=ctrl,
stmt_block=lambda: suzuki_trotter(
pauli_operator=hamiltonian,
evolution_coefficient=1,
order=ORDERS_for_ST[k],
repetitions=REPETITIONS_for_ST[k],
qbv=qbv,
),
)
qmod = create_model(main)
qmod = set_preferences(qmod, preferences=preferences)
qprog = synthesize(qmod)
circuit = QuantumProgram.from_qprog(qprog)
classiq_depths.append(circuit.transpiled_circuit.depth)
classiq_cx_counts.append(circuit.transpiled_circuit.count_ops["cx"])
3. Approximating with the qDRIFT Formula
from classiq import qdrift
for n_qd in N_QDS_for_qDRIFT:
@qfunc
def main() -> None:
qbv = QArray("qbv")
allocate(len(pauli_list[0][0]), qbv)
qdrift(
pauli_operator=hamiltonian,
evolution_coefficient=1,
num_qdrift=n_qd,
qbv=qbv,
)
qmod = create_model(main)
qmod = set_preferences(qmod, preferences=preferences)
qprog = synthesize(qmod)
circuit = QuantumProgram.from_qprog(qprog)
classiq_depths.append(circuit.transpiled_circuit.depth)
classiq_cx_counts.append(circuit.transpiled_circuit.count_ops["cx"])
4. Comparing to Qiskit Implementation
Comments: * Qiskit's Suzuki-Trotter of order 1 is a separate function, called the Lie-Trotter function. * Qiskit's qDRIFT takes a long time to run for unclear problems. Alternatively, a random product formula is implemented, which is equivalent to qDRIFT.
The qiskit data was generated using qiskit version 1.0.0. To run the qiskit code uncomment the commented cells below.
qiskit_cx_counts = [25164, 33508, 41882, 186067, 248232, 3581, 7404]
qiskit_depths = [30627, 42879, 53581, 300924, 402125, 3463, 7141]
# from importlib.metadata import version
# try:
# import qiskit
# if version('qiskit') != "1.0.0":
# !pip uninstall qiskit -y
# !pip install qiskit==1.0.0
# except ImportError:
# !pip install qiskit==1.0.0
# from qiskit.circuit.library import PauliEvolutionGate
# from qiskit.quantum_info import SparsePauliOp
# from qiskit.synthesis import LieTrotter as LieTrotter_qiskit
# qiskit_depths = []
# qiskit_cx_counts = []
# operator = SparsePauliOp.from_list(pauli_list)
# gate = PauliEvolutionGate(operator, 1)
# ## Suzuki-Trotter of order 1 = Lie-Trotter
# from qiskit import transpile
# lt = LieTrotter_qiskit(reps=REPETITIONS_for_ST[0])
# circ = lt.synthesize(gate)
# tqc = transpile(
# circ, basis_gates=["u", "cx"], optimization_level=transpilation_options["qiskit"]
# )
# qiskit_depths.append(tqc.depth())
# qiskit_cx_counts.append(tqc.count_ops()["cx"])
# ## Suzuki-Trotter
# from qiskit.synthesis import SuzukiTrotter as SuzukiTrotter_qiskit
# for k in range(1, 3):
# st = SuzukiTrotter_qiskit(order=ORDERS_for_ST[k], reps=REPETITIONS_for_ST[k])
# circ = st.synthesize(gate)
# tqc = transpile(
# circ,
# basis_gates=["u", "cx"],
# optimization_level=transpilation_options["qiskit"],
# )
# qiskit_depths.append(tqc.depth())
# qiskit_cx_counts.append(tqc.count_ops()["cx"])
# # Controlled dynamics
# lt_ctrl = LieTrotter_qiskit(reps=REPETITIONS_for_ST[0])
# circ = lt_ctrl.synthesize(gate).control(1)
# tqc = transpile(
# circ, basis_gates=["u", "cx"], optimization_level=transpilation_options["qiskit"]
# )
# qiskit_depths.append(tqc.depth())
# qiskit_cx_counts.append(tqc.count_ops()["cx"])
# for k in range(1, 2):
# st_ctrl = SuzukiTrotter_qiskit(order=ORDERS_for_ST[k], reps=REPETITIONS_for_ST[k])
# circ = st_ctrl.synthesize(gate).control(1)
# tqc = transpile(
# circ,
# basis_gates=["u", "cx"],
# optimization_level=transpilation_options["qiskit"],
# )
# qiskit_depths.append(tqc.depth())
# qiskit_cx_counts.append(tqc.count_ops()["cx"])
# ## For qDRIFT, generate a random sequence
# import numpy as np
# def index_channel(n, list_coe):
# """
# This function gets an ordered list of coefficients 'list_coe' and a number of calls 'n' as inputs.
# It returns a random ordered list of size n with elements from list_coe',
# where the probability of choosing the i-th elements is list_coe[i]/sum(list_coe),
# """
# coe = np.array(list_coe) / sum(list_coe)
# c_coe = np.cumsum(coe)
# return np.searchsorted(c_coe, np.random.uniform(size=n))
# assert (
# pauli_list[0][0] == len(pauli_list[0][0]) * "I"
# ), """The Identity term is not the first on the list of Paulis,
# please modify the code accordingly """
# pauli_list_without_id = pauli_list[1::]
# po_coe = [
# np.abs(np.real(p[1])) for p in pauli_list_without_id
# ] # gets absolute value of coefficients
# for n_qd in N_QDS_for_qDRIFT:
# new_indices = index_channel(n_qd, po_coe)
# small_lambda = sum(po_coe)
# randomly_generated_pauli_list = [
# (
# pauli_list_without_id[new_indices[k]][0],
# np.sign(pauli_list_without_id[new_indices[k]][1]) * small_lambda / n_qd,
# )
# for k in range(n_qd)
# ]
# randomly_generated_pauli_list += [pauli_list[0]] # adding the identity
# gate_qd = PauliEvolutionGate(
# SparsePauliOp.from_list(randomly_generated_pauli_list), 1
# )
# qd = LieTrotter_qiskit(reps=1)
# circ = qd.synthesize(gate_qd)
# tqc = transpile(
# circ,
# basis_gates=["u", "cx"],
# optimization_level=transpilation_options["qiskit"],
# )
# qiskit_depths.append(tqc.depth())
# qiskit_cx_counts.append(tqc.count_ops()["cx"])
5. Plotting the Data
print("The cx-counts on Classiq:", classiq_cx_counts)
print("The cx-counts on Qiskit:", qiskit_cx_counts)
The cx-counts on Classiq: [15588, 20798, 25992, 22188, 29598, 3330, 6330]
The cx-counts on Qiskit: [25164, 33508, 41882, 186067, 248232, 3581, 7404]
import matplotlib.pyplot as plt
classiq_color = "#D7F75B"
qiskit_color = "#6FA4FF"
plt.rcParams["font.family"] = "serif"
plt.rc("savefig", dpi=300)
plt.rcParams["axes.linewidth"] = 1
plt.rcParams["xtick.major.size"] = 5
plt.rcParams["xtick.minor.size"] = 5
plt.rcParams["ytick.major.size"] = 5
plt.rcParams["ytick.minor.size"] = 5
plt.semilogy(
qiskit_cx_counts[-2::] + qiskit_cx_counts[0:5],
"s",
label="qiskit",
markerfacecolor=qiskit_color,
markeredgecolor="k",
markersize=6,
markeredgewidth=1.5,
)
plt.semilogy(
classiq_cx_counts[-2::] + classiq_cx_counts[0:5],
"o",
label="classiq",
markerfacecolor=classiq_color,
markeredgecolor="k",
markersize=6.5,
markeredgewidth=1.5,
)
labels = [
"qDRIFT(N=1000)",
"qDRIFT(N=2000)",
"TS$_1$(reps=6)",
"TS$_2$(reps=4)",
"TS$_4$(reps=1)",
"ctrl-TS$_1$(reps=6)",
"ctrl-TS$_2$(reps=4)",
]
plt.xticks([0, 1, 2, 3, 4, 5, 6], labels, rotation=45, fontsize=16, ha="right")
plt.ylabel("CX-counts", fontsize=18)
plt.yticks(fontsize=16)
plt.legend(loc="upper left", fontsize=18, fancybox=True, framealpha=0.5)
<matplotlib.legend.Legend at 0x7f165654b210>
