Using Quara’s standard tomography features from Forest (1-qubit)
Quara supports object conversions for executing standard quantum tomography from Forest SDK. Here we briefly explain how to perform quantum tomography by using Quara and measurements obtained from quantum virtual machine of Forest SDK.
[1]:
import numpy as np
from numpy import pi
# quara
from quara.interface.forest.api import (
generate_preprocess_program,
generate_pauli_strings_from_povm_name,
calc_empi_dist_from_observables
)
from quara.objects.composite_system_typical import generate_composite_system
from quara.objects.povm_typical import generate_povm_from_name
from quara.objects.state_typical import generate_state_from_name
from quara.protocol.qtomography.standard.standard_qst import StandardQst
from quara.protocol.qtomography.standard.standard_qpt import StandardQpt
from quara.protocol.qtomography.standard.linear_estimator import LinearEstimator
# pyquil
from pyquil import get_qc, Program
from pyquil.gates import H, PHASE
from pyquil.paulis import PauliTerm
from pyquil.experiment import (
Experiment,
ExperimentSetting,
zeros_state,
)
from pyquil.operator_estimation import measure_observables
Quantum State Tomography (QST)
First, we define a program for state preparation for tomographic experiment. Then, we create a function obtain_expectations_for_qst() which is a function that performs measuments based on given pauli operators. Pauli operators are operators that can be performed as measurements in Forest SDK. More details will be explained in the next cell. Make sure you have started quilc and qvm on your other consoles. See this
page for more information.
[2]:
# example for 1 qubits system
qc = get_qc("1q-qvm")
qubits = [0]
# define initialization of the quantum system
# creating a state called "a state"
num_shots = 10000
p = Program()
p += H(qubits[0])
p += PHASE(pi/4, qubits[0])
p.wrap_in_numshots_loop(num_shots)
def obtain_expectations_for_qst(qc, program, pauli_strings):
settings = []
for pauli_str in pauli_strings:
out_operator = PauliTerm.from_list(list(zip(pauli_str, qubits)))
settings.append(ExperimentSetting(zeros_state(qubits), out_operator))
tomo_experiment = Experiment(settings, program)
expectations = []
for pauli_str, res in zip(
pauli_strings,
measure_observables(
qc,
tomo_experiment,
),
):
if res.raw_expectation is None:
# This is the result for II...I operator
expectations.append(1.0)
else:
expectations.append(res.raw_expectation)
print(f"Raw expectation of {pauli_str}: {expectations[-1]}")
return expectations
Next, we define a composite system and POVM objects which will make a tomographically complete measuring experiment. Then we calculate a set of desired pauli mesruments by using a pre-difned function named generate_pauli_strings_from_povm_name() in Quara for each POVM. In this example, we create a dictionary for POVM objects and Pauli operator names which both are indexed with povm_name.
[3]:
c_sys = generate_composite_system("qubit", 1)
povms = {}
pauli_strings = {}
povm_names = ["x", "y", "z"]
for povm_name in povm_names:
povms[povm_name] = generate_povm_from_name(povm_name, c_sys)
pauli_strings[povm_name] = generate_pauli_strings_from_povm_name(povm_name)
print(f"Sample of Quara tester object for povm_name: {povm_names[0]}")
print(povms[povm_names[0]])
print("\nPauli strings that will be used in this tomographical experiment:")
pauli_strings
Sample of Quara tester object for povm_name: x
Type:
Povm
Dim:
2
Number of outcomes:
2
Vecs:
[[ 0.70710678 0.70710678 0. 0. ]
[ 0.70710678 -0.70710678 0. 0. ]]
Pauli strings that will be used in this tomographical experiment:
[3]:
{'x': ['X', 'I'], 'y': ['Y', 'I'], 'z': ['Z', 'I']}
Now we perform projective measurements based on Pauli operators and obtain_expectations_for_qst() which we defined in previous cells. Observed expectations are also stored as dictionary and indexed with povm_name.
[4]:
observables = {}
for povm_name in povm_names:
observables[povm_name] = obtain_expectations_for_qst(qc, p, generate_pauli_strings_from_povm_name(povm_name))
print(f'Expectations for {povm_names[0]} POVM: {observables[povm_names[0]]}')
Raw expectation of X: 0.6976
Raw expectation of I: 1.0
Raw expectation of Y: 0.6948
Raw expectation of I: 1.0
Raw expectation of Z: 0.0024
Raw expectation of I: 1.0
Expectations for x POVM: [0.6976, 1.0]
We need empirical distributions that correspond to selected POVMs in order to perform QST. In ideal case of projective measurements, probability distribution of a POVM can be calculated from a set of Pauli observables and coefficient matrix \(A_{\mathrm{coefficient}}\).
Here we create a set of empirical distributions empi_dists, which are pairs of a repetition number and relative frequencies. We calculate empi_dists from measured expectations by using calc_empi_dist_from_observables() function.
[5]:
empi_dists = []
for povm_name in povm_names:
empi_dist = calc_empi_dist_from_observables(observables[povm_name], num_shots, pauli_strings[povm_name], povms[povm_name])
empi_dists.append(empi_dist)
print(f"{povm_name}: {empi_dist}")
x: (10000, array([0.8488, 0.1512]))
y: (10000, array([0.8474, 0.1526]))
z: (10000, array([0.5012, 0.4988]))
We generate a StandardQst object. From now on, the flow of QST will be the same as explained in tutorial_for_standardqtomography.
[6]:
tester_povms = []
for povm_name in povm_names:
tester_povms.append(povms[povm_name])
qst = StandardQst(tester_povms, on_para_eq_constraint=True, schedules="all")
When we choose a linear estimator, the estimate is calculated as follows
[7]:
estimator = LinearEstimator()
result = estimator.calc_estimate(qtomography=qst, empi_dists=empi_dists, is_computation_time_required=True)
estimate = result.estimated_qoperation
print(estimate)
print("is estimate physical? : ", estimate.is_physical())
print("\nEigenvalues are: ", estimate.calc_eigenvalues())
Type:
State
Dim:
2
Vec:
[0.70710678 0.49327769 0.49129779 0.00169706]
is estimate physical? : True
Eigenvalues are: [0.9922901989680475, 0.007709801031952318]
An eigenvalue of the estimated density matrix is negative, which violates the requirement of positive-semidefiniteness on density matrix. This kind of violation can occur when we choose a linear estimator. In order to avoid the problem, we need to perform a constraint optimization at the data-processing.
Quantum Process Tomography
In this section, we consider process tomography on 1-qubit system. First, we define a target_program which is the quantum process we want to investigate. We also define a function obtain_expectations_for_qpt which is the function that calls pyquil backend to obtain observations.
When performing process tomography, we need a pre process for state preparation which will be defined as preprocess_program. After executing preprocess_program, we run the target_program and then we measure the expectations.
[8]:
# example for 1 qubits system
qc = get_qc("1q-qvm")
qubits = [0]
num_shots = 10000
# define a process that we want to estimate
target_program = Program()
target_program += PHASE(pi/4, qubits[0])
def obtain_expectations_for_qpt(qc, num_shots, preprocess_program, program, pauli_strings):
settings = []
for pauli_str in pauli_strings:
out_operator = PauliTerm.from_list(list(zip(pauli_str, qubits)))
settings.append(ExperimentSetting(zeros_state(qubits), out_operator))
tomo_experiment_program = preprocess_program + program
tomo_experiment_program.wrap_in_numshots_loop(num_shots)
tomo_experiment = Experiment(settings, tomo_experiment_program)
expectations = []
for pauli_str, res in zip(
pauli_strings,
measure_observables(
qc,
tomo_experiment,
),
):
if res.raw_expectation is None:
# This is the result for II...I operator
expectations.append(1.0)
else:
expectations.append(res.raw_expectation)
print(f"Raw expectation of {pauli_str}: {expectations[-1]}")
return expectations
In quara, the name of a quantum system is expressed using strings such as x0 or z0_y0 in multiple qubit systems. There are 6 basic quantum state that are widely used in 1 qubit system which are the following:
We create a state object in Quara and then generate programs to realize those states using generate_preprocess_program.
[9]:
c_sys = generate_composite_system("qubit", 1)
# state names of tester states
state_names = ["x0", "y0", "z0", "z1"]
states = {}
preprocess_programs = {}
for state_name in state_names:
states[state_name] = generate_state_from_name(c_sys, state_name)
# create pre-process programs for pyquil that correspond to quara's state
preprocess_programs[state_name] = generate_preprocess_program(qubits, state_name)
print(f"Sample of Quara tester object for state_name: {state_names[0]}")
print(states[state_names[0]])
print("\nPrograms that will be used for state preparation in this tomographical experiment:")
for state_name in state_names:
print(f"{state_name}: {preprocess_programs[state_name]}")
Sample of Quara tester object for state_name: x0
Type:
State
Dim:
2
Vec:
[0.70710678 0.70710678 0. 0. ]
Programs that will be used for state preparation in this tomographical experiment:
x0: H 0
y0: RX(-pi/2) 0
z0:
z1: X 0
Alike quantum state tomography we generate Pauli strings and POVM for measurement.
[10]:
c_sys = generate_composite_system("qubit", 1)
povms = {}
pauli_strings = {}
povm_names = ["x", "y", "z"]
for povm_name in povm_names:
povms[povm_name] = generate_povm_from_name(povm_name, c_sys)
pauli_strings[povm_name] = generate_pauli_strings_from_povm_name(povm_name)
print(f"Sample of Quara tester object for povm_name: {povm_names[0]}")
print(povms[povm_names[0]])
print("\nPauli strings that will be used in this tomographical experiment:")
pauli_strings
Sample of Quara tester object for povm_name: x
Type:
Povm
Dim:
2
Number of outcomes:
2
Vecs:
[[ 0.70710678 0.70710678 0. 0. ]
[ 0.70710678 -0.70710678 0. 0. ]]
Pauli strings that will be used in this tomographical experiment:
[10]:
{'x': ['X', 'I'], 'y': ['Y', 'I'], 'z': ['Z', 'I']}
Now we obtain expectations for each POVM measurement and preprocess states.
[11]:
observables = {}
for state_name in state_names:
observables[state_name] = {}
print(state_name)
for povm_name in povm_names:
print(povm_name)
observable = obtain_expectations_for_qpt(qc, num_shots, preprocess_programs[state_name], target_program, pauli_strings[povm_name])
observables[state_name][povm_name] = observable
x0
x
Raw expectation of X: 0.7016
Raw expectation of I: 1.0
y
Raw expectation of Y: 0.7198
Raw expectation of I: 1.0
z
Raw expectation of Z: 0.0102
Raw expectation of I: 1.0
y0
x
Raw expectation of X: -0.7178
Raw expectation of I: 1.0
y
Raw expectation of Y: 0.7034
Raw expectation of I: 1.0
z
Raw expectation of Z: -0.0034
Raw expectation of I: 1.0
z0
x
Raw expectation of X: 0.0078
Raw expectation of I: 1.0
y
Raw expectation of Y: -0.0018
Raw expectation of I: 1.0
z
Raw expectation of Z: 1.0
Raw expectation of I: 1.0
z1
x
Raw expectation of X: -0.0012
Raw expectation of I: 1.0
y
Raw expectation of Y: 0.0092
Raw expectation of I: 1.0
z
Raw expectation of Z: -1.0
Raw expectation of I: 1.0
Alike state tomography we calculate the empirical distribution using calc_empi_dist_from_observation. Make sure the orders of empi_dists, tester_states, tester_povms match to each other, otherwise the estimation cannot be caculated properly.
[12]:
empi_dists = []
for state_name in state_names:
print(f"empirical distributions for state preparation {state_name}")
for povm_name in povm_names:
empi_dist = calc_empi_dist_from_observables(observables[state_name][povm_name], num_shots, pauli_strings[povm_name], povms[povm_name])
empi_dists.append(empi_dist)
print(f"{povm_name}: {empi_dist}")
tester_states = []
for state_name in state_names:
tester_states.append(states[state_name])
tester_povms = []
for povm_name in povm_names:
tester_povms.append(povms[povm_name])
qpt = StandardQpt(states=tester_states, povms=tester_povms, on_para_eq_constraint=True, schedules="all")
empirical distributions for state preparation x0
x: (10000, array([0.8508, 0.1492]))
y: (10000, array([0.8599, 0.1401]))
z: (10000, array([0.5051, 0.4949]))
empirical distributions for state preparation y0
x: (10000, array([0.1411, 0.8589]))
y: (10000, array([0.8517, 0.1483]))
z: (10000, array([0.4983, 0.5017]))
empirical distributions for state preparation z0
x: (10000, array([0.5039, 0.4961]))
y: (10000, array([0.4991, 0.5009]))
z: (10000, array([1., 0.]))
empirical distributions for state preparation z1
x: (10000, array([0.4994, 0.5006]))
y: (10000, array([0.5046, 0.4954]))
z: (10000, array([1.11022302e-16, 1.00000000e+00]))
The estimate is calculated as follows when linear estimator is chosen.
[13]:
estimator = LinearEstimator()
result = estimator.calc_estimate(qtomography=qpt, empi_dists=empi_dists, is_computation_time_required=True)
estimate = result.estimated_qoperation
print(estimate)
print("is estimate physical? : ", estimate.is_physical())
evals, evecs = np.linalg.eigh(estimate.to_choi_matrix())
print("\nEigenvalues are:", evals)
Type:
Gate
Dim:
2
HS:
[[ 1.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00]
[ 3.30000000e-03 6.98300000e-01 -7.21100000e-01 4.50000000e-03]
[ 3.70000000e-03 7.16100000e-01 6.99700000e-01 -5.50000000e-03]
[-1.11022302e-16 1.02000000e-02 -3.40000000e-03 1.00000000e+00]]
is estimate physical? : False
Eigenvalues are: [-8.42148956e-03 -1.64624349e-03 7.55747516e-03 2.00251026e+00]