Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Single- and two-qubit state and process tomography #246

Merged
merged 22 commits into from
Nov 10, 2024
Merged
Changes from all commits
Commits
Show all changes
22 commits
Select commit Hold shift + click to select a range
0890d33
Add Use Case 4 folder on both state tomography and standard process t…
bornman-nick Sep 4, 2024
f202b9b
Add Use Case 4 folder on both state tomography and standard process t…
bornman-nick Sep 4, 2024
8d0203c
Add license file
bornman-nick Sep 4, 2024
f411b58
Modify License file to BSD-3
bornman-nick Sep 6, 2024
ccd3faf
Update single-qubit-state-tomography.py
bornman-nick Oct 10, 2024
7de5184
exclude helper_function.py from black
yomach Oct 10, 2024
3d7d809
Update helper_functions.py
bornman-nick Oct 10, 2024
5903617
Merge pull request #2 from bornman-nick/bornman-nick-patch-1
bornman-nick Oct 10, 2024
303b0fc
Update pyproject.toml
bornman-nick Oct 10, 2024
60d9595
Update helper_functions.py
bornman-nick Oct 11, 2024
a78a346
Update single-qubit-process-tomography.py
bornman-nick Oct 11, 2024
6fdde33
Update single-qubit-state-tomography.py
bornman-nick Oct 11, 2024
b991f1e
Update two-qubit-process-tomography.py
bornman-nick Oct 11, 2024
8f449ee
Update two-qubit-state-tomography.py
bornman-nick Oct 11, 2024
e1f6bb6
Formatting
MichalGQM Oct 14, 2024
07b85d3
Add README and accompanying images, configuration.py file, and change…
bornman-nick Oct 17, 2024
27f45af
Update README.md
bornman-nick Oct 17, 2024
5060b3b
Merge branch 'main' into main
bornman-nick Oct 17, 2024
7ae6905
Update README.md
bornman-nick Oct 18, 2024
46f752c
Update single-qubit-process-tomography.py
bornman-nick Oct 22, 2024
ec5c5e1
Update two-qubit-state-tomography.py
bornman-nick Oct 22, 2024
80b00f4
Reformat configuration.py and two-qubit-state-tomography.py using black
bornman-nick Oct 31, 2024
File filter

Filter by extension

Filter by extension


Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
Original file line number Diff line number Diff line change
@@ -0,0 +1,11 @@
Copyright 2024. Fermi Research Alliance, LLC.
All rights reserved.

Copyright 2024. Fermi Research Alliance, LLC. This software was produced under U.S. Department of Energy, Office of Science, National Quantum Information Science Research Centers, Superconducting Quantum Materials and Systems (SQMS) Center, under the contract No. DE-AC02-07CH11359, for Fermi National Accelerator Laboratory (Fermilab), which is operated by Fermi Research Alliance, LLC for the U.S. Department of Energy. The U.S. Government has rights to use, reproduce, and distribute this software. NEITHER THE GOVERNMENT NOR FERMI RESEARCH ALLIANCE, LLC MAKES ANY WARRANTY, EXPRESS OR IMPLIED, OR ASSUMES ANY LIABILITY FOR THE USE OF THIS SOFTWARE. If software is modified to produce derivative works, such modified software should be clearly marked, so as not to confuse it with the version available from Fermilab.

Additionally, redistribution and use in source and binary forms, with or without modification, are permitted provided that the following conditions are met:
• Redistributions of source code must retain the above copyright notice, this list of conditions and the following disclaimer.
• Redistributions in binary form must reproduce the above copyright notice, this list of conditions and the following disclaimer in the documentation and/or other materials provided with the distribution.
• Neither the name of Fermi Research Alliance, LLC, Fermi National Accelerator Laboratory, Fermilab, the U.S. Government, nor the names of its contributors may be used to endorse or promote products derived from this software without specific prior written permission.

THIS SOFTWARE IS PROVIDED BY FERMI RESEARCH ALLIANCE, LLC AND CONTRIBUTORS "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL FERMI RESEARCH ALLIANCE, LLC OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
Binary file not shown.

Large diffs are not rendered by default.

Loading
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
Loading
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.

Large diffs are not rendered by default.

Large diffs are not rendered by default.

Loading
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
Loading
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
Loading
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
Original file line number Diff line number Diff line change
@@ -0,0 +1,242 @@
#!/usr/bin/env python

"""
SINGLE QUBIT PROCESS TOMOGRAPHY
The sequence consists of preparing the qubit into one of the six cardinal Bloch sphere
states using calibrated gates from the config, applying the operation/process - typically
ideally a unitary matrix - under investigation, and then measuring the state of the qubit,
by way of the readout resonator, in the X, Y and Z bases (specifically, measuring the
projector of one of the same six cardinal Bloch sphere states). This is repeated for
the full set of input states and measurement projectors, which is tomographically
complete. The output, 'probs', containing the 'Bloch sphere' prepared states and
measurement projectors, is mapped to the Pauli basis (giving a 'measurement vector'),
and the equation "measurement_vector = B \times chi_vector" (where B is a constant)
inverted to give the chi process matrix for the process under
investigation.
Prerequisites:
- Having found the resonance frequency of the resonator coupled to the qubit under study (resonator_spectroscopy)
- Having calibrated qubit pi and pi/2 pulses by running qubit spectroscopy, rabi_chevron, power_rabi and updating the config
- Having calibrated the readout (readout_frequency, amplitude, duration_optimization IQ_blobs) for better SNR, and having
saved the derived readout threshold value in the config
- Set the desired flux bias in the case of flux-tunable qubits
This script implements the logic outlined in the notebook
https://github.com/bornman-nick/quantum-state-and-process-tomography. See
that notebook for a detailed explanation of standard single qubit process
tomography.
Note that this routine carries out standard process tomography on a single qubit.
Given that the measurements recorded are simply whether the final state of the
qubit is in the ground state or not for each state preparation and measurement
setting (instead of also deciding whether the qubit is perhaps in the excited
state and incrementing the counter of a complementary measurement setting),
this routine scales as 6**2, rather than 6*3, as is the case in
most formulations of standard process tomography. However, the routine below
is fast enough that this doesn't matter too much.
"""

from qm.qua import *
from qm import QuantumMachinesManager
from qm import SimulationConfig
from configuration import *
from qualang_tools.results import progress_counter, fetching_tool
import numpy as np
from scipy.linalg import solve

from helper_functions import (
P_Pauli1,
plot_process_tomography1,
map_from_bloch_state_to_pauli_basis1,
func_E1,
)


###################
# The QUA program #
###################

# qubit under test, assuming there are multiple qubits
# with XY line elements "q<qubit>_xy", flux lines
# "q<qubit>_z", and readout resonator elements "rr<qubit>"
qubit = 1

if qubit == 1:
threshold = ge_threshold_q1
elif qubit == 2:
threshold = ge_threshold_q2
else:
raise ValueError(f"Incorrect qubit number chosen")


n_avg = 10_000


# subroutine to prepare desired qubit gate/process
def analysed_process(qubit):
# write whatever QUA code you need, here, in order to perform the
# desired process which we want to subject to tomography.
# For example, to analyse the Y gate:
play("y180", f"q{qubit}_xy")


ideal_gate = func_E1(2)


# subroutine to switch between preparing the tomographically-complete
# set of input states
def prepare_state(i):
with switch_(i):
with case_(0):
wait(pi_len // 4, f"q{qubit}_xy")
with case_(1):
play("x180", f"q{qubit}_xy")
with case_(2):
play("y90", f"q{qubit}_xy")
with case_(3):
play("-y90", f"q{qubit}_xy")
with case_(4):
play("-x90", f"q{qubit}_xy")
with case_(5):
play("x90", f"q{qubit}_xy")


# subroutine to switch between the tomographically-complete
# set of bases in which to measure
def measurement_basis_change(j):
with switch_(j):
with case_(0):
wait(pi_len // 4, f"q{qubit}_xy")
with case_(1):
play("x180", f"q{qubit}_xy")
with case_(2):
play("y90", f"q{qubit}_xy")
with case_(3):
play("-y90", f"q{qubit}_xy")
with case_(4):
play("-x90", f"q{qubit}_xy")
with case_(5):
play("x90", f"q{qubit}_xy")


with program() as single_qubit_process_tomography:

n = declare(int)
n_st = declare_stream()

I = declare(fixed)
Q = declare(fixed)

state = declare(bool) # QUA variable for the measured qubit state
state_st = declare_stream() # Stream for the qubit state

c = declare(int) # QUA variable for switching between state preparation/creations
m = declare(int) # QUA variable for switching between basis projections/measurements

with for_(n, 0, n < n_avg, n + 1): # QUA for_ loop for averaging
with for_(c, 0, c <= 5, c + 1): # QUA for_ loop for switching between state preparations
with for_(m, 0, m <= 5, m + 1): # QUA for_ loop for switching between basis projections/measurements

# prepare qubit in one of six cardinal Bloch sphere states
prepare_state(c)

align()

# the process to be analysed
analysed_process(qubit)

align()

# projective measurement basis change
measurement_basis_change(m)

align()

measure(
"readout",
f"rr{qubit}",
None,
dual_demod.full("rotated_cos", "out1", "rotated_sin", "out2", I),
dual_demod.full("rotated_minus_sin", "out1", "rotated_cos", "out2", Q),
)

align()

# track number of times the qubit is found to be
# in the ground state
assign(state, I < threshold)

save(state, state_st)

wait(thermalization_time * u.ns)

save(n, n_st)

with stream_processing():
n_st.save("iteration")
state_st.boolean_to_int().buffer(6).buffer(6).average().save("probs")


#####################################
# Open Communication with the QOP #
#####################################

qmm = QuantumMachinesManager(host=qop_ip, port=qop_port, cluster_name=cluster_name, octave=octave_config)

###########################
# Run or Simulate Program #
###########################

simulate = False

if simulate:
# Simulates the QUA program for the specified duration
simulation_config = SimulationConfig(duration=10_000) # In clock cycles = 4ns
job = qmm.simulate(config, single_qubit_process_tomography, simulation_config)
job.get_simulated_samples().con1.plot()

else:
# Open the quantum machine
qm = qmm.open_qm(config)
# Send the QUA program to the OPX, which compiles and executes it
job = qm.execute(single_qubit_process_tomography)
# Get results from QUA program
results = fetching_tool(job, data_list=["iteration", "probs"], mode="live")

while results.is_processing():

# Fetch results
iteration, probs = results.fetch_all()

# Progress bar
progress_counter(iteration, n_avg, start_time=results.get_start_time())

# Close the quantum machines at the end in order to put all flux biases to 0 so that the fridge doesn't heat-up
qm.close()

# post-processing
PMatrix = np.array(
[
[
P_Pauli1(
np.floor(v / 4).astype(int),
v % 4,
np.floor(w / 4).astype(int),
w % 4,
)
for w in range(16)
]
for v in range(16)
]
)

pauli_basis_measurements = lambda l, k: map_from_bloch_state_to_pauli_basis1(l, k, probs)

measurement_vector = np.array([pauli_basis_measurements(np.floor(v / 4).astype(int), v % 4) for v in range(16)])

chi_vector = solve(PMatrix, measurement_vector)
chi_matrix = chi_vector.reshape(4, 4)

plot_process_tomography1(chi_vector)
Original file line number Diff line number Diff line change
@@ -0,0 +1,174 @@
#!/usr/bin/env python

"""
SINGLE QUBIT STATE TOMOGRAPHY
The sequence consists of preparing the qubit into a chosen state using calibrated gates from the config, and
measuring the state of the qubit, by way of the readout resonator, in the X, Y and Z bases. The output,
'probs', is a triple of probability 'difference' counts for measuring the qubit in one eigenstate of the
X/Y/Z bases, minus the probability of measuring it in that eigenstate's complement.
These values are scaled into the usual Stokes parameters and used to infer the state of the qubit; see
https://research.physics.illinois.edu/QI/Photonics/tomography-files/tomo_chapter_2004.pdf for further details
Note that this program is similar to
qua-libs/Quantum-Control-Applications/Superconducting/Single-Fixed-Transmon
/19_state_tomography.py, which the author became aware of after having
written the current script
Prerequisites:
- Having found the resonance frequency of the resonator coupled to the qubit under study (resonator_spectroscopy)
- Having calibrated qubit pi and pi/2 pulses by running qubit spectroscopy, rabi_chevron, power_rabi and updating the config
- Having calibrated the readout (readout_frequency, amplitude, duration_optimization IQ_blobs) for better SNR, and having
saved the derived readout threshold value in the config
- Set the desired flux bias in the case of flux-tunable qubits
This script implements the logic outlined in the notebook
https://github.com/bornman-nick/quantum-state-and-process-tomography. See that
notebook for a detailed explanation of standard single qubit state
tomography.
"""

from qm.qua import *
from qm import QuantumMachinesManager
from qm import SimulationConfig
from configuration import *
from qualang_tools.results import progress_counter, fetching_tool
import numpy as np


###################
# The QUA program #
###################

# qubit under test, assuming there are multiple qubits
# with XY line elements "q<qubit>_xy", flux lines
# "q<qubit>_z", and readout resonator elements "rr<qubit>"
qubit = 1

if qubit == 1:
threshold = ge_threshold_q1
elif qubit == 2:
threshold = ge_threshold_q2
else:
raise ValueError(f"Incorrect qubit number chosen")


n_avg = 10_000


# subroutine to prepare desired qubit state
def prepare_state(qubit):
# write whatever QUA code you need in order to create the
# state to perform tomography on, from an initial
# ground state. For example, to create the |1> state
play("y180", f"q{qubit}_xy")


with program() as single_qubit_state_tomography:
n = declare(int) # QUA variable for average loop
n_st = declare_stream() # Stream for the averaging iteration 'n'

state = declare(bool) # QUA variable for the qubit state
state_st = declare_stream() # Stream for the qubit state

p = declare(int) # QUA variable for switching between projections

I = declare(fixed)
Q = declare(fixed)

with for_(n, 0, n < n_avg, n + 1): # QUA for_ loop for averaging
with for_(p, 0, p <= 2, p + 1): # QUA for_ loop for switching between basis changes

prepare_state(qubit)
align()

with switch_(c):
with case_(0): # basis X

# Map the X-component of the Bloch vector onto the Z-axis
# 1/sqrt(2)(|0>+|1>) -> |0>; 1/sqrt(2)(|0>-|1>) -> |1>
play("-y90", f"q{qubit}_xy")

align(f"q{qubit}_xy", f"rr{qubit}")

with case_(1): # basis Y

# Map the Y-component of the Bloch vector onto the Z-axis
# 1/sqrt(2)(|0>+i|1>) -> |0>; 1/sqrt(2)(|0>-i|1>) -> |1>
play("x90", f"q{qubit}_xy")

align(f"q{qubit}_xy", f"rr{qubit}")

with case_(2): # basis Z

align(f"q{qubit}_xy", f"rr{qubit}")

measure(
"readout",
f"rr{qubit}",
None,
dual_demod.full("rotated_cos", "out1", "rotated_sin", "out2", I),
dual_demod.full("rotated_minus_sin", "out1", "rotated_cos", "out2", Q),
)

# True if qubit state is |1>, False if |0>

assign(state, I > threshold)

wait(thermalization_time * u.ns)

save(state, state_st)

save(n, n_st)

with stream_processing():
n_st.save("iteration")
state_st.boolean_to_int().buffer(3).average().save("probs")


#####################################
# Open Communication with the QOP #
#####################################
qmm = QuantumMachinesManager(host=qop_ip, port=qop_port, cluster_name=cluster_name, octave=octave_config)

###########################
# Run or Simulate Program #
###########################
simulate = False

if simulate:
# Simulates the QUA program for the specified duration
simulation_config = SimulationConfig(duration=10_000) # In clock cycles = 4ns
job = qmm.simulate(config, single_qubit_state_tomography, simulation_config)
job.get_simulated_samples().con1.plot()

else:
# Open the quantum machine
qm = qmm.open_qm(config)
# Send the QUA program to the OPX, which compiles and executes it
job = qm.execute(single_qubit_state_tomography)
# Get results from QUA program
results = fetching_tool(job, data_list=["probs", "iteration"], mode="live")

while results.is_processing():
# Fetch results
probs, iteration = results.fetch_all()
# Progress bar
progress_counter(iteration, n_avg, start_time=results.get_start_time())

# Converts the (0,1) -> |g>,|e> convention, arising from the I>threshold
# assignment, to (1,-1) -> |g>,|e>, which aligns with the Stokes parameter
# definitions from the projector probabilities
prob = -2 * (probs - 0.5)

# Close the quantum machines at the end in order to put all flux biases to 0 so that the fridge doesn't heat-up
qm.close()

# Reconstruct the density matrix
I = np.array([[1, 0], [0, 1]])
sigma_x = np.array([[0, 1], [1, 0]])
sigma_y = np.array([[0, -1j], [1j, 0]])
sigma_z = np.array([[1, 0], [0, -1]])

rho = 0.5 * (I + prob[0] * sigma_x + prob[1] * sigma_y + prob[2] * sigma_z)
print(f"The density matrix is:\n{rho}")
Original file line number Diff line number Diff line change
@@ -0,0 +1,331 @@
#!/usr/bin/env python

"""
TWO QUBIT PROCESS TOMOGRAPHY
The sequence consists of preparing the qubits each into one of their six cardinal Bloch sphere
states using calibrated gates from the config, applying the operation/process - typically
ideally a unitary matrix - under investigation, and then measuring the state of the qubits,
by way of the readout resonator, in their X, Y and Z bases (specifically, measuring the
projector of one of the same six cardinal Bloch sphere states). This is repeated for
the full set of input states and measurement projectors, which is tomographically
complete. The output, 'probs', containing the 'Bloch sphere' prepared states and
measurement projectors, is mapped to the Pauli basis (giving a 'measurement vector'),
and the equation "measurement_vector = B \times chi_vector" (where B is a constant)
inverted to give the chi process matrix for the process under
investigation.
Computing PMatrix from scratch in the two-qubit process case takes a while,
so included with this project is a serialised file of Pmatrix
Prerequisites:
- Having found the resonance frequency of the resonator coupled to the qubits under study (resonator_spectroscopy)
- Having calibrated qubits' pi and pi/2 pulses by running qubit spectroscopy, rabi_chevron, power_rabi and updating the config
- Having calibrated the readout (readout_frequency, amplitude, duration_optimization IQ_blobs) for better SNR, and having
saved the derived readout threshold values in the config
- Set the desired flux biases in the case of flux-tunable qubits
This script implements the logic outlined in the notebook
https://github.com/bornman-nick/quantum-state-and-process-tomography. See that
notebook for a detailed explanation of standard two qubit process tomography
Note that this routine carries out standard process tomography for two qubits.
Given that the measurements recorded are simply whether the final state of the
qubits are in their ground state or not for each state preparation and measurement
setting (instead of also deciding whether the qubits are perhaps in the excited
states and incrementing the counter of a complementary measurement setting),
this routine scales as 6**4, rather than (6**2)*(3**2), as is the case in
most formulations of standard process tomography. However, the routine below
is fast enough that this doesn't matter too much.
"""

from qm.qua import *
from qm import QuantumMachinesManager
from qm import SimulationConfig
from configuration import *
from qualang_tools.results import progress_counter, fetching_tool
import numpy as np
from scipy.linalg import solve

import os
import pickle

from helper_functions import (
P_Pauli2,
plot_process_tomography2,
map_from_bloch_state_to_pauli_basis2,
func_F2,
)


###################
# The QUA program #
###################

# load PMatrix if pickle file is present
# note: this takes a while to compute if you do not load the serialised
# PMatrix file
file_location = "<path to PMatrix2 file>" # ./PMatrix2.pkl

if os.path.isfile(file_location):
with open(file_location, "rb") as file:
PMatrix = pickle.load(file)
else:
# computing this takes a while
PMatrix = np.array(
[
[
P_Pauli2(
np.floor(v / 16).astype(int),
v % 16,
np.floor(w / 16).astype(int),
w % 16,
)
for w in range(16**2)
]
for v in range(16**2)
]
)

# qubits under test, assuming there are multiple qubits
# with XY line elements "q<qubit>_xy", flux lines
# "q<qubit>_z", and readout resonator elements "rr<qubit>"
qubit1 = 1
qubit2 = 2

if qubit1 == 1:
threshold1 = ge_threshold_q1
elif qubit1 == 2:
threshold1 = ge_threshold_q2
else:
raise ValueError(f"Incorrect qubit1 number chosen")

if qubit2 == 1:
threshold2 = ge_threshold_q1
elif qubit2 == 2:
threshold2 = ge_threshold_q2
else:
raise ValueError(f"Incorrect qubit2 number chosen")

if qubit1 == qubit2:
raise ValueError(f"The value of qubit1 cannot equal that of qubit2")


n_avg = 5_000


# subroutine to prepare desired gate/process
def analysed_process(qubit1, qubit2):
# write whatever QUA code you need, here, in order to perform the
# desired process which we want to subject to tomography.
# For example, to analyse the X gate on qubit1 and the -Y90
# gate on qubit2:
play("x180", f"q{qubit1}_xy")
play("-y90", f"q{qubit2}_xy")


ideal_gate = np.kron(func_F2(1), func_F2(3))


# subroutine to switch between preparing the tomographically-complete
# set of input states
def prepare_states(i, j, qubit1, qubit2):

# Prepare Bloch state i on qubit1
with switch_(i):
with case_(0):
wait(pi_len // 4, f"q{qubit1}_xy")
with case_(1):
play("x180", f"q{qubit1}_xy")
with case_(2):
play("y90", f"q{qubit1}_xy")
with case_(3):
play("-y90", f"q{qubit1}_xy")
with case_(4):
play("-x90", f"q{qubit1}_xy")
with case_(5):
play("x90", f"q{qubit1}_xy")

# Prepare Bloch state j on qubit2
with switch_(j):
with case_(0):
wait(pi_len // 4, f"q{qubit2}_xy")
with case_(1):
play("x180", f"q{qubit2}_xy")
with case_(2):
play("y90", f"q{qubit2}_xy")
with case_(3):
play("-y90", f"q{qubit2}_xy")
with case_(4):
play("-x90", f"q{qubit2}_xy")
with case_(5):
play("x90", f"q{qubit2}_xy")


# subroutine to switch between the tomographically-complete
# set of bases in which to measure
def measurement_basis_change(k, l, qubit1, qubit2):

# Change basis with operation k on qubit1
with switch_(k):
with case_(0):
wait(pi_len // 4, f"q{qubit2}_xy")
with case_(1):
play("x180", f"q{qubit1}_xy")
with case_(2):
play("y90", f"q{qubit1}_xy")
with case_(3):
play("-y90", f"q{qubit1}_xy")
with case_(4):
play("-x90", f"q{qubit1}_xy")
with case_(5):
play("x90", f"q{qubit1}_xy")

# Change basis with operation l on qubit2
with switch_(l):
with case_(0):
wait(pi_len // 4, f"q{qubit2}_xy")
with case_(1):
play("x180", f"q{qubit2}_xy")
with case_(2):
play("y90", f"q{qubit2}_xy")
with case_(3):
play("-y90", f"q{qubit2}_xy")
with case_(4):
play("-x90", f"q{qubit2}_xy")
with case_(5):
play("x90", f"q{qubit2}_xy")


with program() as two_qubit_process_tomography:

n = declare(int)
n_st = declare_stream()

I1 = declare(fixed)
Q1 = declare(fixed)
I2 = declare(fixed)
Q2 = declare(fixed)

state = declare(bool) # QUA variable for the measured qubits state
state_st = declare_stream() # Stream for the qubits state

c1 = declare(int) # QUA variable for switching between state preparation/creations on qubit 1
c2 = declare(int) # QUA variable for switching between state preparation/creations on qubit 2
m1 = declare(int) # QUA variable for switching between Bloch basis projections/measurements on qubit 1
m2 = declare(int) # QUA variable for switching between Bloch basis projections/measurements on qubit 2

with for_(n, 0, n < n_avg, n + 1): # QUA for_ loop for averaging
with for_(c1, 0, c1 <= 5, c1 + 1): # QUA for_ loop for switching between state preparations on qubit 1
with for_(c2, 0, c2 <= 5, c2 + 1): # QUA for_ loop for switching between state preparations on qubit 2
with for_(
m1, 0, m1 <= 5, m1 + 1
): # QUA for_ loop for switching between Bloch basis projections/measurements on qubit 1
with for_(
m2, 0, m2 <= 5, m2 + 1
): # QUA for_ loop for switching between Bloch basis projections/measurements on qubit 2

reset_frame(f"q{qubit1}_xy")
reset_frame(f"q{qubit2}_xy")

reset_phase(f"q{qubit1}_xy")
reset_phase(f"q{qubit2}_xy")

# prepare qubit1 and qubit 2 in one of six Bloch sphere states
prepare_states(c1, c2, qubit1, qubit2)

align()

# apply the process to be analysed
analysed_process(qubit1, qubit2)

aling()

# projective measurement basis change
measurement_basis_change(m1, m2, qubit1, qubit2)

align()

measure(
"readout",
f"rr{qubit1}",
None,
dual_demod.full("rotated_cos", "out1", "rotated_sin", "out2", I1),
dual_demod.full("rotated_minus_sin", "out1", "rotated_cos", "out2", Q1),
)

measure(
"readout",
f"rr{qubit2}",
None,
dual_demod.full("rotated_cos", "out1", "rotated_sin", "out2", I2),
dual_demod.full("rotated_minus_sin", "out1", "rotated_cos", "out2", Q2),
)

align()

# track the number of clicks when both qubit 1 and qubit 2 are
# in their ground states
with if_((I1 < threshold1) & (I2 < threshold2)):
assign(state, True)
with else_():
assign(state, False)

save(state, state_st)

wait(thermalization_time * u.ns, f"rr{qubit1}", f"rr{qubit2}")

save(n, n_st)

with stream_processing():
n_st.save("iteration")
state_st.boolean_to_int().buffer(6).buffer(6).buffer(6).buffer(6).average().save("probs")


#####################################
# Open Communication with the QOP #
#####################################

qmm = QuantumMachinesManager(host=qop_ip, port=qop_port, cluster_name=cluster_name, octave=octave_config)

###########################
# Run or Simulate Program #
###########################

simulate = False

if simulate:
# Simulates the QUA program for the specified duration
simulation_config = SimulationConfig(duration=10_000) # In clock cycles = 4ns
job = qmm.simulate(config, two_qubit_process_tomography, simulation_config)
job.get_simulated_samples().con1.plot()

else:
# Open the quantum machine
qm = qmm.open_qm(config)
# Send the QUA program to the OPX, which compiles and executes it
job = qm.execute(two_qubit_process_tomography)
# Get results from QUA program
results = fetching_tool(job, data_list=["iteration", "probs"], mode="live")

while results.is_processing():

# Fetch results
iteration, probs = results.fetch_all()

# Progress bar
progress_counter(iteration, n_avg, start_time=results.get_start_time())

# Close the quantum machines at the end in order to put all flux biases to 0 so that the fridge doesn't heat-up
qm.close()

# post-processing
pauli_basis_measurements = lambda q, n: map_from_bloch_state_to_pauli_basis2(q, n, probs)

measurement_vector = np.array(
[pauli_basis_measurements(np.floor(v / 16).astype(int), v % 16) for v in range(16**2)]
)

chi_vector = solve(PMatrix, measurement_vector)
chi_matrix = chi_vector.reshape(16, 16)

plot_process_tomography2(chi_vector)
Original file line number Diff line number Diff line change
@@ -0,0 +1,349 @@
#!/usr/bin/env python

"""
TWO QUBIT STATE TOMOGRAPHY
The sequence consists of preparing the qubits into a chosen state using
calibrated gates from the config, and measuring the state of the qubits, by way
of the readout resonator, in their X, Y and Z bases. The output, 'probs', is a
list of 15 probability 'difference' counts corresponding with the 15 Stokes
parameters for two qubit state tomography. These values are used to infer the
two-qubit state; see
https://research.physics.illinois.edu/QI/Photonics/tomography-files/tomo_chapter_2004.pdf
for further details.
Prerequisites:
- Having found the resonance frequencies of the resonators coupled to the
qubits under study (resonator_spectroscopy)
- Having calibrated qubits' pi and pi/2 pulses by running qubit
spectroscopy, rabi_chevron, power_rabi and updating the config
- Having calibrated the readout (readout_frequency, amplitude,
duration_optimization IQ_blobs) for each qubit
for better SNR, and having saved the derived readout threshold values in
the config
- Set the desired flux biases in the case of flux-tunable qubits
This script implements the logic outlined in the notebook
https://github.com/bornman-nick/quantum-state-and-process-tomography. See that
notebook for a detailed explanation of standard two-qubit state tomography
"""

import numpy as np
from scipy.linalg import sqrtm

from qm.qua import *
from qm import QuantumMachinesManager
from qm import SimulationConfig
from configuration import *
from qualang_tools.results import progress_counter, fetching_tool

from helper_functions import rotated_multiplexed_state_discrimination


###################
# The QUA program #
###################

# qubits under test, assuming there are multiple qubits
# with XY line elements "q<qubit>_xy", flux lines
# "q<qubit>_z", and readout resonator elements "rr<qubit>"

qubit1 = 1
qubit2 = 2

if qubit1 == 1:
threshold1 = ge_threshold_q1
elif qubit1 == 2:
threshold1 = ge_threshold_q2
else:
raise ValueError(f"Incorrect qubit1 number chosen")

if qubit2 == 1:
threshold2 = ge_threshold_q1
elif qubit2 == 2:
threshold2 = ge_threshold_q2
else:
raise ValueError(f"Incorrect qubit2 number chosen")

if qubit1 == qubit2:
raise ValueError(f"The value of qubit1 cannot equal that of qubit2")


n_avg = 10_000


# subroutine to prepare desired qubit state
def prepare_state(qubit1, qubit2):
# write whatever QUA code you need in order to create the
# state to perform tomography on, from an initial
# ground state. For example, to create the |1> 1/sqrt(2)(|0>+i|1>)
# state
play("y180", f"q{qubit1}_xy")
play("-x90", f"q{qubit2}_xy")


# matrix representation of the ideal composite qubit state from above
# (|1> 1/sqrt(2)(|0>+i|1>)) (<1| 1/sqrt(2)(<0|-i<1|))
ideal_state = np.array([[0, 0, 0, 0], [0, 0, 0, 0], [0, 0, 1 / 2, 1j / 2], [0, 0, -1j / 2, 1 / 2]])


with program() as two_qubit_state_tomography:

n = declare(int)
n_st = declare_stream()
I = [declare(fixed) for _ in range(2)]
Q = [declare(fixed) for _ in range(2)]
# I_st = [declare_stream() for _ in range(2)]
# Q_st = [declare_stream() for _ in range(2)]

states = [
declare(bool),
declare(bool),
] # QUA variable for the measured qubit states
# states_st = [declare_stream(), declare_stream()] # Stream for the qubits states

p00 = declare(int) # variable to track number of |0>|0> measurements
p01 = declare(int) # variable to track number of |0>|1> measurements
p10 = declare(int) # variable to track number of |1>|0> measurements
p11 = declare(int) # variable to track number of |1>|1> measurements

prob_vec_results_st = declare_stream() # Stream for the probability vector - average of states vector

c = declare(int) # QUA variable for switching between projections

with for_(n, 0, n < n_avg, n + 1): # QUA for_ loop for averaging
with for_(c, 1, c <= 15, c + 1): # QUA for_ loop for switching between projections

prepare_state(qubit1, qubit2)
align()

assign(p00, 0)
assign(p01, 0)
assign(p10, 0)
assign(p11, 0)

with switch_(c):
with case_(1): # projection along Z1, X2

play("-y90", f"q{qubit2}_xy")
# Stokes parameter for this case: P00 + P10 - P01 - P11

with case_(2): # projection along Z1, Y2

play("x90", f"q{qubit2}_xy")
# Stokes parameter for this case: P00 + P10 - P01 - P11

with case_(3): # projection along Z1, Z2

wait(pi_len // 4, f"q{qubit1}_xy", f"q{qubit2}_xy")
# Stokes parameter for this case: P00 + P10 - P01 - P11

with case_(4): # projection along X1, Z2

play("-y90", f"q{qubit1}_xy")
# Stokes parameter for this case: P00 + P01 - P10 - P11

with case_(5): # projection along X1, X2

play("-y90", f"q{qubit1}_xy")
play("-y90", f"q{qubit2}_xy")
# Stokes parameter for this case: P00 + P11 - P10 - P01

with case_(6): # projection along X1, Y2

play("-y90", f"q{qubit1}_xy")
play("x90", f"q{qubit2}_xy")
# Stokes parameter for this case: P00 + P11 - P10 - P01

with case_(7): # projection along X1, Z2

play("-y90", f"q{qubit1}_xy")
# Stokes parameter for this case: P00 + P11 - P10 - P01

with case_(8): # projection along Y1, Z2

play("x90", f"q{qubit1}_xy")
# Stokes parameter for this case: P00 + P01 - P10 - P11

with case_(9): # projection along Y1, X2

play("x90", f"q{qubit1}_xy")
play("-y90", f"q{qubit2}_xy")
# Stokes parameter for this case: P00 + P11 - P10 - P01

with case_(10): # projection along Y1, Y2

play("x90", f"q{qubit1}_xy")
play("x90", f"q{qubit2}_xy")
# Stokes parameter for this case: P00 + P11 - P10 - P01

with case_(11): # projection along Y1, Z2

play("x90", f"q{qubit1}_xy")
# Stokes parameter for this case: P00 + P11 - P10 - P01

with case_(12): # projection along Z1, Z2

wait(pi_len // 4, f"q{qubit1}_xy", f"q{qubit2}_xy")
# Stokes parameter for this case: P00 + P01 - P10 - P11

with case_(13): # projection along Z1, X2

play("-y90", f"q{qubit2}_xy")
# Stokes parameter for this case: P00 + P11 - P10 - P01

with case_(14): # projection along Z1, Y2

play("x90", f"q{qubit2}_xy")
# Stokes parameter for this case: P00 + P11 - P10 - P01

with case_(15): # projection along Z1, Z2

wait(pi_len // 4, f"q{qubit1}_xy", f"q{qubit2}_xy")
# Stokes parameter for this case: P00 + P11 - P10 - P01

align()

rotated_multiplexed_state_discrimination(
I,
None,
Q,
None,
states,
None,
[qubit1, qubit2],
[threshold1, threshold2],
)

align()

# populate correct variable depending on which qubit states
# were measured
with if_(states[0]):
with if_(states[1]):
assign(p11, 1) # |1>|1> was measured
with else_():
assign(p10, 1) # |1>|0> was measured
with else_():
with if_(states[1]):
assign(p01, 1) # |0>|1> was measured
with else_():
assign(p00, 1) # |0>|0> was measured

# stream all four possible values - only a single one of these
# four should be 1, the rest, 0
save(p00, prob_vec_results_st)
save(p01, prob_vec_results_st)
save(p10, prob_vec_results_st)
save(p11, prob_vec_results_st)

wait(thermalization_time * u.ns)

save(n, n_st)

with stream_processing():
n_st.save("iteration")
prob_vec_results_st.buffer(4).buffer(15).average().save("probs")


#####################################
# Open Communication with the QOP #
#####################################
qmm = QuantumMachinesManager(host=qop_ip, port=qop_port, cluster_name=cluster_name, octave=octave_config)

###########################
# Run or Simulate Program #
###########################
simulate = False

if simulate:
# Simulates the QUA program for the specified duration
simulation_config = SimulationConfig(duration=10_000) # In clock cycles = 4ns
job = qmm.simulate(config, two_qubit_state_tomography, simulation_config)
job.get_simulated_samples().con1.plot()

else:
# Open the quantum machine
qm = qmm.open_qm(config)
# Send the QUA program to the OPX, which compiles and executes it
job = qm.execute(two_qubit_state_tomography)
# Get results from QUA program
results = fetching_tool(job, data_list=["probs", "iteration"], mode="live")

while results.is_processing():

# Fetch results
iteration, probs = results.fetch_all()

# Progress bar
progress_counter(iteration, n_avg, start_time=results.get_start_time())

# Close the quantum machines at the end in order to put all flux biases to
# 0 so that the fridge doesn't heat-up
qm.close()

# Use probs to reconstruct density matrix

# initialise vector to contain the 15 stokes parameters
stokes = np.zeros(15)

# For cases 1, 2,and 3 - Stokes parameter is P00 + P10 - P01 - P11
# For cases 4, 8 and 12 - Stokes parameter is P00 + P01 - P10 - P11
# For cases 5, 6, 7, 9, 10, 11, 13, 14, 15 - P00 + P11 - P10 - P01
for i in range(1, 16):
if i in [1, 2, 3]:
stokes[i - 1] = probs[i - 1][0] + probs[i - 1][2] - probs[i - 1][1] - probs[i - 1][3]
elif i in [4, 8, 12]:
stokes[i - 1] = probs[i - 1][0] + probs[i - 1][1] - probs[i - 1][2] - probs[i - 1][3]
elif i in [5, 6, 7, 9, 10, 11, 13, 14, 15]:
stokes[i - 1] = probs[i - 1][0] + probs[i - 1][3] - probs[i - 1][1] - probs[i - 1][2]

# Derive the density matrix
I = np.array([[1, 0], [0, 1]])
sigma_x = np.array([[0, 1], [1, 0]])
sigma_y = np.array([[0, -1j], [1j, 0]])
sigma_z = np.array([[1, 0], [0, -1]])

# Density matrix Pauli operator basis
II = np.kron(I, I)
IX = np.kron(I, sigma_x)
IY = np.kron(I, sigma_y)
IZ = np.kron(I, sigma_z)
XI = np.kron(sigma_x, I)
XX = np.kron(sigma_x, sigma_x)
XY = np.kron(sigma_x, sigma_y)
XZ = np.kron(sigma_x, sigma_z)
YI = np.kron(sigma_y, I)
YX = np.kron(sigma_y, sigma_x)
YY = np.kron(sigma_y, sigma_y)
YZ = np.kron(sigma_y, sigma_z)
ZI = np.kron(sigma_z, I)
ZX = np.kron(sigma_z, sigma_x)
ZY = np.kron(sigma_z, sigma_y)
ZZ = np.kron(sigma_z, sigma_z)

rho = 0.25 * (
II
+ stokes[0] * IX
+ stokes[1] * IY
+ stokes[2] * IZ
+ stokes[3] * XI
+ stokes[4] * XX
+ stokes[5] * XY
+ stokes[6] * XZ
+ stokes[7] * YI
+ stokes[8] * YX
+ stokes[9] * YY
+ stokes[10] * YZ
+ stokes[11] * ZI
+ stokes[12] * ZX
+ stokes[13] * ZY
+ stokes[14] * ZZ
)

print(f"The density matrix is:\n{np.round(rho, decimals=3)}")

sqrt_ideal_state = sqrtm(ideal_state)

state_fidelity = (np.abs(sqrtm(sqrt_ideal_state @ rho @ sqrt_ideal_state).trace())) ** 2

print(f"The state fidelity is: {np.round(state_fidelity, decimals=4)}")
1 change: 1 addition & 0 deletions pyproject.toml
Original file line number Diff line number Diff line change
@@ -17,6 +17,7 @@ flake8 = "^3.9.1"

[tool.black]
line-length = 120
#force-exclude = '''Two-Flux-Tunable-Transmons/Use Case 4 - Single- and Two-Qubit State and Process Tomography/helper_functions.py'''
#include = '(examples)/.+\.pyi?$'
#include = '(Quantum Control Applications)/.+\.pyi?$'
#include = '(Tutorials)/.+\.pyi?$'