Rank - 1 decomposition to visualize Thorp shuffling graph

(* swap p'th pair of the list *)
swap[l_,p_Integer]:=Module[{ll=l},ll[[{2p-1,2p}]]=ll[[{2p,2p-1}]];ll];
(* swap i'th pair for every i where p[i]=True *)
swap[l_,p_List]:=Fold[swap[#1,#2]&,l,Position[p,True]//Flatten]
(* generate all permutations resulting from one shuffle of given permutation *)
thorp[l_]:=Module[{k,ll},
k=Length[l]/2;
ll=Riffle[l[[;;k]],l[[k+1;;]]];
FromDigits[swap[ll,#]]&/@Tuples[{False,True},k]
]
(* generate all the arrows of the thorp graph *)
thorpEdges=Thread[FromDigits[#]->thorp[#]]&/@Permutations[Range[4]]//Flatten;

(* converts list of edges in permutation graph to adjacency matrix, canonical permutation numbering *)
edges2mat[edges_,size_]:=Module[{perms=Permutations[Range[size]]},
Array[If[MemberQ[edges,FromDigits[perms[[#1]]]->FromDigits[perms[[#2]]]],1,0]&,{Length[perms],Length[perms]}]
];

(* adjacency matrix of Thorp's shuffling graph *)
thorpMat=edges2mat[thorpEdges,4];
GraphPlot[thorpEdges,DirectedEdges->True,VertexLabeling->True]

GraphPlot3D[thorpMat,EdgeRenderingFunction->(Cylinder[#1,.04]&),VertexRenderingFunction->({Black,Sphere[#1,.1]}&),PlotStyle->Directive[Specularity[White,20]]]

r=10;
anim=Table[Show[g,ViewPoint->{r Cos[t],r Sin[t],2},ViewAngle->Pi/30,Boxed->False],{t,0,2Pi,2Pi/40}];
Export["C:\\research\\ideas\\shuffling\\thorp.gif",anim]

(* rank-1 decomposition of the graph *)
{u,d,v}={#1,Diagonal[#2],#3}&@@(SingularValueDecomposition[thorpMat]);
rank1mats=Table[d[[i]]*Outer[Times,u[[All,i]],v[[All,i]]],{i,1,6}];
MatrixPlot[Total[rank1mats]]
MatrixPlot/@rank1mats//Partition[#,3]&//GraphicsGrid

(* converts adjacency matrix of a permutation graph to list of edges, restoring permutation labels *)
mat2edges[mat_,size_]:=Module[{perms,edges,mapping},
perms=FromDigits[#]&/@Permutations[Range[size]];
edges=ArrayRules[mat]//Cases[#,p:(x_->1)->x]&//Rule@@@#&;
mapping=Thread[Range[Length[perms]]->perms];
edges/.mapping
];
GraphPlot[mat2edges[#,4],VertexLabeling->True,DirectedEdges->True]&/@rank1mats//Partition[#,3]&//GraphicsGrid

Merging nodes to see general structure

(* find indices of nodes that have positive out degree *)
outNodes[mat_]:=MemberQ[#,1]&/@mat//Position[#,True]&//Flatten
partitions=outNodes[#]&/@rank1mats;
(* takes adjacency matrix of a graph, partition of indices, merges those partitions into single nodes *)
mergeGraph[adjMat_,partitions_]:=Module[{l=Length[partitions]},
Table[Unitize[Total[adjMat[[partitions[[i]],partitions[[j]]]],2]],{i,1,l},{j,1,l}]
];
mergeGraph[thorpMat,partitions]//GraphPlot[#,DirectedEdges->True]&

p2=Permutations[partitions][[66]];
rules[i_]:=Thread[p2[[i]]->Table[{Cos[i 2Pi/6],Sin[i 2Pi/6]}+2/20{Cos[j 2Pi/4],Sin[j  2Pi/4]},{j,1,4}]];
allrules=rules[#]&/@Range[6]//Flatten;
GraphPlot[thorpMat,DirectedEdges->True,VertexLabeling->True,VertexCoordinateRules->allrules]