Note
This is the documentation for the current state of the development branch of Qiskit Experiments. The documentation or APIs here can change prior to being released.
Quantum State Tomography¶
Quantum tomography is an experimental procedure to reconstruct a description of part of a quantum system from the measurement outcomes of a specific set of experiments. In particular, quantum state tomography reconstructs the density matrix of a quantum state by preparing the state many times and measuring them in a tomographically complete basis of measurement operators.
Note
This tutorial requires the qiskit-aer and qiskit-ibm-runtime
packages to run simulations. You can install them with python -m pip
install qiskit-aer qiskit-ibm-runtime
.
We first initialize a simulator to run the experiments on.
from qiskit_aer import AerSimulator
from qiskit_ibm_runtime.fake_provider import FakePerth
backend = AerSimulator.from_backend(FakePerth())
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.
import qiskit
from qiskit_experiments.framework import ParallelExperiment
from qiskit_experiments.library import StateTomography
# 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)
AnalysisResult
- name: state
- value: DensityMatrix([[ 0.48307292+0.j , -0.0069987 +0.00683594j,
0.00325521+0.00130208j, -0.00488281-0.43408203j],
[-0.0069987 -0.00683594j, 0.01334635+0.j ,
-0.01757812+0.00341797j, -0.01041667+0.00227865j],
[ 0.00325521-0.00130208j, -0.01757812-0.00341797j,
0.0296224 +0.j , -0.00992839-0.00878906j],
[-0.00488281+0.43408203j, -0.01041667-0.00227865j,
-0.00992839+0.00878906j, 0.47395833+0.j ]],
dims=(2, 2))
- quality: unknown
- extra: <9 items>
- device_components: ['Q0', 'Q1']
- verified: False
AnalysisResult
- name: state_fidelity
- value: 0.9125976562499989
- quality: unknown
- extra: <9 items>
- device_components: ['Q0', 'Q1']
- verified: False
AnalysisResult
- name: positive
- value: True
- quality: unknown
- extra: <9 items>
- 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([[ 0.48307292+0.j , -0.0069987 +0.00683594j,
0.00325521+0.00130208j, -0.00488281-0.43408203j],
[-0.0069987 -0.00683594j, 0.01334635+0.j ,
-0.01757812+0.00341797j, -0.01041667+0.00227865j],
[ 0.00325521-0.00130208j, -0.01757812-0.00341797j,
0.0296224 +0.j , -0.00992839-0.00878906j],
[-0.00488281+0.43408203j, -0.01041667-0.00227865j,
-0.00992839+0.00878906j, 0.47395833+0.j ]],
dims=(2, 2))
We can also visualize the density matrix:
from qiskit.visualization import plot_state_city
plot_state_city(qstdata1.analysis_results("state").value, title='Density Matrix')
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.91260
Additional state metadata¶
Additional data is stored in the tomography under the
"state_metadata"
field. This includes
eigvals
: the eigenvalues of the fitted statetrace
: the trace of the fitted statepositive
: Whether the eigenvalues are all non-negativepositive_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.0000000000000013,
'eigvals': array([9.13011298e-01, 4.81233589e-02, 3.85285588e-02, 3.36783977e-04]),
'raw_eigvals': array([9.13011298e-01, 4.81233589e-02, 3.85285588e-02, 3.36783977e-04]),
'rescaled_psd': False,
'fitter_metadata': {'fitter': 'linear_inversion',
'fitter_time': 0.004160642623901367},
'conditional_probability': 1.0,
'positive': True,
'experiment': 'StateTomography',
'run_time': None}
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
AnalysisResult
- name: state
- value: DensityMatrix([[ 0.4385565 +0.00000000e+00j, 0.04403168-6.63251038e-02j,
-0.07404532+1.00399458e-01j, 0.02235449-4.14912645e-01j],
[ 0.04403168+6.63251038e-02j, 0.06766099+1.73472348e-18j,
-0.05555338+7.04805084e-03j, 0.10247157-8.57795571e-03j],
[-0.07404532-1.00399458e-01j, -0.05555338-7.04805084e-03j,
0.05712557+0.00000000e+00j, -0.11740077+4.08010434e-02j],
[ 0.02235449+4.14912645e-01j, 0.10247157+8.57795571e-03j,
-0.11740077-4.08010434e-02j, 0.43665694+0.00000000e+00j]],
dims=(2, 2))
- quality: unknown
- extra: <9 items>
- device_components: ['Q0', 'Q1']
- verified: False
{'trace': 1.000000000000001,
'eigvals': array([0.91125604, 0.08874396, 0. , 0. ]),
'raw_eigvals': array([ 1.01321646e+00, 1.90704372e-01, -2.82380580e-04, -2.03638449e-01]),
'rescaled_psd': True,
'fitter_metadata': {'fitter': 'linear_inversion',
'fitter_time': 0.003038644790649414},
'conditional_probability': 1.0,
'positive': True,
'experiment': 'StateTomography',
'run_time': None}
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-positive reconstructed state. This state is rescaled to
be positive-semidefinite (PSD) by computing its eigen-decomposition and
rescaling its eigenvalues using the approach from Ref. [1].
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")
AnalysisResult
- name: state
- value: DensityMatrix([[ 0.47815411+0.j , -0.00664353+0.00996264j,
0.00425638+0.0070674j , 0.0130099 -0.43995094j],
[-0.00664353-0.00996264j, 0.0342797 +0.j ,
0.00397864+0.00143443j, 0.00120946-0.00607383j],
[ 0.00425638-0.0070674j , 0.00397864-0.00143443j,
0.021964 +0.j , 0.00752466-0.00672796j],
[ 0.0130099 +0.43995094j, 0.00120946+0.00607383j,
0.00752466+0.00672796j, 0.46560219+0.j ]],
dims=(2, 2))
- quality: unknown
- extra: <9 items>
- device_components: ['Q0', 'Q1']
- verified: False
extra:
- trace: 1.0000000014581047
- eigvals: [0.9121827 0.04817207 0.02622239 0.01342285]
- raw_eigvals: [0.9121827 0.04817207 0.02622239 0.01342285]
- rescaled_psd: False
- fitter_metadata: {'fitter': 'cvxpy_gaussian_lstsq', 'cvxpy_solver': 'SCS', 'cvxpy_status': ['optimal'], 'psd_constraint': True, 'trace_preserving': True, 'fitter_time': 0.021170616149902344}
- conditional_probability: 1.0
- positive: True
- experiment: StateTomography
- run_time: None
Parallel Tomography Experiment¶
We can also use the 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, physical_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)
AnalysisResult
- name: state
- value: DensityMatrix([[ 0.97460938+0.j , -0.02832031+0.00390625j],
[-0.02832031-0.00390625j, 0.02539063+0.j ]],
dims=(2,))
- quality: unknown
- extra: <9 items>
- device_components: ['Q0']
- verified: False
AnalysisResult
- name: state_fidelity
- value: 0.9746093750000002
- quality: unknown
- extra: <9 items>
- device_components: ['Q0']
- verified: False
AnalysisResult
- name: positive
- value: True
- quality: unknown
- extra: <9 items>
- device_components: ['Q0']
- verified: False
AnalysisResult
- name: state
- value: DensityMatrix([[ 0.81445313+0.j , -0.02832031+0.34375j],
[-0.02832031-0.34375j, 0.18554688+0.j ]],
dims=(2,))
- quality: unknown
- extra: <9 items>
- device_components: ['Q1']
- verified: False
AnalysisResult
- name: state_fidelity
- value: 0.965419893085677
- quality: unknown
- extra: <9 items>
- device_components: ['Q1']
- verified: False
AnalysisResult
- name: positive
- value: True
- quality: unknown
- extra: <9 items>
- device_components: ['Q1']
- verified: False
AnalysisResult
- name: state
- value: DensityMatrix([[0.49511719+0.j , 0.00585938+0.46289062j],
[0.00585938-0.46289062j, 0.50488281+0.j ]],
dims=(2,))
- quality: unknown
- extra: <9 items>
- device_components: ['Q2']
- verified: False
AnalysisResult
- name: state_fidelity
- value: 0.9628906250000007
- quality: unknown
- extra: <9 items>
- device_components: ['Q2']
- verified: False
AnalysisResult
- name: positive
- value: True
- quality: unknown
- extra: <9 items>
- device_components: ['Q2']
- verified: False
AnalysisResult
- name: state
- value: DensityMatrix([[0.16894531+0.j , 0.01660156+0.34570312j],
[0.01660156-0.34570312j, 0.83105469+0.j ]],
dims=(2,))
- quality: unknown
- extra: <9 items>
- device_components: ['Q3']
- verified: False
AnalysisResult
- name: state_fidelity
- value: 0.9785400384397234
- quality: unknown
- extra: <9 items>
- device_components: ['Q3']
- verified: False
AnalysisResult
- name: positive
- value: True
- quality: unknown
- extra: <9 items>
- device_components: ['Q3']
- verified: False
AnalysisResult
- name: state
- value: DensityMatrix([[ 0.02441406+0.j , -0.00195312+0.01855469j],
[-0.00195312-0.01855469j, 0.97558594+0.j ]],
dims=(2,))
- quality: unknown
- extra: <9 items>
- device_components: ['Q4']
- verified: False
AnalysisResult
- name: state_fidelity
- value: 0.9755859375000002
- quality: unknown
- extra: <9 items>
- device_components: ['Q4']
- verified: False
AnalysisResult
- name: positive
- value: True
- quality: unknown
- extra: <9 items>
- device_components: ['Q4']
- verified: False
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))
References¶
See also¶
API documentation:
StateTomography