# How to write the $iSWAP$ unitary as a linear combination of tensor products between 1-qubit gates? [duplicate]

As far as I understood, it should always be possible to decompose any $$n$$-qubits unitary $$W$$ into a linear combination of tensor products between $$n$$ single-qubit unitaries $$U_i$$: $$W = \sum_k \lambda_k \left( \bigotimes_{i=1}^{n} U_i \right)$$

For example, the unitary matrix of the $$SWAP$$ gate can be written as follows: $$SWAP = \frac{1}{2} \left( I \otimes I + X \otimes X + Y \otimes Y + Z \otimes Z\right)$$

How can I write the unitary of the $$iSWAP$$ gate in the same form? And, more in general, is there a mechanical procedure to find such a decomposition for any given $$W$$?

• Your EDIT is another question. I would recommend to open new thread and post link to this one in order to preserve requirement one question per thread. Mar 24 at 5:51

ISWAP =
+ II * (0.5+0j)
+ XX * 0.5j
+ YY * 0.5j
+ ZZ * (0.5+0j)


Output is from this code:

iswap = np.array([
[1, 0, 0, 0],
[0, 0, 1j, 0],
[0, 1j, 0, 0],
[0, 0, 0, 1],
])
iswap_terms = matrix_to_pauli_terms(iswap)
print("ISWAP =")
for iswap_pauli_string, iswap_coefficient in iswap_terms.items():
print("    +", iswap_pauli_string, '*', repr(iswap_coefficient))
np.testing.assert_array_equal(pauli_terms_to_matrix(iswap_terms), iswap)


Using these utilities, which are not optimized but work on any matrix. The complexity here is $$16^n$$ where $$n$$ is the number of qubits. I think it should be possible to achieve $$n 4^n$$ instead:

import itertools
from typing import Dict, Iterable

import numpy as np

i = np.eye(2)
x = np.array([
[0, 1],
[1, 0],
])
z = np.array([
[1, 0],
[0, -1],
])
y = np.array([
[0, -1j],
[1j, 0],
])
paulis = {"I": i, "X": x, "Y": y, "Z": z}

def pauli_string_to_matrix(pauli_string: Iterable[str]) -> np.ndarray:
t = np.eye(1)
for p in pauli_string:
t = np.kron(t, paulis[p])
return t

def matrix_to_pauli_terms(matrix: np.ndarray) -> Dict[str, complex]:
w, h = matrix.shape
n = w.bit_length() - 1
assert w == h and w == 2**n  # check matrix is square with power of 2 size

terms: Dict[str, complex] = {}
for pauli_tuple in itertools.product("IXYZ", repeat=n):
pauli_string = ''.join(pauli_tuple)
pauli_matrix = pauli_string_to_matrix(pauli_string)
coefficient = complex(np.dot(matrix.flat, np.conj(pauli_matrix.flat))) / 2**n
if coefficient:
terms[pauli_string] = coefficient
if not terms:
terms["I" * n] = 0
return terms

def pauli_terms_to_matrix(terms: Dict[str, complex]) -> np.ndarray:
assert terms
n = len(next(iter(terms.keys())))
total = np.zeros(shape=(2**n, 2**n), dtype=np.complex64)
for pauli_string, coefficient in terms.items():
total += pauli_string_to_matrix(pauli_string) * coefficient


Update: here's the fast version that takes $$O(n4^n)$$ time to compute the terms

def matrix_to_pauli_terms_fast(matrix: np.ndarray) -> Dict[str, complex]:
w, h = matrix.shape
n = w.bit_length() - 1
N = 2**n
assert w == h == N  # check matrix is square with power of 2 size

# Permute by xoring row coordinate into column coordinate
term_matrix = np.empty(shape=matrix.shape, dtype=np.complex64)
term_matrix[0, :] = matrix[0, :]
indices = np.array(range(N))
for k in range(1, N):
indices ^= k ^ (k - 1)
term_matrix[k, :] = matrix[k, indices]

# Hadamard transform the columns and account for scalar phase from Y.
term_matrix.shape = (2,) * (2 * n)
for k in range(n):
index: List[Union[slice, int]] = [slice(None)] * (2 * n)
index[k] = 0
a = tuple(index)
index[k] = 1
b = tuple(index)
term_matrix[a] += term_matrix[b]
term_matrix[b] *= -2
term_matrix[b] += term_matrix[a]
# Scalar phase from Y.
index[k + n] = 1
term_matrix[tuple(index)] *= 1j
term_matrix /= 2**n
term_matrix.shape = (N, N)

# Convert from dense matrix representation to sparse dict representation.
terms: Dict[str, complex] = {}
for pauli_tuple in itertools.product("IXYZ", repeat=n):
pauli_string = ''.join(pauli_tuple)
xk = sum(2 **k * (pauli_string[n-k-1] in 'XY') for k in range(n))
zk = sum(2 **k * (pauli_string[n-k-1] in 'YZ') for k in range(n))
coefficient = term_matrix[zk, xk]
if coefficient:
terms[pauli_string] = coefficient
if not terms:
terms["I" * n] = 0

return terms


The funny thing is that this fancy looking transformation of the matrix actually corresponds to a very familiar quantum circuit. Flattening the matrix into vector, what's happening is equivalent to applying the CNOT-and-H prefix of a bunch of bell basis measurements! So basically you convert from matrix to Pauli terms by Bell basis measurement. • +1 - something, something, teach a man about fish. Mar 24 at 19:54

You can check that $$\mathrm{iSWAP}=\frac12(I\otimes I+iX\otimes X+iY\otimes Y + Z\otimes Z)$$