Quantum State Tomography¶
import qiskit
from qiskit_experiments.framework import ParallelExperiment
from qiskit_experiments.library import StateTomography
# For simulation
from qiskit.providers.aer import AerSimulator
from qiskit.test.mock import FakeParis
# Noisy simulator backend
backend = AerSimulator.from_backend(FakeParis())
State Tomography Experiment¶
To run a state tomography experiment we initialize the experiment with a
circuit to prepare the state to be measured. We can also pass in an
Operator, or a Statevector to describe the preparation circuit.
# Run experiments
# GHZ State preparation circuit
nq = 2
qc_ghz = qiskit.QuantumCircuit(nq)
qc_ghz.h(0)
qc_ghz.s(0)
for i in range(1, nq):
qc_ghz.cx(0, i)
# QST Experiment
qstexp1 = StateTomography(qc_ghz)
qstdata1 = qstexp1.run(backend, seed_simulation=100).block_for_results()
# Print results
for result in qstdata1.analysis_results():
print(result)
DbAnalysisResultV1
- name: state
- value: DensityMatrix([[ 4.92710916e-01+0.00000000e+00j,
-6.36699653e-03-6.63575796e-03j,
3.50422342e-02-1.88269472e-02j,
4.91481964e-03-4.42405733e-01j],
[-6.36699653e-03+6.63575796e-03j,
3.08516467e-02+0.00000000e+00j,
1.02374555e-02+4.38415725e-03j,
-2.67891920e-04+4.04268830e-03j],
[ 3.50422342e-02+1.88269472e-02j,
1.02374555e-02-4.38415725e-03j,
2.60556731e-02-5.42101086e-20j,
2.14873699e-02-5.83766965e-04j],
[ 4.91481964e-03+4.42405733e-01j,
-2.67891920e-04-4.04268830e-03j,
2.14873699e-02+5.83766965e-04j,
4.50381764e-01+8.67361738e-19j]],
dims=(2, 2))
- extra: <5 items>
- device_components: ['Q0', 'Q1']
- verified: False
DbAnalysisResultV1
- name: state_fidelity
- value: 0.9139520730831032
- device_components: ['Q0', 'Q1']
- verified: False
DbAnalysisResultV1
- name: positive
- value: True
- device_components: ['Q0', 'Q1']
- verified: False
Tomography Results¶
The main result for tomography is the fitted state, which is stored as a
DensityMatrix object:
state_result = qstdata1.analysis_results("state")
print(state_result.value)
DensityMatrix([[ 4.92710916e-01+0.00000000e+00j,
-6.36699653e-03-6.63575796e-03j,
3.50422342e-02-1.88269472e-02j,
4.91481964e-03-4.42405733e-01j],
[-6.36699653e-03+6.63575796e-03j,
3.08516467e-02+0.00000000e+00j,
1.02374555e-02+4.38415725e-03j,
-2.67891920e-04+4.04268830e-03j],
[ 3.50422342e-02+1.88269472e-02j,
1.02374555e-02-4.38415725e-03j,
2.60556731e-02-5.42101086e-20j,
2.14873699e-02-5.83766965e-04j],
[ 4.91481964e-03+4.42405733e-01j,
-2.67891920e-04-4.04268830e-03j,
2.14873699e-02+5.83766965e-04j,
4.50381764e-01+8.67361738e-19j]],
dims=(2, 2))
The state fidelity of the fitted state with the ideal state prepared by
the input circuit is stored in the "state_fidelity" result field.
Note that if the input circuit contained any measurements the ideal
state cannot be automatically generated and this field will be set to
None.
fid_result = qstdata1.analysis_results("state_fidelity")
print("State Fidelity = {:.5f}".format(fid_result.value))
State Fidelity = 0.91395
Additional state metadata¶
Additional data is stored in the tomography under the
"state_metadata" field. This includes - eigvals: the eigenvalues
of the fitted state - trace: the trace of the fitted state -
positive: Whether the eigenvalues are all non-negative -
positive_delta: the deviation from positivity given by 1-norm of
negative eigenvalues.
If trace rescaling was performed this dictionary will also contain a
raw_trace field containing the trace before rescaling. Futhermore,
if the state was rescaled to be positive or trace 1 an additional field
raw_eigvals will contain the state eigenvalues before rescaling was
performed.
state_result.extra
{'trace': 1.0000000000000016,
'eigvals': array([0.91623978, 0.05346034, 0.03029988, 0. ]),
'raw_eigvals': array([ 0.91683962, 0.05406018, 0.03089973, -0.00179953]),
'rescaled_psd': True,
'fitter_metadata': {'fitter': 'linear_inversion',
'fitter_time': 0.00407099723815918}}
To see the effect of rescaling we can perform a “bad” fit with very low counts
# QST Experiment
bad_data = qstexp1.run(backend, shots=10, seed_simulation=100).block_for_results()
bad_state_result = bad_data.analysis_results("state")
# Print result
print(bad_state_result)
# Show extra data
bad_state_result.extra
DbAnalysisResultV1
- name: state
- value: DensityMatrix([[ 0.49835129+0.00000000e+00j, -0.08599238-2.64882188e-02j,
-0.06151686-1.00775529e-02j, -0.17450078-4.48480505e-01j],
[-0.08599238+2.64882188e-02j, 0.02279176+0.00000000e+00j,
0.01667578-1.65429957e-03j, 0.05089237+6.94258072e-02j],
[-0.06151686+1.00775529e-02j, 0.01667578+1.65429957e-03j,
0.01246366-1.08420217e-19j, 0.02800529+5.28834159e-02j],
[-0.17450078+4.48480505e-01j, 0.05089237-6.94258072e-02j,
0.02800529-5.28834159e-02j, 0.4663933 +2.77555756e-17j]],
dims=(2, 2))
- extra: <5 items>
- device_components: ['Q0', 'Q1']
- verified: False
{'trace': 1.0000000000000018,
'eigvals': array([0.98797027, 0.01202973, 0. , 0. ]),
'raw_eigvals': array([ 1.06313803, 0.08719748, -0.04398389, -0.10635162]),
'rescaled_psd': True,
'fitter_metadata': {'fitter': 'linear_inversion',
'fitter_time': 0.002368927001953125}}
Tomography Fitters¶
The default fitters is linear_inversion, which reconstructs the
state using dual basis of the tomography basis. This will typically
result in a non-postive reconstructed state. This state is rescaled to
be postive-semidfinite (PSD) by computing its eigen-decomposition and
rescaling its eigenvalues using the approach from *J Smolin, JM
Gambetta, G Smith, Phys. Rev. Lett. 108, 070502 (2012), open
access.
There are several other fitters are included (See API documentation for
details). For example if cvxpy is installed we can use the
cvxpy_gaussian_lstsq fitter which allows constraining the fit to be
PSD without requiring rescaling.
try:
import cvxpy
# Set analysis option for cvxpy fitter
qstexp1.analysis.set_options(fitter='cvxpy_gaussian_lstsq')
# Re-run experiment
qstdata2 = qstexp1.run(backend, seed_simulation=100).block_for_results()
state_result2 = qstdata2.analysis_results("state")
print(state_result2)
print("\nextra:")
for key, val in state_result2.extra.items():
print(f"- {key}: {val}")
except ModuleNotFoundError:
print("CVXPY is not installed")
DbAnalysisResultV1
- name: state
- value: DensityMatrix([[ 4.73669164e-01+0.00000000e+00j,
1.10952631e-02+1.36892391e-03j,
1.28787403e-03-2.85440631e-03j,
-7.02449084e-05-4.42530533e-01j],
[ 1.10952631e-02-1.36892391e-03j,
2.69467318e-02+0.00000000e+00j,
1.10698978e-02-5.95263273e-05j,
8.75263019e-03-2.16792231e-02j],
[ 1.28787403e-03+2.85440631e-03j,
1.10698978e-02+5.95263273e-05j,
3.47894260e-02+0.00000000e+00j,
2.09348786e-02-6.31235420e-03j],
[-7.02449084e-05+4.42530533e-01j,
8.75263019e-03+2.16792231e-02j,
2.09348786e-02+6.31235420e-03j,
4.64594679e-01+0.00000000e+00j]],
dims=(2, 2))
- extra: <5 items>
- device_components: ['Q0', 'Q1']
- verified: False
extra:
- trace: 1
- eigvals: [0.91267941 0.05246603 0.02154755 0.01330701]
- raw_eigvals: [0.91267941 0.05246603 0.02154755 0.01330701]
- rescaled_psd: False
- fitter_metadata: {'cvxpy_solver': 'SCS', 'cvxpy_status': 'optimal', 'fitter': 'cvxpy_gaussian_lstsq', 'fitter_time': 0.02846813201904297}
Parallel Tomography Experiment¶
We can also use the qiskit_experiments.ParallelExperiment class to
run subsystem tomography on multiple qubits in parallel.
For example if we want to perform 1-qubit QST on several qubits at once:
from math import pi
num_qubits = 5
gates = [qiskit.circuit.library.RXGate(i * pi / (num_qubits - 1))
for i in range(num_qubits)]
subexps = [
StateTomography(gate, qubits=[i])
for i, gate in enumerate(gates)
]
parexp = ParallelExperiment(subexps)
pardata = parexp.run(backend, seed_simulation=100).block_for_results()
for result in pardata.analysis_results():
print(result)
View component experiment analysis results
for i, expdata in enumerate(pardata.child_data()):
state_result_i = expdata.analysis_results("state")
fid_result_i = expdata.analysis_results("state_fidelity")
print(f'\nPARALLEL EXP {i}')
print("State Fidelity: {:.5f}".format(fid_result_i.value))
print("State: {}".format(state_result_i.value))
PARALLEL EXP 0
State Fidelity: 0.99121
State: DensityMatrix([[ 0.99121094+0.j , -0.00585938+0.01171875j],
[-0.00585938-0.01171875j, 0.00878906+0.j ]],
dims=(2,))
PARALLEL EXP 1
State Fidelity: 0.98752
State: DensityMatrix([[0.85351563+0.j , 0.00976562+0.3359375j],
[0.00976562-0.3359375j, 0.14648438+0.j ]],
dims=(2,))
PARALLEL EXP 2
State Fidelity: 0.96973
State: DensityMatrix([[0.50585938+0.j , 0.02832031+0.46972656j],
[0.02832031-0.46972656j, 0.49414063+0.j ]],
dims=(2,))
PARALLEL EXP 3
State Fidelity: 0.98199
State: DensityMatrix([[ 0.16503906+0.j , -0.00390625+0.34667969j],
[-0.00390625-0.34667969j, 0.83496094+0.j ]],
dims=(2,))
PARALLEL EXP 4
State Fidelity: 0.97168
State: DensityMatrix([[0.02832031+0.j , 0.00976562-0.00390625j],
[0.00976562+0.00390625j, 0.97167969+0.j ]],
dims=(2,))
import qiskit.tools.jupyter
%qiskit_copyright
This code is a part of Qiskit
© Copyright IBM 2017, 2022.
This code is licensed under the Apache License, Version 2.0. You may
obtain a copy of this license in the LICENSE.txt file in the root directory
of this source tree or at http://www.apache.org/licenses/LICENSE-2.0.
Any modifications or derivative works of this code must retain this
copyright notice, and modified files need to carry a notice indicating
that they have been altered from the originals.