I was trying to obtain the quasi-probability decomposition of the CNOT gate by using the information in this paper.

The authors give us the example for the CZ gate (Figure 2, i.e. the one below).


The problem is that by using that decomposition I get $\begin{bmatrix} 1-i & 0 & 0 & 0 \\ 0 & 1+i& 0 & 0 \\ 0 & 0 & 1+i & 0 \\ 0 & 0 & 0 & -1-i \end{bmatrix}$ instead of the CZ gate that is $\begin{bmatrix} 1 & 0 & 0 & 0 \\ 0 & 1& 0 & 0 \\ 0 & 0 & 1 & 0 \\ 0 & 0 & 0 & -1 \end{bmatrix}$.

I'm not even sure on how to convert this information to obtain the CNOT gate, but first I need to get the CZ gate with that calculation. Is there anybody that can help me obtaining the correct result?

I used a small python script to obtain that matrix:

import numpy as np

i = np.eye(2)

z = np.array([[1,0],[0,-1]])

t1  = 1/2*np.kron(1/np.sqrt(2)*(i+1j*z), 1/np.sqrt(2)*(i+1j*z))

t2  = 1/2*np.kron(1/np.sqrt(2)*(i-1j*z), 1/np.sqrt(2)*(i-1j*z))

t3 = 0
for a1 in [-1,1]:
   for a2 in [-1,1]:
       t3 += -1/2 *a1*a2*( np.kron((i+a1*z)/2, np.cos((a2+1)*np.pi/4)*i+1j*np.sin((a2+1)*np.pi/4)*z) + np.kron(np.cos((a1+1)*np.pi/4)*i+1j*np.sin((a1+1)*np.pi/4)*z, (i+a2*z)/2) )


EDIT: I'm starting to think that the "=" in that image is not totally true. If I use that decomposition the expectation values of all observables are the same of the original circuit ones, but the final state changes. Since they focus on the expectation value, probably they mean that they are equal in that way.

  • 1
    $\begingroup$ I've done the computations by hand and your script gives the right result. I encourage you to contact the authors to ask them the question directly. Either you're right and they can correct the article, or they will explain what you've misunderstood. Once that done, don't hesitate to post the answer here accordingly! I've quickly read the paper and they justify that this decomposition gives the $CZ$ gate up to a global phase. However, that would be the case if the first coefficient was $1+\mathrm{i}$ instead of $1-\mathrm{i}$ (and if there was a $\frac{1}{\sqrt{2}}$ factor on the matrix. $\endgroup$
    – Tristan Nemoz
    Sep 28, 2022 at 7:46
  • 1
    $\begingroup$ I could obviously be wrong though, but I think the authors would be able to clarify this for you! More often than not, paper authors are happy to answer questions on their work (except when it means they have to correct it but hey, test your luck :D) $\endgroup$
    – Tristan Nemoz
    Sep 28, 2022 at 7:47
  • 1
    $\begingroup$ @TristanNemoz Thank you for your suggestion. I have sent the email! $\endgroup$
    – stopper
    Sep 28, 2022 at 9:02

3 Answers 3


The channel representation of a quantum gate $U$ essentially corresponds to the action of a gate on the density matrix $\rho$, as in $\mathcal{S}(U) = U \rho U^{\dagger}$, instead of its action on the pure state wave function $|\psi \rangle$, as in $U |\psi \rangle$. The conventional description of quantum circuits, particularly in the absence of mid-circuit measurements, is the latter, as we can regard the state of the quantum register as a pure state $|\psi\rangle$ acted upon by unitary operators $\{U_i\}$ (i.e., quantum gates), which maintain the purity of the state. The channel representation is the appropriate description of a circuit where some of its operations are non-unitary (e.g., mid-circuit measurements), in which case the state of the quantum register can only be fully described in terms of a density operator $\rho$ acted upon by superoperators $\mathcal{S}(U)$.

On page 3 of the paper, Mitarai and Fujii explain why they adopted this superoperator decomposition instead of the standard operator decomposition.

It is noteworthy that as we perform decompositions of a superoperator rather than an operator such as $U$ itself, the method becomes friendly for a realistic quantum device. A direct decomposition of $U$ into some simple operators {$V_i$}, i.e. $U = \sum_i c_i V_i$, can also be utilized for the same task; however, as expectation values are calculated as $\langle 0|U^{\dagger} O U|0\rangle$ where $|0\rangle$ is an initial state, this approach requires us to evaluate $\sum_{i,j} c_i c^{*}_j \langle 0|V_{j}^{\dagger} O V_i|0\rangle$, which are rather hard for the NISQ devices. This fact demonstrates the advantage of using the above formalism. The tensor network representation of the superoperator formalism allows us to graphically understand the decompositions.

To rephrase the quote, it is advantageous to use the superoperator formalism in this context because it allows us to express the expectation value that we wish to compute as a linear combination of expectation values with the symmetric form $\langle 0|P^{\dagger} O P|0\rangle$, while the operator formalism would yield some expectation values $\langle 0|Q^{\dagger} O P|0\rangle$ with $P \neq Q$. The former can be straightforwardly computed, while the latter typically require more sophisticated methods that roughly double the circuit depth (cf., e.g., Fig. 9 in R. Somma et al., Phys. Rev. A 65, 042323 (2002))

The point about this paper by Mitarai and Fujii is not to express a two-qubit gate as a linear combination of tensor products of single-qubit gates. In fact, for the sort of two-qubit gate they consider, $e^{i \theta A_1 \otimes A_2}$, with $A_1^2 = \mathbb{1}$ and $A_2^2 = \mathbb{1}$, there is one such trivial decomposition: $e^{i \theta A_1 \otimes A_2} = \cos \theta \; \mathbb{1} \otimes \mathbb{1} + i \sin \theta \; A_1 \otimes A_2$. Hence, as you rightly pointed out in your edit, the equalities in Fig. 1 of the paper do not refer to the $4 \times 4$ matrix representations of the two-qubit operations but rather to the expectation values computed with such two-qubit operations. Or, equivalently, in the words of one of the authors, the equality holds when we consider the operations as quantum channels.

Regarding your ultimate goal of obtaining a quasiprobabilistic decomposition of the CNOT gate, it is straightforward to adapt the expression provided by Mitarai and Fujii for the CZ gate using the identity $CNOT = (\mathbb{1} \otimes H) CZ (\mathbb{1} \otimes H)$, where $H$ is the Hadamard gate. The analogue of Fig. 2 for the CNOT gate is as follows.

enter image description here

Basically, all the Z operations are replaced by X operations in the target qubit only. You can easily check that this linear combination yields exactly the same expectation values as the original CNOT gate. However, its $4 \times 4$ matrix representation is

$\begin{pmatrix} 1 & -i & 0 & 0 \\ -i & 1 & 0 & 0 \\ 0 & 0 & 0 & 1+i \\ 0 & 0 & 1+i & 0 \end{pmatrix}$,

which is clearly not the CNOT gate.

  • $\begingroup$ Thank you for your clarification! I've also came up with the same decomposition for the CNOT gate. I would like to ask you another question if possible: in the decomposition we use the projective measurement (I+Z)/2, and if I use exactly this operation it seems that everything works well. Instead, if I try to really implement that operation by projecting the state into |0> I do not obtain the correct result for certain expectation values. The difference between the 2 methods is that in the first one I have a non-normalized state. How do I implement such a thing in a real quantum device? $\endgroup$
    – stopper
    Sep 30, 2022 at 10:20
  • 1
    $\begingroup$ You can fix the second approach if you multiply the projective measurement by the probability of such measurement outcome being obtained, as noted in Eq. (8) of the paper. On an actual quantum device, whenever (1 + $\alpha$ A)/2 appears in the middle of the circuit, we measure the operator A. This yields one of two outcomes, +1 or -1, which define the $\alpha$ associated with this particular circuit execution. For one such execution, no multiplication by the probability is made; this is only done at the end, when we have the statistics from the multiple circuit executions. $\endgroup$
    – bm442
    Sep 30, 2022 at 11:47
  • $\begingroup$ Thank you again for your answer, now it is finally clear. I have just one more question to proceed with my work: in this paper( arxiv.org/pdf/1712.09271.pdf ) and in other ones they introduce the C factor to evaluate the sampling overhead (C^(2n)) of the decomposition. It is the sum of the modulus of the coefficients of the decomposition. In the abstract of the paper of the 1st question, they say that the sampling overhead scales as 9^n, where "n" is the number of non-local gates that we decompose. Why 9? In the sum we have 10 terms with coefficient 1/2. So it should be 25(5^2). $\endgroup$
    – stopper
    Sep 30, 2022 at 13:48
  • 1
    $\begingroup$ In the decomposition derived by Mitarai and Fujii, there are only six terms (cf. Fig. 1 or Fig. 2), each with prefactor with absolute value $1/2$, so $C = 6 \times 1/2 = 3$, in which case the sampling overhead for $n$ such two-qubit operations is $3^{2n} = 9^n$, in agreement with Endo, Benjamin and Li. A more detailed explanation of this scaling is provided by Mitarai and Fujii in the proof of Theorem 5 (cf. page 6 of main text and Appendix E). $\endgroup$
    – bm442
    Oct 3, 2022 at 16:31
  • $\begingroup$ Thank you for the answer. So they consider the term inside the sum symbol as a unique term repeated 4 times (for the 4 values of alpha)? I thought I should have considered each term inside the sum separately, and in this way I obtain 10 terms in total. Thank you for the clarification $\endgroup$
    – stopper
    Oct 4, 2022 at 9:59

I contacted the author of the paper, and he told me that the equality holds when we consider the operations as quantum channels, which leads to the equal expectation values. So the channel representation of CZ gate and the decomposed one should be the same.

I still don't know what is a channel representation of a quantum gate.

  • $\begingroup$ Ask a new question giving this context to get some visibility on it! :) $\endgroup$
    – Tristan Nemoz
    Sep 28, 2022 at 17:02

Just for completeness I'll add how to apply this CZ decomposition. One must apply each summand individually to the state and add them afterwards. Let's assume we want to apply the decomposed CZ to an arbitrary state $\newcommand{\ket}[1]{\vert#1\rangle}\newcommand{\bra}[1]{\langle#1\vert}\ket{\psi}=\left[\begin{matrix}a&b&c&d\end{matrix}\right]^\top$, $\rho=\ket{\psi}\bra{\psi}$. Then we have to individually compute

$$\frac{1}{2}e^{i\pi Z/4}\otimes e^{i\pi Z/4}\ket{\psi}$$ $$\frac{1}{2}e^{-i\pi Z/4}\otimes e^{-i\pi Z/4}\ket{\psi}$$

and also each individual summand in this sum

$$-\frac{1}{2}\sum_{\alpha\in\{\pm 1\}^2}\alpha_1\alpha_2\left((I+\alpha_1 Z)/2\otimes e^{i(\alpha_2 + 1)\pi Z/4}\ket{\psi} + e^{i(\alpha_1 + 1)\pi Z/4}\otimes (I+\alpha_2 Z)/2\ket{\psi}\right)$$

After that we need to compute the mixed density matrix of all the states which is the same for the normal CZ operator

$$CZ~\rho~CZ^\dagger = \left[\begin{matrix}a \overline{a} & a \overline{b} & a \overline{c} & - a \overline{d}\\b \overline{a} & b \overline{b} & b \overline{c} & - b \overline{d}\\c \overline{a} & c \overline{b} & c \overline{c} & - c \overline{d}\\- d \overline{a} & - d \overline{b} & - d \overline{c} & d \overline{d}\end{matrix}\right]$$

Here is a python implementation, since I didn't bother to write out all summands above.

from sympy import *
from sympy.physics.quantum import TensorProduct
from sympy.physics.quantum.gate import Z

z = Z().get_target_matrix() # define the Z gate
cz = Matrix([[1, 0, 0, 0], [0, 1, 0, 0], [0, 0, 1, 0], [0, 0, 0, -1]]) # define the CZ gate
state = Matrix(symbols("a b c d")) # define the state
rho = state * state.H # define the state's density matrix

# apply each summand to the state's density matrix
t1 = TensorProduct(Matrix(exp(I * pi * z / 4)), Matrix(exp(I * pi * z / 4)))
s1 = t1 * rho * t1.H / 2
t2 = TensorProduct(Matrix(exp(-I * pi * z / 4)), Matrix(exp(-I * pi * z / 4)))
s2 = t2 * rho * t2.H / 2
a1, a2 = symbols("a1 a2")
t3 = TensorProduct((eye(2) + a1 * z) / 2, Matrix(exp(I * (a2 + 1) * pi * z / 4)))
t4 = TensorProduct(Matrix(exp(I * (a1 + 1) * pi * z / 4)), (eye(2) + a2 * z) / 2)
s3 = t3.subs({a1: -1, a2: -1}) * rho * t3.subs({a1: -1, a2: -1}).H / 2
s4 = t4.subs({a1: -1, a2: -1}) * rho * t4.subs({a1: -1, a2: -1}).H / 2
s5 = -t3.subs({a1: -1, a2: 1}) * rho * t3.subs({a1: -1, a2: 1}).H / 2
s6 = -t4.subs({a1: -1, a2: 1}) * rho * t4.subs({a1: -1, a2: 1}).H / 2
s7 = -t3.subs({a1: 1, a2: -1}) * rho * t3.subs({a1: 1, a2: -1}).H / 2
s8 = -t4.subs({a1: 1, a2: -1}) * rho * t4.subs({a1: 1, a2: -1}).H / 2
s9 = t3.subs({a1: 1, a2: 1}) * rho * t3.subs({a1: 1, a2: 1}).H / 2
s10 = t4.subs({a1: 1, a2: 1}) * rho * t4.subs({a1: 1, a2: 1}).H / 2

normal_cz = simplify(cz * rho * cz.H) # density matrix of CZ applied to the state
decomp_cz = simplify(s1 + s2 - (s3 + s4 + s5 + s6 + s7 + s8 + s9 + s10)) # add the summands
assert decomp_cz.equals(normal_cz) # check if the two results are equal

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.