Allow QuantumCircuit to contain non-CPTP channels

[kevinsung]

What should we add?

It's useful to be able store non-CPTP channels in a QuantumCircuit. The density matrix can still be evolved, it just won't be a normalized state. Apparently, the error is thrown upon conversion to an Instruction, so we should just make it so the error is not thrown there.

Code:

import numpy as np
from qiskit import QuantumCircuit
from qiskit.quantum_info import PTM

circuit = QuantumCircuit(2)
ptm = PTM(0.5 * np.eye(16))

circuit.append(ptm, [0, 1])

Result:

[jakelishman]

What situations did you have in mind where you would want to be able to store a non-CPTP channel in a circuit? Going a shade further: how much of the CPTP condition is necessary to relax for your use-cases - could we still limit operations to being CP?

We don't have clear documentation the limitations of Instruction, but the line in QuantumChannel.to_instruction that throws the error you're seeing is clearly deliberate from when it was added in #2173, which makes me a bit nervous that there will be downstream consumers of QuantumCircuit that assume that at any program execution point, there can be a valid (normalised) quantum state associated with it.

I imagine that the only place that really allows/interprets quantum channels is Aer already, so the biggest concerns are that we wouldn't be breaking any assumptions that Aer makes. I can already see a couple of assumptions baked into QuantumChannel.to_instruction that are based on the CPTP nature - for example, if the channel comprises a single Kraus operator, then it outputs a UnitaryGate, which ofc isn't valid if the channel isn't TP. If we're to continue outputting Kraus operators, I think we at least must enforce the CP requirement (and also, just to be clear: in your example, the PTM wouldn't actually get stored as a PTM, it'd get stored as a list of Kraus operators anyway).

[kevinsung]

In my use case, the PTM is obtained from process tomography, and due to statistical and experimental error, it may fail to be TP. I ran some numerical experiments and it seemed to always be CP, but I'm not sure this will always be the case.

Here is the code for my experiment:

from qiskit_aer import AerSimulator, noise
from qiskit.circuit import QuantumCircuit, QuantumRegister
from qiskit.circuit.library import CZGate
from qiskit_experiments.library import ProcessTomography

basis_gates = ["cz", "id", "rz", "sx", "x"]

depolarizing_error = 0.4
amplitude_damping_error = 0.3
readout_error = 0.2
shots=100

noise_model = noise.NoiseModel(basis_gates=basis_gates)
noise_model.add_all_qubit_quantum_error(
    noise.depolarizing_error(depolarizing_error, 2), "cz"
)
noise_model.add_all_qubit_quantum_error(
    noise.amplitude_damping_error(amplitude_damping_error).tensor(
        noise.amplitude_damping_error(amplitude_damping_error)
    ),
    "cz",
)
noise_model.add_all_qubit_readout_error(
    noise.ReadoutError(
        [
            [1 - readout_error, readout_error],
            [readout_error, 1 - readout_error],
        ]
    )
)

n_qubits = 10
coupling_map = [[i, i + 1] for i in range(n_qubits - 1)]
backend = AerSimulator(
    n_qubits=n_qubits,
    basis_gates=basis_gates,
    coupling_map=coupling_map,
    noise_model=noise_model,
    seed_simulator=1234,
)

register = QuantumRegister(2, name="q")
target_circuit = QuantumCircuit(register)
a, b = register
target_circuit.append(CZGate(), [a, b])
physical_qubits = [0, 1]

experiment = ProcessTomography(
    target_circuit, backend=backend, physical_qubits=physical_qubits
)
data = experiment.run(backend, shots=shots)
data.block_for_results()
choi = data.analysis_results("state").value

print(f"CP: {choi.is_cp()}")
print(f"TP: {choi.is_tp()}")

CP: True
TP: False

[jakelishman]

Coming back to this after being on holiday: I think if we can continue with the CP requirement, it's probably ok, since QuantumChannel.to_instruction uses Kraus operators. I do wonder a little if a better root solution, though, would be to have some method on QuantumChannel that scales the channel to make it TP (or some other "nearest" TP channel for some definition of nearest) might be more appropriate - letting QuantumCircuit store non-TP channels with the goal of making things silently easier might just be shifting numerical-stability problems further down the line, like a simulation now returning a non-normalised state unexpectedly.

[kevinsung]

This paper discusses CPTP projection of a map: https://journals.aps.org/pra/abstract/10.1103/PhysRevA.98.062336. I believe @ebennewitz and/or @neversakura coded it up for a separate project.

[ebennewitz]

Yea, I coded up the CPTP projection (as outlined in section III: Composite Projection and in Algorithm 1) for another project. I'd be happy to share my code! The authors of the above paper also put out a matlab implementation which was a useful reference: https://github.com/geoknee/CPTPprojection/tree/master

``