Simulating the dynamics of the molecule#
Quantum computing offers a promising avenue for simulating complex molecular systems, potentially revolutionizing fields like drug discovery, materials science, and chemical engineering. At the heart of this capability lies Hamiltonian simulation, a powerful technique that allows us to model the quantum behavior of molecules with exponentially less resources compared to classical computers. In this tutorial, we’ll explore how to simulate molecules using Hamiltonian simulation techniques implemented in Qrisp. Since version 0.5, Qrisp provides a user-friendly interface for simulating quantum many body systems, making it a suitable tool for tackling complex quantum chemistry problems. We’ll cover the following key aspects:
By the end of this tutorial, you’ll have a solid understanding of how to leverage Qrisp’s capabilities to perform molecular simulations on quantum computers. Let’s dive in and discover how Qrisp can help us unlock the quantum secrets of molecules!
The fundamentals of Hamiltonian simulation and the second quantization#
To understand what “second quantization” means it might be helpful to start with the first quantization. The idea behind this is that physical systems containing only a single particale are described by a Hilbert space
An arbitrary state
What these particular

This plot represents the real part of the amplitude of one of the basis vectors (the one with
The second quantization now describes how systems including multiple particles are treated. For this the so called ladder operators are introduced, that is, for each basis state
For this reason, this operator is called the “creator”. Similarly there is also the “annihilator”
The ladder operators therefore determine the intrinsic interaction behavior of the particles. Observations from experiments show that there are two types of particles whose ladder operators follow different sets of laws.
For Bosons the ladder operators obey commutator laws:
For Fermions the ladder operators obey anti-commutator laws:
Note that the fermionic laws imply
Within Qrisp it is currently only possible to model fermions, which is for many applications in chemistry the more important case. A modelling framework for bosons will follow in a future release. To start building a fermionic operator, we import the functions c
and a
for creators and annihilators.
from qrisp.operators import c, a
O = a(0)*c(1) + a(1)*a(2)
print(O)
# Yields: a0*c1 + a1*a2
To learn more how to build and manipulate these expressions, please look at the documentation page of FermionicOperator. For instance, the hermitian conjugate can be computed using the .dagger
method.
print(O.dagger())
# Yields: a1*c0 + c2*c1
To apply the Pauli exclusion principle but also other anti-commutation laws for simplifaction, you can call the reduce
method.
O = a(0)*a(0) + a(1)*a(2) - a(2) * a(1)
print(O.reduce())
#Yields 2*a1*a2
The Jordan-Wigner embedding#
A natural question that comes up is how to represent the ladder operators and the corresponding states on a quantum computer. The most established way to do this is to use the Jordan-Wigner embedding (even though there are several interesting alternatives). The Jordan-Wigner embedding identifies each ladder term with an operator that acts on a qubit space:
Where
O_fermionic = a(4)
O_qubit = O_fermionic.to_qubit_operator(mapping_type = "jordan_wigner")
print(O_qubit)
# Yields: Z_0*Z_1*Z_2*Z_3*A_4
This gives us an instance of the QubitOperator class. What is the difference to a FermionicOperator? While FermionicOperators model the more abstract fermion space, qubit operators represent operators on the qubit space
Dynamics#
Both boson and fermion systems evolve under the Schrödinger equation:
Where
For bosonic systems, the Hamiltonian can only be a linear combination of products of the bosonic ladder operators. The equivalent holds for fermionic systems.
Where all
The particular values of the coefficients (like
Loading molecular data into Qrisp#
If you don’t feel like solving integrals right now, we’ve got you covered! Qrisp has a convenient interface to PySCF, which loads the molecular data directly as FermionicOperator. For that you need PySCF installed (pip install pyscf
). If you’re on Windows you might need to do some WSL gymnastics.
from pyscf import gto
from qrisp.operators import FermionicOperator
mol = gto.M(atom = '''H 0 0 0; H 0 0 0.74''', basis = 'sto-3g')
H_ferm = FermionicOperator.from_pyscf(mol)
print(H_ferm)
This snippet uses the .from_pyscf
method to load the FermionicOperator representing the orbitals of the Dihydrogen molecule
-0.181210462015197*a0*a1*c2*c3 + 0.181210462015197*a0*c1*c2*a3
- 1.25330978664598*c0*a0 + 0.674755926814448*c0*a0*c1*a1
+ 0.482500939335616*c0*a0*c2*a2 + 0.663711401350814*c0*a0*c3*a3
+ 0.181210462015197*c0*a1*a2*c3 - 0.181210462015197*c0*c1*a2*a3
- 1.25330978664598*c1*a1 + 0.663711401350814*c1*a1*c2*a2
+ 0.482500939335616*c1*a1*c3*a3 - 0.475068848772178*c2*a2
+ 0.697651504490463*c2*a2*c3*a3 - 0.475068848772178*c3*a3
Or if preferred, the Jordan-Wigner embedding:
H_qubit = H_ferm.to_qubit_operator()
print(H_qubit)
0.181210462015197*A_0*A_1*C_2*C_3 - 0.181210462015197*A_0*C_1*C_2*A_3
- 0.181210462015197*C_0*A_1*A_2*C_3 + 0.181210462015197*C_0*C_1*A_2*A_3
- 1.25330978664598*P^0_0 + 0.674755926814448*P^0_0*P^0_1
+ 0.482500939335616*P^0_0*P^0_2 + 0.663711401350814*P^0_0*P^0_3
- 1.25330978664598*P^0_1 + 0.663711401350814*P^0_1*P^0_2
+ 0.482500939335616*P^0_1*P^0_3 - 0.475068848772178*P^0_2
+ 0.697651504490463*P^0_2*P^0_3 - 0.475068848772178*P^0_3
Performing Hamiltonian simulation#
To perform Hamiltonian simulation, we use the .trotterization
method, which gives us a Python function that performs a simulation of the Hamiltonian on a QuantumVariable.
from qrisp import QuantumVariable
electron_state = QuantumVariable(4)
electron_state[:] = {"1100": 2**-0.5, "0001": 2**-0.5}
This snippet initializes the state
U = H_ferm.trotterization()
U(electron_state, t = 100, steps = 20)
This snippets simulates the Dihydrogen molecule for
Finally, we want to extract some physical quantity from our simulation. Our quantity of choice is the particle number operator:
In python code:
N = sum(c(i)*a(i) for i in range(4))
For each state .get_measurement
method.
expectation_value = N.get_measurement(electron_state, precision = 0.01)
print(expectation_value)
# Yields: 1.50440973329083
We see that the expectaction value is (almost) 1.5 because
print(N.get_measurement(electron_state, precision = 0.0001))
# Yields: 1.4999799028124874
This concludes our little tutorial on Hamiltonian simulation. We hope you could learn something and feel motivated to explore more systems and techniques! Make sure to also check out the Molecular Potential Energy Curves example to learn how to compute the ground state for both types of operators!
For your reference, we give the full code below:
# Loading molecular data
from qrisp.operators import a,c, FermionicOperator
from pyscf import gto
mol = gto.M(atom = '''H 0 0 0; H 0 0 0.74''', basis = 'sto-3g')
H_ferm = FermionicOperator.from_pyscf(mol)
# Initializing the quantum state
from qrisp import QuantumVariable
electron_state = QuantumVariable(4)
electron_state[:] = {"1100": 2**-0.5, "0001": 2**-0.5}
# Performing the simulation
U = H_ferm.trotterization()
U(electron_state, t = 100, steps = 20)
# Evaluating the number operator
N = sum(c(i)*a(i) for i in range(4))
expectation_value = N.get_measurement(electron_state, precision = 0.01)
print(expectation_value)