Re: Counting

*To*: mathgroup at smc.vnet.net*Subject*: [mg114979] Re: Counting*From*: Leonid Shifrin <lshifr at gmail.com>*Date*: Sun, 26 Dec 2010 04:01:06 -0500 (EST)*References*: <201012200540.AAA22760@smc.vnet.net>

Yaroslav, Looks like I have produced a chain of incorrect solutions. Here is a modified one, which I hope is correct. At least, it is more clear what is going on there. The idea is that I permute all the list of permutations with a single group element (transformation) at once, this can be done very efficiently with the list[[All, permutation]] idiom. Then, I use an auxiliary function positionsOfSame, which finds pairs of positions of same permutations in the original list and the permuted one, so for example {{1,3},{2,4}} would mean that first permutation when transformed became third one in the original list, and second became fourth. Then, I use these pairs as a graph specification (edges), to find connected components - they will give the positions in the original list of all permutations which are equivalent to each other under a given transformation. I pick the first element (position) in each component, and delete from the permutation list all elements (permutations) which are at remaining positions in these connected components. This is done for each group element (transformation), one after another, with Fold. The result now seems to be independent on the order in which the transformations are applied. Starting with the longest cycles is an optimization, since the connected components are larger then (of the size of more or less the length of the cycle), which means that long permutation list is quickly reduced in size and other transformations work with smaller lists and are then faster. Here is the code: (* Find connected components and a list of positions in the permutation list, at which positions to delete elements *) Clear[getTree, listSplit, gatherBy, getComponents, getPositionsToDelete]; getTree = Compile[{{pairs, _Integer, 2}}, Module[{t = 0, i = 0, j = 0, xl = 0, yl = 0, k = 0, len = Max[pairs], dad = Table[0, {Max[pairs]}]}, For[k = 1, k <= Length[pairs], k++, xl = i = pairs[[k, 1]]; yl = j = pairs[[k, 2]]; If[xl == yl, Continue[]]; While[dad[[i]] > 0, i = dad[[i]]]; While[dad[[j]] > 0, j = dad[[j]]]; While[dad[[xl]] > 0, t = xl; xl = dad[[xl]]; dad[[t]] = i]; While[dad[[yl]] > 0, t = yl; yl = dad[[yl]]; dad[[t]] = j]; If[i != j, If[dad[[j]] <= dad[[i]], dad[[j]] += dad[[i]] - 1; dad[[i]] = j;,(*else*)(dad[[i]] += dad[[j]] - 1); dad[[j]] = i];];]; For[k = 1, k <= len, k++, If[dad[[k]] <= 0, dad[[k]] = k]]; dad]]; listSplit[x_List, lengths_List] := MapThread[Take[x, {##}] &, {Most[#], Rest[#] - 1}] &@ Accumulate[Prepend[lengths, 1]]; gatherBy[lst_List, flst_List] := listSplit[lst[[Ordering[flst]]], (Sort@Tally[flst])[[All, 2]]]; getComponents[pairs_] := With[{dad = getTree[pairs]}, gatherBy[Range[Length[dad]], FixedPoint[dad[[#]] &, dad]]]; (* Delete all elements except first, in each component *) getPositionsToDelete[pairs_] := Flatten[#[[All, 2 ;; -1]]] &@getComponents[pairs]; (* Functions to find correspondence between positions of same permutations in the original and transformed list of permutations *) Clear[less]; less = Compile[{{fst, _Integer, 1}, {sec, _Integer, 1}}, Module[{i = 1}, While[fst[[i]] == sec[[i]], i++]; fst[[i]] < sec[[i]]]]; Clear[positionsOfSameInSorted ]; positionsOfSameInSorted = Compile[{{fsorted, _Integer, 2}, {ssorted, _Integer, 2}}, Module[{i, j, pctr, positions = Table[{0, 0}, {Max[Length[fsorted], Length[ssorted]]}]}, For[i = 1; j = 1; pctr = 0, i <= Length[fsorted] && j <= Length[ssorted], If[fsorted[[i]] == ssorted[[j]], positions[[++pctr]] = {i, j}; i++; j++, (* else *) If[less[fsorted[[i]], ssorted[[j]]], i++, (* else *) j++ ]]]; Take[positions, pctr] ], CompilationOptions -> {"InlineCompiledFunctions" -> True, "InlineExternalDefinitions" -> True}]; Clear[positionsOfSame]; positionsOfSame = Compile[{{fst, _Integer, 2}, {sec, _Integer, 2}}, Module[{ord1, ord2, ps, psrt}, ord1 = Ordering[fst]; ord2 = Ordering[sec]; psrt = positionsOfSameInSorted[fst[[ord1]], sec[[ord2]]]; If[Length[psrt] != 0, ps = Transpose[ psrt]; Transpose[{ord1[[ps[[1]]]], ord2[[ps[[2]]]]}], (* else *) {{0}} ] ], CompilationOptions -> {"InlineCompiledFunctions" -> True, "InlineExternalDefinitions" -> True}]; (* Main functions *) Clear[reducePermutationsC]; reducePermutationsC = Compile[{{list, _Integer, 2}, {permutation, _Integer, 1}}, Module[{perms = list, pos = Table[-1, {Length[list]}], permuted = list, psame = Table[{-1, -1}, {Length[list]}]}, permuted = perms[[All, permutation]]; psame = positionsOfSame[perms , permuted]; If[psame != {{0}}, pos = getPositionsToDelete[psame]; If[pos != {}, perms = Delete[perms, List /@ pos]]]; perms], {{getPositionsToDelete[_], _Integer, 1}}, CompilationOptions -> {"InlineCompiledFunctions" -> True, "InlineExternalDefinitions" -> True}, CompilationTarget -> "C"]; Clear[nonequivalentPermutations]; nonequivalentPermutations[list_, group_] := With[{perms = Permutations[list], actions = DeleteCases[GroupElements@group, Cycles[{}]]}, With[{rlen = Range[Length[perms[[1]]]]}, Fold[reducePermutationsC, perms, Map[Permute[rlen, #] &, SortBy[actions, Length@#[[1]] &]]]]]; The syntax is now the same as in your solution. Here are examples of use: In[179]:= nonequivalentPermutations[{2,1,1,0},DihedralGroup[4]] Out[179]= {{2,1,1,0},{2,1,0,1}} In[180]:= nonequivalentPermutations[{2,2,2,2,2,2,2,1,1,1,1,1,1,0,0,0},DihedralGroup[16]]//Short//Timing Out[180]= {5.437,{{2,2,2,2,2,2,2,1,1,1,1,1,1,0,0,0},<<30098>>,{2,1,2,0,2,1,2,0,2,1,1,2,1,1,2,0}}} The resulting number of permutations is somewhat close to an estimate (30030) given to you by someone at StackOverflow. The connected component functionality is my port of implementation of weight-balancing path-compression union-find algorithm taken straight from the Sedgewick's book on algorithms in C, with some minor changes as needed for adoption to Mathematica. I also heavily optimized it for speed. The reason I did not compile it to C is that the bottleneck there is actually in the gatherBy function (rather than getTree), which can not alas be Compiled since it produces non-tensor structures (list of sublists of non-constant length). If it could, we could perhaps get another 30-40 % speed increase for the entire solution. But even this "home-brewed" gatherBy is about 4 times faster than the built-in GatherBy for this particular problem. Anyway, I really hope that now this is a correct solution, and this is then my last post on this topic :) Hope this helps. Regards, Leonid On Mon, Dec 20, 2010 at 8:40 AM, Yaroslav Bulatov <yaroslavvb at gmail.com>wrote: > I'd like to count the number of permutations of {2, 2, 2, 2, 2, 2, 2, > 1, 1, 0, 0, 0, 0, 0, 0, 0} that are not equivalent under the symmetry > of DihedralGroup[16]. In other words, count the ways of assigning > those integers to vertices of a 4 dimensional cube. > > This takes about a minute in another systme using "OrbitsDomain" command. > My > Mathematica approach is below, however it doesn't finish within 10 > minutes, any advice how to make it tractable? > > nonequivalentPermutations[lst_, group_] := ( > removeEquivalent[{}] := {}; > removeEquivalent[list_] := ( > Sow[First[list]]; > equivalents = Permute[First[list], #] & /@ GroupElements[group]; > DeleteCases[list, Alternatives @@ equivalents] > ); > > reaped = Reap@FixedPoint[removeEquivalent, Permutations@lst]; > reaped[[2, 1]] // Length > ); > nonequivalentPermutations[{2, 2, 2, 2, 2, 2, 2, 1, 1, 0, 0, 0, 0, 0, > 0, 0}, DihedralGroup[16]] > >

**References**:**Counting***From:*Yaroslav Bulatov <yaroslavvb@gmail.com>