# How to represent as a matrix the CSWAP on non-adjacent qubits?

I was reading this post recently on the CNOT gate between non-adjacent qubits in a 3 qubit system.

And the accepted answer generalizes to general controlled unitaries, saying that any controlled operation is really just

$$CU = | 0 \rangle \langle 0 | I + | 1 \rangle \langle 1 | U$$

Where $$I$$ is the identity matrix with the same dimension as $$U$$.

And the author later explains how this generalizes to non-adjacent qubits.

But in the case of CSWAP, this isn't entirely clear. This is because the SWAP gate itself (the unitary) can act on non-adjacent qubits.

So how is the topic of non-adjacent qubits dealt with in the SWAP case?

It's pretty much the same thing.

Say you have an $$n$$-qubit register, and some $$U$$ acting (nontrivially) on $$k qubits. Say you want to write the action of the associated controlled operation $$CU$$, where the control qubit is the $$i$$-th, and the target qubits have labels $$J\equiv (j_1,...,j_k)$$ (and thus obviously we also assume $$i\neq j_k$$ for all $$k$$). An easy approach is to first assume $$i=1$$ and $$j_1=2,j_2=3,..., j_k=k+1$$, that is, use the first qubit as control, and the other qubits in the natural order. With this choice we simply have $$CU = |0\rangle\!\langle0| \otimes I_{n-1} + |1\rangle\!\langle1|\otimes U\otimes I_{n-k-1},$$ which you presumably know how to compute. You can now use a permutation matrix to rearrange all qubits in the order you want them. So using the assumptions above for control and target qubits, we get $$CU_{i,J} = P (|0\rangle\!\langle0| \otimes I_{n-1} + |1\rangle\!\langle1|\otimes U\otimes I_{n-k-1}) P, \\ P|\ell_1,...,\ell_n\rangle = |\pi(\ell_1),...,\pi(\ell_n)\rangle,$$ with $$\pi\in S_n$$ the permutation such that $$\pi(1)=i$$, $$\pi(2)=j_1,...,\pi(k+1)=j_k$$, and acts arbitrarily on the rest of the inputs.

Formalism aside, the action of the permutation matrices amounts to simply shuffling around the qubits to place control and targets where suitable. Taking the CSWAP as a concrete example, we start writing a decomposition for it assuming first qubit is control, and the targets are in natural succession: $$\mathrm{CSWAP} = |0\rangle\!\langle0|\otimes I_2\otimes I_{n-3}+ |1\rangle\!\langle1| \otimes\mathrm{SWAP}\otimes I_{n-3}.$$ It's also useful to write the SWAP itself as a sum of local operations: $$\mathrm{SWAP} = \mathbb{P}_{00} + E_{01}\otimes E_{10} + E_{10}\otimes E_{01} + \mathbb{P}_{11}, \qquad \mathbb{P}_{00}\equiv |00\rangle\!\langle 00|, E_{01}\equiv |0\rangle\!\langle 1|.$$ Now permute the spaces as suitable to get the correct control and targets. For example, if we have $$n=4$$, with the fourth qubit as control and the first two qubits as targets, the action of the permutation matrices will give you $$I_3\otimes |0\rangle\!\langle0| + \mathrm{SWAP}\otimes I\otimes |1\rangle\!\langle1|,$$ from which you easily get the matrix representation. Or if control is (again) the fourth qubit, but now the targets are first and third qubits, then you'll get $$I_3\otimes|0\rangle\!\langle 0| + \mathrm{SWAP}_{1,3}\otimes|1\rangle\!\langle1|,$$ where $$\mathrm{SWAP}_{1,3}$$ is the three-qubit unitary swapping first and third qubit, which we can write explicitly taking the decomposition given above for the SWAP, and permitting the spaces: $$\mathrm{SWAP}_{1,3} = |0\rangle\!\langle0|\otimes I\otimes |0\rangle\!\langle0| + E_{01}\otimes I \otimes E_{10} + E_{10} \otimes I \otimes E_{01} + |1\rangle\!\langle1| \otimes I \otimes |1\rangle\!\langle1|.$$ Practically speaking, in most cases it's simpler to just write the permutation matrices and apply them to the gate.

As per how you'd actually implement this, there's of course several approaches. Here's a Mathematica snippet I wrote many years ago to achieve this. Looking back at this, it's probably that you can make this more efficient by multiplying by permutation matrices as discussed in the post above, rather than shuffling indices as done in the code; still, it works in at least simple scenarios, so I haven't had the incentive to change it:

QBasePermutation[
matrix_?MatrixQ,
basisLengths : {__Integer},
newIndices : {__Integer}
] := Module[{matrixAsTP, transposedTP},
(* Convert matrix to TensorProduct structure *)
matrixAsTP = ArrayReshape[matrix, Join[#, #] & @ basisLengths];
(* Properly transpose the indices *)
transposedTP = Transpose[matrixAsTP,
Join[newIndices, newIndices + Length @ basisLengths]
];
(* Convert back into matrix structure *)
Flatten[transposedTP,
{#, # + Length @ basisLengths}& @ Range @ Length @ basisLengths
]
];

ProjectionMatrix[dim_Integer, y_Integer, x_Integer] := Normal @ SparseArray[
{{y, x} -> 1},
{dim, dim}
];

ProjectionMatrix[dim_Integer, x_Integer] := ProjectionMatrix[dim, x, x];

KP[x_] := x;
KP[x___] := KroneckerProduct @ x;

QControlledGate[numQubits_Integer,
controlQubit_Integer,
targetQubits_List,
gateMatrix_] /; (
And[
1 <= controlQubit <= numQubits,
Sequence @@ Thread[1 <= targetQubits <= numQubits],
]
) := Module[{gate},
gate = Plus[
KP[ProjectionMatrix[2, 1, 1], IdentityMatrix[2 ^ Length @ targetQubits]],
KP[ProjectionMatrix[2, 2, 2], gateMatrix]
];
gate = KP[gate, IdentityMatrix[2 ^ (numQubits - Length @ targetQubits - 1)]];

QBasePermutation[gate,
ConstantArray[2, numQubits],
{controlQubit, Sequence @@ targetQubits,
Sequence @@ Complement[
Range @ numQubits, Append[targetQubits, controlQubit]
]
}
]
];

QControlledGate[
numQubits_Integer,
controlQubit_Integer,
targetQubits_Integer,
gateMatrix_
] := QControlledGate[numQubits, controlQubit, {targetQubits}, gateMatrix];

swap = {{1, 0, 0, 0},
{0, 0, 1, 0},
{0, 1, 0, 0},
{0, 0, 0, 1}};