Quantum gates and computation.

My notations for this notebook : ket[{0}] and ket[{1}] stand for Dirac’s ket representation, |0> and |1>.  Their matrix representation are |0>={1,0} and |1>={0,1}.
More generally : ket[{0,1,0}] stands for Dirac’s ket representation, |010>.

Ket to Matrix represention :

BasisKetToMat[ket_]:=ReplacePart[Table[0,{k,2^Length[ket[[1]]]}],1,1+FromDigits[ket[[1]],2]]

KTM[x_  ket[y_]]:=x BasisKetToMat[ket[y]]

KTM[ket[y_]]:= BasisKetToMat[ket[y]]

KetToMat[ket_]:=Distribute[KTM[ket]]

(*3 Examples : *)     {KetToMat[  ket[{1}]],KetToMat[  ket[{0,1}]],KetToMat[  ket[{0,1}]+(1/(Sqrt[2])) ket[{1,0}]+a  ket[{1,1}]]}

Matrix to Ket  representation :

MatToKet[list_]:=Sum[{list[[k]] ket[IntegerDigits[Position[UnitVector[Length[list],k],1][[1,1]]-1,2,Log[2,Length[UnitVector[Length[list],k]]]]] },{k,Length[list]}][[1]]

(*3 Examples : *)    {MatToKet[{0,1}],MatToKet[{0,1,0,0}],MatToKet[{0,1,(1/(Sqrt[2])),a}]}

External product |ket><bra| via matrix representation :

x_ ⊗ y_:=ArrayFlatten[Outer[Times,x,y]]

(* Example : *)    (Subscript[D, i_]):= KetToMat[ket[{i}]]⊗KetToMat[ket[{i}]]          (* Equivalent to |1><1| in matrix form *)

{(Subscript[D, 0]),(Subscript[D, 1]),(Subscript[D, 0])//MatrixForm ,(Subscript[D, 1])//MatrixForm }                (* Representations of '0' or '1' detectors *)

{(Subscript[D, 0]).{1,0},(Subscript[D, 1]).{1,0},(Subscript[D, 0]).{0,1},(Subscript[D, 1]).{0,1}}      (* Confirms that (Subscript[D, i]) detects |i> with probability one *)

Tensor product |ket1>|ket2> via matrix representation :

Flatten[KetToMat[ket[{1}]]⊗KetToMat[ket[{1}]]]     (*The same "Outer" plus "Flatten" does the expected job : there is little syntactic diference between both products *)

MatToKet[Flatten[KetToMat[ket[{1}]]⊗KetToMat[ket[{1}]]]]        (* In agreement with  |1>|1> = |11>  *)

Some useful gates acting on a single qubit : Hadamard (H), Dephasing (F), Not (M) :

H=(1/(Sqrt[2]))(({{1, 1}, {1, -1}}));Φ[φ_]:=(({{1, 0}, {0, (E^(I φ))}}));M=(({{0, 1}, {1, 0}}));

Example : a balanced Mach - Zehnder interferometer :

H.M.H//MatrixForm

Pauli and rotation matrices :

(Subscript[σ, x])=(({{0, 1}, {1, 0}}));(Subscript[σ, y])=(({{0, -I}, {I, 0}}));(Subscript[σ, z])=(({{1, 0}, {0, -1}}));

Rx[θ_]:=MatrixExp[I θ (Subscript[σ, x])];Ry[θ_]:=MatrixExp[I θ (Subscript[σ, y])];Rz[θ_]:=MatrixExp[I θ (Subscript[σ, z])];

Preparation of a single qubit in an arbitrary state (starting with | 0 >= {1, 0}) ; U does the job up to an inessential phase :

U[φ_,θ_]:= Φ[φ+π/2].H.Φ[θ].H

U[φ,θ]==(E^(I θ/2)) (({{Cos[θ/2], -I Sin[θ/2]}, {(E^(I φ))  Sin[θ/2], I (E^(I φ)) Cos[θ/2]}}))//Simplify

Every gate can be reduced to H  and F (H  and F are universal for a single qubit) :

Rx[θ]==(E^(I θ)) H.Φ[-2θ].H//Simplify

Ry[θ]==(E^(I θ)) Φ[π/2].H.Φ[-2θ].H.Φ[-π/2]//Simplify

{Rx[θ],Ry[θ]}

Altering the polarization state of a single photon with a plate (Subscript[L, β,δ]) :

L[β_,δ_]:= (({{(E^(I δ/2))((Cos[β])^2)+(E^(-I δ/2))((Sin[β])^2), I Sin[2β] Sin[δ/2]}, {I Sin[2β] Sin[δ/2], (E^(-I δ/2))((Cos[β])^2)+(E^(I δ/2))((Sin[β])^2)}}))

Typical examples (up to an inessential phase) :

L[0,δ]//MatrixForm           (*A dephasing gate*)

L[π/4,π/2]//Simplify//MatrixForm           (*Sqrt of Not gate which transforms linear polarization to circular polarization*)

L[π/4,π/2].L[π/4,π/2]//ExpToTrig//MatrixForm        (*Not gate, up to factor i, results of the squaring of the former gate*)

L[π/8,π]//Simplify//MatrixForm              (*Hadamard gate*)

L[π/4,π]//Simplify//MatrixForm           (*Not gate*)

Special case : rotating the linear polarization of a photon (Zeno effect) :

(Subscript[D, 0]).{1,0}

Ry[θ]//MatrixForm

Abs[(Subscript[D, 0]).MatrixPower[Ry[ (π/(2n))],n].{1,0}//FullSimplify//PowerExpand]

Abs[(Subscript[D, 0]).MatrixPower[(Subscript[D, 0]).Ry[ (π/(2n))],n].{1,0}//FullSimplify//PowerExpand]

Gates acting on registers (at least two qubits). This needs tensor products  :

IdentityMatrix[2]⊗ {{a,b},{d,e}}//MatrixForm

The Controlled-Not gate (C-Not) :

CNot=(Subscript[D, 0])⊗IdentityMatrix[2]+(Subscript[D, 1])⊗M;

CNot//MatrixForm

The Swap gate (Sw) :

Sw=CNot.(IdentityMatrix[2]⊗(Subscript[D, 0])+M⊗(Subscript[D, 1])).CNot

The Controlled-Controlled-Not gate (CC-Not) :

CCNot=((Subscript[D, 0])⊗(Subscript[D, 0]))⊗IdentityMatrix[2]+((Subscript[D, 0])⊗(Subscript[D, 1]))⊗IdentityMatrix[2]+((Subscript[D, 1])⊗(Subscript[D, 0]))⊗IdentityMatrix[2]+((Subscript[D, 1])⊗(Subscript[D, 1]))⊗M;

CCNot//MatrixForm

Preparation of a 3 qubits register starting with |000> :

((Ry[-π/4]⊗Ry[-π/4])⊗Ry[-ArcTan[2]]).{1,0,0,0,0,0,0,0}//MatrixForm           (*Encodes the sequence {1,2,1,2,1,2,1,2} for future use*)

MatToKet[((Ry[-π/4]⊗Ry[-π/4])⊗Ry[-ArcTan[2]]).{1,0,0,0,0,0,0,0}]

(∑_(k=0)^7)(1+Mod[k,2])(u^k)                      (*The factorization of the associated generating function detects if the register is separable*)

Factor[(∑_(k=0)^7)(1+Mod[k,2])(u^k)]               (*The register is of course separable*)

Sum[(u^k)(1+Mod[k,4]),{k,0,15}]

Factor[Sum[(u^k)(1+Mod[k,4]),{k,0,15}]]             (*The corresponding register is not fully separable*)

Detection projectors for a particular qubit in the register (two qubits in the example) :

{(Subscript[D, 0])⊗(Subscript[D, 0])//MatrixForm,(Subscript[D, 0])⊗(Subscript[D, 1])//MatrixForm,(Subscript[D, 1])⊗(Subscript[D, 0])//MatrixForm,(Subscript[D, 1])⊗(Subscript[D, 1])//MatrixForm}

Controlled-U gates (C-SqrtNot in the example) :

SqrtNot=(E^(-I π/4)) L[π/4,π/2]//ExpToTrig

SqrtNot==(E^(-I π/4))Φ[π/2].H.Φ[π/2]//Simplify

CSqrtNot=(Subscript[D, 0])⊗IdentityMatrix[2]+(Subscript[D, 1])⊗SqrtNot

C-Phase and C-H gates :

CΦ[φ_]:=(Subscript[D, 0])⊗IdentityMatrix[2]+(Subscript[D, 1])⊗Φ[φ]         (*Controlled Phase gate*)

CΦ[φ]//MatrixForm

CH=(Subscript[D, 0])⊗IdentityMatrix[2]+(Subscript[D, 1])⊗H //MatrixForm         (*Controlled Hadamard gate*)

Reduction of the C-Not gate to Hadamard and C-Phase gates :

CNot==(IdentityMatrix[2]⊗H).CΦ[π/2].CΦ[π/2].(IdentityMatrix[2]⊗H)//Simplify

Reduction of the CC-Not (3 qubits) to Hadamard, CNot and C-Phase gates :

CCNot==(IdentityMatrix[4]⊗H).((Subscript[D, 0])⊗IdentityMatrix[4]+((Subscript[D, 1])⊗IdentityMatrix[2])⊗Φ[π/2]).(CNot⊗IdentityMatrix[2]).(IdentityMatrix[2]⊗CΦ[-π/2]).(CNot⊗IdentityMatrix[2]).(IdentityMatrix[2]⊗CΦ[π/2]).(IdentityMatrix[4]⊗H)

CCNot + CNot gates acting as a binary  half-adder

halfadder[ket[{x_,y_,0}]]:=MatToKet[(CNot⊗IdentityMatrix[2]).CCNot.KetToMat[ket[{x,y,0}]]]

{halfadder[ket[{0,0,0}]],halfadder[ket[{0,1,0}]],halfadder[ket[{1,0,0}]],halfadder[ket[{1,1,0}]]}

Identification algorithm among 4 (out of 16 existing) operators (reversibly) calculating Boolean functions in two variables (following the quantum scheme : |a>|b>|0> -> |a>|b>|f(a,b)>)

and=CCNot;xlessthany=CCNot.(IdentityMatrix[2]⊗CNot);xgreaterthany=CCNot.((Subscript[D, 1])⊗IdentityMatrix[4]+((Subscript[D, 1])⊗IdentityMatrix[2])⊗M);nor=(IdentityMatrix[4]⊗M).((Subscript[D, 0])⊗IdentityMatrix[4]+((Subscript[D, 1])⊗IdentityMatrix[2])⊗M).CCNot.(IdentityMatrix[2]⊗CNot);

Table[MatToKet[nor.KetToMat[ket[{i,j,0}]] ],{i,0,1},{j,0,1}]     (*The result translated with kets : the third qubit is the Nor of the two first*)

Deutsch[op_]=((H⊗H)⊗IdentityMatrix[2]).(CΦ[π]⊗IdentityMatrix[2]).((Φ[π]⊗Φ[π])⊗IdentityMatrix[2]).((H⊗H)⊗IdentityMatrix[2]).op.((H⊗H)⊗H);       (*Identification algorithm in one pass between 4 boolean predefined functions*)

Deutsch[and].KetToMat[ket[{0,0,1}]]     (*The 2 first qubits of the result,|11>(1/(Sqrt[2]))(|0>-|1>), is the signature of and*)

Deutsch[nor].KetToMat[ket[{0,0,1}]]        (*The 2 first qubits of the result,|00>(1/(Sqrt[2]))(|0>-|1>), is the signature of nor*)

Deutsch[xlessthany].KetToMat[ket[{0,0,1}]]     (*The 2 first qubits of the result,|01>(1/(Sqrt[2]))(|0>-|1>), is the signature of a<b*)

Deutsch[xgreaterthany].KetToMat[ket[{0,0,1}]]     (*The 2 first qubits of the result,|10>(1/(Sqrt[2]))(|0>-|1>), is the signature of a>b*)

Grover algorithm   (Lorsque n=2 qubits sont en jeu une seule itération suffit à localiser l’élément intrus (incarné par un changement de signe : -1 dans une bande de 1 (Une porte CCNot trafiquée fait l’affaire)) sinon il en faut Sqrt(n).

2(Subscript[D, 0])-IdentityMatrix[2]

2(Subscript[D, 1])-IdentityMatrix[2]

Grover[op_]=((H⊗H)⊗IdentityMatrix[2]).((2(Subscript[D, 0])⊗(Subscript[D, 0])-IdentityMatrix[4])⊗IdentityMatrix[2]).((H⊗H)⊗IdentityMatrix[2]).op.((H⊗H)⊗H);

CCNot=((Subscript[D, 0])⊗(Subscript[D, 0]))⊗IdentityMatrix[2]+((Subscript[D, 0])⊗(Subscript[D, 1]))⊗IdentityMatrix[2]+((Subscript[D, 1])⊗(Subscript[D, 0]))⊗IdentityMatrix[2]+((Subscript[D, 1])⊗(Subscript[D, 1]))⊗M;   (*Pour rappel*)

Grover[CCNot].KetToMat[ket[{0,0,1}]]      (*The 2 first qubits of the result,|11>(1/(Sqrt[2]))(|0>-|1>), indicate which item produces a change of sign (11 in the case of the CCNot gate)*)

Grover[((Subscript[D, 0])⊗(Subscript[D, 0]))⊗IdentityMatrix[2]+((Subscript[D, 0])⊗(Subscript[D, 1]))⊗IdentityMatrix[2]+((Subscript[D, 1])⊗(Subscript[D, 0]))⊗M+((Subscript[D, 1])⊗(Subscript[D, 1]))⊗IdentityMatrix[2]].KetToMat[ket[{0,0,1}]]

Grover[((Subscript[D, 0])⊗(Subscript[D, 0]))⊗IdentityMatrix[2]+((Subscript[D, 0])⊗(Subscript[D, 1]))⊗M+((Subscript[D, 1])⊗(Subscript[D, 0]))⊗IdentityMatrix[2]+((Subscript[D, 1])⊗(Subscript[D, 1]))⊗IdentityMatrix[2]].KetToMat[ket[{0,0,1}]]

Grover[((Subscript[D, 0])⊗(Subscript[D, 0]))⊗M+((Subscript[D, 0])⊗(Subscript[D, 1]))⊗IdentityMatrix[2]+((Subscript[D, 1])⊗(Subscript[D, 0]))⊗IdentityMatrix[2]+((Subscript[D, 1])⊗(Subscript[D, 1]))⊗IdentityMatrix[2]].KetToMat[ket[{0,0,1}]]

Shor algorithm

Quantum circuit for the Discrete Fourier Transform :

(Subscript[DFT, 2])=Sw.(IdentityMatrix[2]⊗H) .CΦ[π/2].(H⊗IdentityMatrix[2]);     (*valid for a two-qubits register : be careful with the order*)

(Subscript[DFT, 2])==(1/(Sqrt[(2^2)]))Table[Exp[((2I π)/(2^2))j k],{j,0,(2^2)-1},{k,0,(2^2)-1}]

Probabilities of measuring one of the 4 possible output-states (in the order 00 01 10 11) :

resu=(Subscript[DFT, 2]).KetToMat[(Subscript[s, 0])  ket[{0,0}] +(Subscript[s, 1]) ket[{0,1}]+(Subscript[s, 2])  ket[{1,0}]+(Subscript[s, 3])  ket[{1,1}]]

resu Conjugate[resu]//ComplexExpand//Simplify

(Subscript[DFT, 3])=((Sw⊗IdentityMatrix[2]).(IdentityMatrix[2]⊗Sw).(Sw⊗IdentityMatrix[2])).(IdentityMatrix[4]⊗H).(IdentityMatrix[2]⊗CΦ[π/2]).((Subscript[D, 0])⊗IdentityMatrix[4]+((Subscript[D, 1])⊗IdentityMatrix[2])⊗Φ[π/4] ).((IdentityMatrix[2]⊗H)⊗IdentityMatrix[2]).(CΦ[π/2]⊗IdentityMatrix[2]).(H⊗IdentityMatrix[4]);   (*valid for a three-qubits register : be careful with the order*)

(Subscript[DFT, 3])==(1/(Sqrt[(2^3)]))Table[Exp[((2I π)/(2^3))j k],{j,0,(2^3)-1},{k,0,(2^3)-1}]

resu=(Subscript[DFT, 3]).KetToMat[(Subscript[s, 0])  ket[{0,0,0}] +(Subscript[s, 1])  ket[{0,0,1}] +(Subscript[s, 2])  ket[{0,1,0}] +(Subscript[s, 3])  ket[{0,1,1}] +(Subscript[s, 4])  ket[{1,0,0}] +(Subscript[s, 5])  ket[{1,0,1}] +(Subscript[s, 6])  ket[{1,1,0}] +(Subscript[s, 7])  ket[{1,1,1}]]//Factor

Probabilities of measuring one of the 8 possible output-states (in the order 000 001 010 011 100 101 110 111) :

resu Conjugate[resu]//ComplexExpand//Simplify

ListPlot[Chop[resu ComplexExpand[Conjugate[resu]]/.{(Subscript[s, 0]) ->(1/(Sqrt[33])),(Subscript[s, 1]) ->(2/(Sqrt[33])),(Subscript[s, 2]) ->(3/(Sqrt[33])),(Subscript[s, 3]) ->(1/(Sqrt[33])),(Subscript[s, 4]) ->(2/(Sqrt[33])),(Subscript[s, 5]) ->(3/(Sqrt[33])),(Subscript[s, 6]) ->(1/(Sqrt[33])),(Subscript[s, 7]) ->(2/(Sqrt[33]))}//N],Frame->True,Filling->Bottom,FillingStyle->Directive[Thick]]

Chop[resu ComplexExpand[Conjugate[resu]]/.{(Subscript[s, 0]) ->(1/(Sqrt[20])),(Subscript[s, 1]) ->(2/(Sqrt[20])),(Subscript[s, 2]) ->(1/(Sqrt[20])),(Subscript[s, 3]) ->(2/(Sqrt[20])),(Subscript[s, 4]) ->(1/(Sqrt[20])),(Subscript[s, 5]) ->(2/(Sqrt[20])),(Subscript[s, 6]) ->(1/(Sqrt[20])),(Subscript[s, 7]) ->(2/(Sqrt[20]))}//N]

ListPlot[Chop[resu ComplexExpand[Conjugate[resu]]/.{(Subscript[s, 0]) ->(1/(Sqrt[20])),(Subscript[s, 1]) ->(2/(Sqrt[20])),(Subscript[s, 2]) ->(1/(Sqrt[20])),(Subscript[s, 3]) ->(2/(Sqrt[20])),(Subscript[s, 4]) ->(1/(Sqrt[20])),(Subscript[s, 5]) ->(2/(Sqrt[20])),(Subscript[s, 6]) ->(1/(Sqrt[20])),(Subscript[s, 7]) ->(2/(Sqrt[20]))}//N],PlotRange->All,Frame->True,Filling->Bottom,FillingStyle->Directive[Thick]]

Total[{0.8522727272727273`,0.005087673297377347`,0.018939393939393936`,0.047942629732925686`,0.003787878787878788`,0.047942629732925665`,0.018939393939393936`,0.005087673297377346`}]

Random Quantum Walk

(*The following recursion computes the Feynman path-sum corresponding to the transit of a photon '0' in a triangular network of half-silvered plates. Additional phases are added so that the only alteration is a sign inversion when a plate is crossed in the South-North direction*)

c[0]=f[0,0];c[i_]:=c[i-1]/.f[x_,c_]->((f[x-1,0]+(((-1))^c)f[x+1,1])/(Sqrt[2]))//Simplify

c[1]

c[2]

n=32;

(*Two additional functions are necessary in order to translate the results as detectors numbers*)

arg[a_ f_[x_,y_]]:=(n-1)y+(x+n)/2

treat[z_]:=If[ArrayQ[z],((n-1)z[[2]]+(z[[1]]+n)/2),z]

coeff=Table[Expand[c[n]][[k,1]],{k,2n}]

det=Table[treat[arg[2^n Expand[c[n]][[k]]]],{k,2n}]

probas=Flatten[Table[coeff[[Position[det,k][[1]]]],{k,0,2n-1}]]^2

ListPlot[probas,PlotJoined->True,PlotStyle->{Thick,Red},Frame->True,PlotRange->All]

f1=ListPlot[{1,7,9,-5,3,-3,3,1,1,5,1,-3,3,-1,-5,-1}^2/256,PlotJoined->True,PlotStyle->{Thick,Red},Frame->True]

f2=ListPlot[Table[Binomial[15,k]/2^15,{k,0,15}],PlotJoined->True,PlotStyle->{Thick,Blue},Frame->True]

Sum[Binomial[15,k]/2^15,{k,0,15}]

Show[{f1,f2}]