Process Tomography

Process tomography allows for the characterisation of a quantum process \(\mathcal{E}\), where \(\mathcal{E}(\rho) = \rho^\prime\). This is very useful when assessing the exact operation of different quantum gates and can be used for the assessment of gate fidelities.

This notebook demonstates the utilisation of the different process tomography algorithms included within Lightworks. In each the choi matrix is found and used as a characterisation of the quantum process which occured.

import matplotlib.pyplot as plt
import numpy as np

import lightworks as lw
from lightworks import emulator, qubit
from lightworks.tomography import (
    LIProcessTomography,
    MLEProcessTomography,
    choi_from_unitary,
)

Before starting, a general function is defined to quickly perform the plotting of choi matrices. This takes a complex matrix and plots the real and imaginary parts separately.

def plot_choi_matrix(choi: np.ndarray) -> tuple:
    """
    General function for plotting a choi matrix. It will split up and plot
    the real and imaginary components using a common colorbar between them.
    """
    # Create required plotting data
    xx, yy = np.meshgrid(range(choi.shape[0]), range(choi.shape[1]))
    x, y = xx.ravel(), yy.ravel()
    width = depth = 0.7
    h1 = np.real(choi).flatten()
    h2 = np.imag(choi).flatten()
    base = np.zeros_like(h1)

    # Create figure
    fig = plt.figure(figsize=(11, 5))
    ax1 = fig.add_subplot(121, projection="3d")
    ax2 = fig.add_subplot(122, projection="3d")

    # Plot bar charts
    ax1.bar3d(x, y, base, width, depth, h1, shade=True)
    ax1.set_title("Re(C)")
    ax2.bar3d(x, y, base, width, depth, h2, shade=True)
    ax2.set_title("Im(C)")
    for ax in [ax1, ax2]:
        ax.set_zlim(-1, 1)

    return fig, (ax1, ax2)

Setup

To perform process tomography, the circuit to analyse needs to be provided. This should be a Lightworks circuit object, where the number of available modes (excluding heralded) is equal to the number of qubits.

In this case, we’ll examine the post-selected CNOT gate, which acts across two qubits. The post-selection requires measurement of 1 photon across each pair of qubits modes, this is also configured below.

n_qubits = 2

cnot = qubit.CNOT()
cnot.display()

post_select = lw.PostSelection()
post_select.add((0, 1), 1)
post_select.add((2, 3), 1)

There is then a choice of algorithms that can be used for the calculation of the choi matrix for the process:

1. Linear inversion - this is a simpler algorithm with less measurements, but can produce non-physical process matrices (particularily in noisy scenarios).

2. Maximum likelihood estimation - this uses a gradient-descent based optimisation, along with a set of projections to enforce conditions for a physical process matrix.

Linear Inversion

This method is achieved with the LIProcessTomography object. To initialise this, the number of qubits & the circuit is provided. The required experiments are then generated with get_experiments, which is provided as a list of python ProcessTomographyExperiment objects.

li_tomo = LIProcessTomography(n_qubits, cnot)
experiments = li_tomo.get_experiments()

The experiments are then run on the target backend. For process tomography, each experiments specifies a circuit and input state to use, these should both be provided to the Sampler (or whichever task you choose).

backend = emulator.Backend("slos")

# Generate results and return
results = []
for exp in experiments:
    sampler = lw.Sampler(
        exp.circuit,
        exp.input_state,
        50000,
        post_selection=post_select,
        random_seed=98,
    )
    results.append(backend.run(sampler))

Finally, the choi matrix is calculated by providing this result set to the process method.

choi = li_tomo.process(results)

Once calculated, we plot this to view the structure of the density matrix.

plot_choi_matrix(choi)
plt.show()

The choi_from_untiary function can then be used to find the expected unitary for a CNOT gate and the fidelity calculated with the relevant fidelity method. As expected this is near to 100%, with the variance resulting mostly from sample noise.

choi_exp = choi_from_unitary(
    [[1, 0, 0, 0], [0, 1, 0, 0], [0, 0, 0, 1], [0, 0, 1, 0]]
)
print(f"Fidelity = {li_tomo.fidelity(choi_exp) * 100:.3f} %")

Fidelity = 99.997 %

Projection to a physical process

As with state tomography, it is also possible for the produced choi matrix to be unphysical. For the choi matrix the following is required:

1. It is completely positive, this is equivalent to having no negative eigenvalues.

2. It is trace preserving, or equivalently \(\text{Tr}_{\mathcal{H}_{out}}(C)=\mathbb{I}\).

There are two different ways to deal with this in Lightworks, both of which are detailed below.

Maximum Likelihood Estimation

The first option is the pgdB algorithm from [KBLG18], which is used for the projection and fitting of the choi matrix through maximum likelihood estimation. This is a series of projections which minimises the variance between the fitted and measured processes. There is a downside however, in that it requires \(6^n\) input states rather than the \(4^n\) required for normal process tomography.

This method is implemented using the MLEProcessTomography object and utilised in the same way as above. However, instead of looping through the experiments, we use all_circuits and all_inputs to get a list of the circuits and input states directly, which are compiled within a batch task.

# Create tomography object and get experiments
mle_tomo = MLEProcessTomography(n_qubits, cnot)
experiments = mle_tomo.get_experiments()

backend = emulator.Backend("slos")

# Generate results and return
batch = lw.Batch(
    lw.Sampler,
    [experiments.all_circuits, experiments.all_inputs, [10000]],
    {"post_selection": [post_select], "random_seed": [98]},
)
results = backend.run(batch)

# Find choi matrix
choi = mle_tomo.process(results)

Again, the fidelity can be calculated using the expected matrix and the results plotted to confirm the structure remains consistent.

print(f"Fidelity = {mle_tomo.fidelity(choi_exp) * 100:.4f} %")

plot_choi_matrix(choi)
plt.show()

Fidelity = 100.0035 %

Imperfect Experiment

As with state tomography, it is also possible to introduce imperfections into the experimental simulation, for a better understanding of how these errors affect fidelity.

For this, we’ll use a class to create a parametrised experiment with optional single photon source noise and then use this for running the experiments.

class SamplerExperiment:
    """
    Runs experiment using the emulator Sampler for a given n_qubits with an
    optional imperfect single photon source.
    """

    def __init__(
        self, n_qubits: int, source: emulator.Source | None = None
    ) -> None:
        self.n_qubits = n_qubits
        self.source = source
        self.n_samples = 50000
        self.random_seed = 123
        self.backend = emulator.Backend("slos")

    def run(self, experiments) -> list:
        """
        Generalised version of experiment function above, designed for any
        number of qubits. It is assumes the provided circuits contain dual-rail
        encoded qubits across pairs of adjacent modes.
        """
        # Post-select on 1 photon across each pair of qubit modes
        post_select = lw.PostSelection()
        for i in range(self.n_qubits):
            post_select.add((2 * i, 2 * i + 1), 1)

        # Use a batch to run all tasks at once
        batch = lw.Batch(
            lw.Sampler,
            task_args=[
                experiments.all_circuits,
                experiments.all_inputs,
                [self.n_samples],
            ],
            task_kwargs={
                "source": [self.source],
                "post_selection": [post_select],
                "random_seed": [self.random_seed],
            },
        )
        return self.backend.run(batch)

The single photon source parameters are then configured and assigned below, setting the source indistinguishability to 95% and purity to 99%.

source = emulator.Source(indistinguishability=0.95, purity=0.99)

imperfect_exp = SamplerExperiment(2, source=source)

Once the tomography is performed, as expected we see a drop in fidelity and the presence of additional erroneous elements in the choi matrix.

mle_tomo_2 = MLEProcessTomography(n_qubits, cnot)
experiments = mle_tomo_2.get_experiments()

results = imperfect_exp.run(experiments)
choi = mle_tomo_2.process(results)

print(f"Fidelity = {mle_tomo_2.fidelity(choi_exp) * 100:.2f} %")

plot_choi_matrix(choi)
plt.show()

Fidelity = 85.25 %

Direct choi projection

The alternative is to use linear inversion tomography and then perform a projection step at the end. In Lightworks, the CPTP projection subroutine from MLE tomography (Algorithm 1 in [KBLG18]) is used.

To examine the effect of this, the experiments are run on a noisy system and then the choi matrix calculated without projection.

li_tomo = LIProcessTomography(n_qubits, cnot)
experiments = li_tomo.get_experiments()

results = imperfect_exp.run(experiments)

choi = li_tomo.process(results)

Then, to check if the calculated matrix is physical, the eigenvalues can be found directly. In this case it can be seen how a significant number of these are negative.

np.sort(np.real(np.linalg.eigvals(choi)))

array([-6.03772803e-02, -3.91692590e-02, -2.01208228e-02, -8.08525147e-03,
       -4.99687703e-03, -1.16218470e-03, -2.25089136e-04,  2.00668272e-04,
        5.33464843e-04,  2.70653003e-03,  1.37126584e-02,  2.60911814e-02,
        4.34609807e-02,  1.18602889e-01,  2.45125605e-01,  3.68370279e+00])

Alternatively, the check eigenvalues method can be used to perform a state tomography on each input and find the eigenvalues of these density matrices. This is not directly equivalent to above, but should allow for an indication of potential issues and also allow for the identification of any problematic inputs.

li_tomo.check_eigenvalues(results)

{'Z+,Z+': array([-8.27398291e-04, -6.66382273e-07,  8.20000000e-04,  1.00000806e+00]),
 'Z+,Z-': array([-8.27398291e-04, -6.66382273e-07,  8.20000000e-04,  1.00000806e+00]),
 'Z+,X+': array([-8.27398291e-04, -6.66382273e-07,  8.20000000e-04,  1.00000806e+00]),
 'Z+,Y+': array([-8.27398291e-04, -6.66382273e-07,  8.20000000e-04,  1.00000806e+00]),
 'Z-,Z+': array([-3.48642964e-05, -2.39445246e-06,  8.64395747e-02,  9.13597684e-01]),
 'Z-,Z-': array([-3.48642964e-05, -2.39445246e-06,  8.64395747e-02,  9.13597684e-01]),
 'Z-,X+': array([-8.27398291e-04, -6.66382273e-07,  8.20000000e-04,  1.00000806e+00]),
 'Z-,Y+': array([-2.89682582e-05, -1.25454177e-06,  8.64296910e-02,  9.13600532e-01]),
 'X+,Z+': array([-0.00162737,  0.00108129,  0.07161104,  0.92893504]),
 'X+,Z-': array([-4.75020216e-04,  2.98610722e-04,  7.12395729e-02,  9.28936837e-01]),
 'X+,X+': array([-1.21858635e-05, -6.76875063e-07,  8.68024090e-02,  9.13210454e-01]),
 'X+,Y+': array([-8.28289905e-04,  2.78897208e-04,  7.16137855e-02,  9.28935607e-01]),
 'Y+,Z+': array([-1.03924709e-03,  4.95328673e-04,  7.16072486e-02,  9.28936670e-01]),
 'Y+,Z-': array([-1.07988172e-03,  8.92478638e-04,  7.12524403e-02,  9.28934963e-01]),
 'Y+,X+': array([-1.21858635e-05, -6.76875063e-07,  8.68024090e-02,  9.13210454e-01]),
 'Y+,Y+': array([-1.38458136e-03,  8.44230721e-04,  7.16050764e-02,  9.28935274e-01])}

The choi matrix can then be projected by adding project_to_physical = True to the original process call.

choi = li_tomo.process(results, project_to_physical=True)

Again, the eigenvalues can be examined, in this case it can be seen that while a few are technically negative, they are effectively zero and therefore not an issue.

np.sort(np.real(np.linalg.eigvals(choi)))

array([-2.27992457e-17, -3.13751246e-18, -2.16073520e-18, -8.36340788e-19,
       -4.64465885e-20,  1.32871001e-21,  4.00708695e-19,  4.95347843e-18,
        5.03172483e-18,  5.38552379e-18,  4.53125866e-16,  8.96280143e-03,
        2.89721701e-02,  9.04046604e-02,  2.20360468e-01,  3.65379082e+00])

This new matrix can then be plotted and the fidelity checked to confirm it seems reasonable.

print(f"Fidelity = {li_tomo.fidelity(choi_exp) * 100:.2f} %")

plot_choi_matrix(choi)
plt.show()

Fidelity = 91.22 %