The aim of this tutorial is to introduce some of the basic and more advanced features of the ProjectQ simulator. Please note that all the simulator features can be found in our code documentation.
Contents:
Our simulator can be used to simulate any circuit model quantum algorithm. This requires storing the state, also called wavefunction, of all qubits which requires storing 2n complex values (each of size 16 bytes) for an n-qubit algorithm. This can get very expensive in terms of required RAM:
| Number of qubits n | Required RAM to store wavefunction |
|---|---|
| 10 | 16 KByte |
| 20 | 16 MByte |
| 30 | 16 GByte |
| 31 | 32 GByte |
| 32 | 64 GByte |
| 40 | 16 TByte |
| 45 | 512 TByte (world's largest quantum computer simulation) |
The number of qubits you can simulate with the ProjectQ simulator is only limited by the amount of RAM in your notebook or workstation. Applying quantum gates is expensive as we have to potentially update each individual value in the full wavefunction. Therefore, we have implemented a high-performance simulator which is significantly faster than other simulators (see our papers for a detailed comparison [1], [2]). The performance of such simulators is hardware-dependent and therefore we have decided to provide 4 different versions. A simulator implemented in C++ which uses multi-threading (OpenMP) and instrinsics, a C++ simulator which uses intrinsics, a C++ simulator, and a slower simulator which only uses Python. During the installation we try to install the fastest of these options given your hardware and available compiler.
Our simulator is simultaneously also a quantum emulator. This concept was first introduced by us in [2]. A quantum emulator takes classical shortcuts when performing the simulation and hence can be orders of magnitude faster. E.g., for simulating Shor's algorithm, we only require to store the wavefunction of n+1 qubits while the algorithm on a real quantum computer would require 2n+O(1) qubits. Using these emulation capabilities, we can easily emulate Shor's algorithm for factoring, e.g., 4028033 = 2003 · 2011 on a notebook [1].
Please follow the installation instructions in the docs. The Python interface to all our different simulators (C++ or native Python) is identical. The different versions only differ in performance. If you have a C++ compiler installed on your system, the setup will try to install the faster C++ simulator. To figure out which simulator is installed just execute the following code after installing ProjectQ:
In [1]:
import projectq
eng = projectq.MainEngine() # This loads the simulator as it is the default backend
If you now see the following message
(Note: This is the (slow) Python simulator.)you are using the slow Python simulator and we would recommend to reinstall ProjectQ with the C++ simulator following the installation instructions. If this message doesn't show up, then you are using one of the fast C++ simulator versions. Which one exactly depends on your compiler and hardware.
To write a quantum program, we need to create a compiler called MainEngine and provide a backend for which the compiler should compile the quantum program. In this tutorial we are focused on the simulator as a backend. We can create a compiler with a simulator backend by importing the simulator class and creating a simulator instance which is passed to the MainEngine:
In [2]:
from projectq.backends import Simulator
eng = projectq.MainEngine(backend=Simulator())
Note that the MainEngine contains the simulator as the default backend, hence one can equivalently just do:
In [3]:
eng = projectq.MainEngine()
In order to simulate the probabilistic measurement process, the simulator internally requires a random number generator. When creating a simulator instance, one can provide a random seed in order to create reproducible results:
In [4]:
eng = projectq.MainEngine(backend=Simulator(rnd_seed=10))
Let's write a simple test program which allocates a qubit in state |0>, performs a Hadamard gate which puts it into a superposition 1/sqrt(2) * ( |0> + |1>) and then measure the qubit:
In [5]:
import projectq
from projectq.ops import Measure, H
eng = projectq.MainEngine()
qubit = eng.allocate_qubit() # Allocate a qubit from the compiler (MainEngine object)
H | qubit
Measure | qubit # Measures the qubit in the Z basis and collapses it into either |0> or |1>
eng.flush() # Tells the compiler (MainEninge) compile all above gates and send it to the simulator backend
print(int(qubit)) # The measurement result can be accessed by converting the qubit object to a bool or int
This program randomly outputs 0 or 1. Note that the measurement does not return for example a probability of measuring 0 or 1 (see below how this could be achieved). The reason for this is that a program written in our high-level syntax should be independent of the backend. In other words, the code can be executed either by our simulator or by exchanging only the MainEngine's backend by a real device which cannot return a probability but only one value. See the other examples on GitHub of how to execute code on the IBM quantum chip.
Note that we have used eng.flush() which tells the compiler to send all the instructions to the backend. In a simplified version, when the Python interpreter executes a gate (e.g. the above lines with H, or Measure), this gate is sent to the compiler (MainEngine), where it is cached. Compilation and optimization of cached gates happens irregularly, e.g., an optimizer in the compiler might wait for 5 gates before it starts the optimization process. Therefore, if we require the result of a quantum program, we need to call eng.flush() which tells the compiler to compile and send all cached gates to the backend for execution. eng.flush() is therefore necessary when accessing the measurement result by converting the qubit object to an int or bool. Or when using the advanced features below, where we want to access properties of the wavefunction at a specific point in the quantum algorithm.
The above statement is not entirely correct as for example there is no eng.flush() strictly required before accessing the measurement result. The reason is that the measurement gate in our compiler is not cached but is sent directly to the local simulator backend because it would allow for performance optimizations by shrinking the wavefunction. However, this is not a feature which your code should/can rely on and therefore you should always use eng.flush().
When a qubit goes out of scope, it gets deallocated automatically. If the backend is a simulator, it checks that the qubit was in a classical state and otherwise it raises an exception. This is an important debugging feature as in many quantum algorithms, ancilla qubits are used for intermediate results and then "uncomputed" back into state |0>. If such ancilla qubits now go out of scope, the simulator throws an error if they are not in either state |0> or |1>, as this is most likely a bug in the user code:
In [6]:
def test_program(eng):
# Do something
ancilla = eng.allocate_qubit()
H | ancilla
# Use the ancilla for something
# Here the ancilla is not reset to a classical state but still in a superposition and will go out of scope
test_program(eng)
If you are using a qubit as an ancilla which should have been reset, this is a great feature which automatically checks the correctness of the uncomputation if the simulator is used as a backend. Should you wish to deallocate qubits which might be in a superposition, always apply a measurement gate in order to avoid this error message:
In [7]:
from projectq.ops import All
def test_program_2(eng):
# Do something
ancillas = eng.allocate_qureg(3) # allocates a quantum register with 3 qubits
All(H) | ancillas # applies a Hadamard gate to each of the 3 ancilla qubits
All(Measure) | ancillas # Measures all ancilla qubits such that we don't get
#an error message when they are deallocated
# Here the ancillas will go out of scope but because of the measurement, they are in a classical state
test_program_2(eng)
Here we will show some features which are unique to a simulator backend which has access to the full wavefunction. Note that using these features in your code will prohibit the code to run on a real quantum device. Therefore instead of, e.g., using the feature to ininitialize the wavefunction into a specific state, you could execute a small quantum circuit which prepares the desired state and hence the code can be run on both the simulator and on actual hardware.
For details on the simulator please see the code documentation.
In order to use the features of the simulator backend, we need to have a reference to access it. This can be achieved by creating a simulator instance and keeping a reference to it before passing it to the MainEngine:
In [8]:
sim = Simulator(rnd_seed=5) # We can now always access the simulator via the "sim" variable
eng = projectq.MainEngine(backend=sim)
Alternatively, one can access the simulator by accessing the backend of the compiler (MainEngine):
In [9]:
assert id(sim) == id(eng.backend)
In [10]:
from projectq.ops import H, Measure
eng = projectq.MainEngine()
qubit = eng.allocate_qubit()
qubit2 = eng.allocate_qubit()
eng.flush() # sends the allocation of the two qubits to the simulator (only needed to showcase the stuff below)
H | qubit # Put this qubit into a superposition
# qubit is a list with one qubit, qubit2 is another list with one qubit object
# qubit + qubit2 creates a list containing both qubits. Equivalently, one could write [qubit[0], qubit2[0]]
# get_amplitude requires that one provides a list/qureg of all qubits such that it can determine the order
amp_before = eng.backend.get_amplitude('00', qubit + qubit2)
# Amplitude will be 1 as Hadamard gate is not yet executed on the simulator backend
# We forgot the eng.flush()!
print("Amplitude saved in amp_before: {}".format(amp_before))
eng.flush() # Makes sure that all the gates are sent to the backend and executed
amp_after = eng.backend.get_amplitude('00', qubit + qubit2)
# Amplitude will be 1/sqrt(2) as Hadamard gate was executed on the simulator backend
print("Amplitude saved in amp_after: {}".format(amp_after))
# To avoid triggering the warning of deallocating qubits which are in a superposition
Measure | qubit
Measure | qubit2
NOTE: One always has to call eng.flush() before accessing the amplitude as otherwise some of the gates might not have been sent and executed on the simulator. Also don't forget in such an example to measure all the qubits in the end in order to avoid the above mentioned debugging error message of deallocating qubits which are not in a classical state.
In [11]:
import projectq
from projectq.ops import H, Measure, CNOT, All
eng = projectq.MainEngine()
qureg = eng.allocate_qureg(2)
H | qureg[0]
eng.flush()
prob11 = eng.backend.get_probability('11', qureg)
prob10 = eng.backend.get_probability('10', qureg)
prob01 = eng.backend.get_probability('01', qureg)
prob00 = eng.backend.get_probability('00', qureg)
prob_second_0 = eng.backend.get_probability('0', [qureg[1]])
print("Probability to measure 11: {}".format(prob11))
print("Probability to measure 00: {}".format(prob00))
print("Probability to measure 01: {}".format(prob01))
print("Probability to measure 10: {}".format(prob10))
print("Probability that second qubit is in state 0: {}".format(prob_second_0))
All(Measure) | qureg
In [12]:
import projectq
from projectq.ops import H, Measure, CNOT, All, QubitOperator, X, Y
eng = projectq.MainEngine()
qureg = eng.allocate_qureg(3)
X | qureg[0]
H | qureg[1]
eng.flush()
op0 = QubitOperator('Z0') # Z applied to qureg[0] tensor identity on qureg[1], qureg[2]
expectation = eng.backend.get_expectation_value(op0, qureg)
print("Expectation value = <Psi|Z0|Psi> = {}".format(expectation))
op_sum = QubitOperator('Z0 X1') - 0.5 * QubitOperator('X1')
expectation2 = eng.backend.get_expectation_value(op_sum, qureg)
print("Expectation value = <Psi|-0.5 X1 + 1.0 Z0 X1|Psi> = {}".format(expectation2))
All(Measure) | qureg # To avoid error message of deallocating qubits in a superposition
For debugging purposes, one might want to check the algorithm for cases where an intermediate measurement outcome was, e.g., 1. Instead of running many simulations and post selecting only those with the desired intermediate measurement outcome, our simulator allows to force a specific measurement outcome. Note that this is only possible if the desired state has non-zero amplitude, otherwise the simulator will throw an error.
In [13]:
import projectq
from projectq.ops import H, Measure
eng = projectq.MainEngine()
qureg = eng.allocate_qureg(2)
# Create an entangled state:
H | qureg[0]
CNOT | (qureg[0], qureg[1])
# qureg is now in state 1/sqrt(2) * (|00> + |11>)
Measure | qureg[0]
Measure | qureg[1]
eng.flush() # required such that all above gates are executed before accessing the measurement result
print("First qubit measured in state: {} and second qubit in state: {}".format(int(qureg[0]), int(qureg[1])))
Running the above circuit will either produce both qubits in state 0 or both qubits in state 1. Suppose I want to check the outcome if the first qubit was measured in state 0. This can be achieve by telling the simulator backend to collapse the wavefunction for the first qubit to be in state 0:
In [14]:
import projectq
from projectq.ops import H, CNOT
eng = projectq.MainEngine()
qureg = eng.allocate_qureg(2)
# Create an entangled state:
H | qureg[0]
CNOT | (qureg[0], qureg[1])
# qureg is now in state 1/sqrt(2) * (|00> + |11>)
eng.flush() # required such that all above gates are executed before collapsing the wavefunction
# We want to check what happens to the second qubit if the first qubit (qureg[0]) is measured to be 0
eng.backend.collapse_wavefunction([qureg[0]], [0])
# Check the probability that the second qubit is measured to be 0:
prob_0 = eng.backend.get_probability('0', [qureg[1]])
print("After forcing a measurement outcome of the first qubit to be 0, \n"
"the second qubit is in state 0 with probability: {}".format(prob_0))
It is possible to set the state of the simulator to any arbitrary wavefunction. In a first step one needs to allocate all the required qubits (don't forget to call eng.flush()), and then one can use this method to set the wavefunction. Note that this only works if the provided wavefunction is the wavefunction of all allocated qubits. In addition, the wavefunction needs to be normalized. Here is an example:
In [15]:
import math
import projectq
from projectq.ops import H
eng = projectq.MainEngine()
qureg = eng.allocate_qureg(2)
eng.flush()
eng.backend.set_wavefunction([1/math.sqrt(2), 1/math.sqrt(2), 0, 0], qureg)
H | qureg[0]
# At this point both qubits are back in the state 00 and hence there will be no exception thrown
# when the qureg is deallocated
Cheat is the original method to access and manipulate the full wavefunction. Calling cheat with the C++ simulator returns a copy of the full wavefunction plus the mapping of which qubit is at which bit position. The Python simulator returns a reference. If possible we are planning to change the C++ simulator to also return a reference which currently is not possible due to the python export. Please keep this difference in mind when writing code. If you require a copy, it is safest to make a copy of the objects returned by the cheat method.
When qubits are allocated in the code, each of the qubits gets a unique integer id. This id is important in order to understand the wavefunction returned by cheat. The wavefunction is a numpy array of length 2n, where n is the number of qubits. Which bitlocation a specific qubit in the wavefunction has is not predefined (e.g. by the order of qubit allocation) but is rather chosen depending on the compiler optimizations and the simulator. Therefore, cheat also returns a dictionary containing the mapping of qubit id to bit location in the wavefunction. Here is a small example:
In [16]:
import copy
import projectq
from projectq.ops import Rx, Measure, All
eng = projectq.MainEngine()
qubit1 = eng.allocate_qubit()
qubit2 = eng.allocate_qubit()
Rx(0.2) | qubit1
Rx(0.4) | qubit2
eng.flush() # In order to have all the above gates sent to the simulator and executed
# We save a copy of the wavefunction at this point in the algorithm. In order to make sure we get a copy
# also if the Python simulator is used, one should make a deepcopy:
mapping, wavefunction = copy.deepcopy(eng.backend.cheat())
print("The full wavefunction is: {}".format(wavefunction))
# Note: qubit1 is a qureg of length 1, i.e. a list containing one qubit objects, therefore the
# unique qubit id can be accessed via qubit1[0].id
print("qubit1 has bit-location {}".format(mapping[qubit1[0].id]))
print("qubit2 has bit-location {}".format(mapping[qubit2[0].id]))
# Suppose we want to know the amplitude of the qubit1 in state 0 and qubit2 in state 1:
state = 0 + (1 << mapping[qubit2[0].id])
print("Amplitude of state qubit1 in state 0 and qubit2 in state 1: {}".format(wavefunction[state]))
# If one only wants to access one (or a few) amplitudes, get_amplitude provides an easier interface:
amplitude = eng.backend.get_amplitude('01', qubit1 + qubit2)
print("Accessing same amplitude but using get_amplitude instead: {}".format(amplitude))
All(Measure) | qubit1 + qubit2 # In order to not deallocate a qubit in a superposition state
In [17]:
import projectq
from projectq.backends import Simulator
eng = projectq.MainEngine(backend=Simulator(), engine_list=[]) # Removes the default compiler engines
As noted in the code documentation, one can play around with the number of threads in order to increase the simulation speed. Execute the following statements in the terminal before running ProjectQ:
export OMP_NUM_THREADS=2
export OMP_PROC_BIND=spreadA good setting is to set the number of threads to the number of physical cores on your system.
The simulator has a feature called "gate fusion" in which case it combines smaller gates into larger ones in order to increase the speed of the simulation. If the simulator is faster with or without gate fusion depends on many parameters. By default it is currently turned off but one can turn it on and compare by:
In [18]:
import projectq
from projectq.backends import Simulator
eng = projectq.MainEngine(backend=Simulator(gate_fusion=True)) # Enables gate fusion
We'd like to refer interested readers to our paper on the world's largest and fastest quantum computer simulation for more details on how to optimize the speed of a quantum simulator.
In [ ]: