82 views

Skip to first unread message

Dec 20, 2010, 12:40:14 AM12/20/10

to

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.

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]]

Dec 25, 2010, 2:34:54 AM12/25/10

to

Hi Yaroslav,

this is a follow - up to my previous post. I think that the code I just

posted grossly overestimates

the number of non-equivalent permutations, because I forgot about the

non-abelian nature of

permutations and my method explores only a part of the group space. I did

not figure out the

rigorous way to solve this problem (perhaps one should look more carefully

at the group structure),

but we can at least modify the code to give an upper bound on the number of

non-equivalent

permutations:

1. Modify the code for getUniquePermutatations, to apply the longest cycles

first

Clear[getUniquePermutatationsAlt];

getUniquePermutatationsAlt[permutations_, actions : {__Cycles}] /;

Length[permutations] > 0 :=

With[{rlen = Range[Length[permutations[[1]]]]},

Fold[reducePermutationsC, permutations,

Map[Permute[rlen, #] &, SortBy[actions, Length@#[[1]] &]]]]

2. Use this function inside FixedPoint:

In[80]:= FixedPoint[ getUniquePermutatationsAlt[#,actions]&,perms]//Timin=

g

Out[80]=

{5.016,{{0,0,0,0,0,0,0,2,1,2,2,2,2,1,2,2},{0,0,0,0,0,0,0,2,1,2,2,2,1,2,2,2}=

,{0,0,0,0,0,0,0,2,1,2,2,1,2,2,2,2},{0,0,0,0,0,0,0,2,1,2,1,2,2,2,2,2},{0,0,0=

,0,0,0,0,2,1,1,2,2,2,2,2,2},{0,0,0,0,0,0,0,1,2,2,2,2,2,2,2,1},{0,0,0,0,0,0,=

0,1,2,2,2,2,2,2,1,2},{0,0,0,0,0,0,0,1,2,2,2,2,2,1,2,2},{0,0,0,0,0,0,0,1,2,2=

,2,2,1,2,2,2},{0,0,0,0,0,0,0,1,2,2,2,1,2,2,2,2},{0,0,0,0,0,0,0,1,2,2,1,2,2,=

2,2,2},{0,0,0,0,0,0,0,1,2,1,2,2,2,2,2,2},{0,0,0,0,0,0,0,1,1,2,2,2,2,2,2,2}}=

}

I don't have a proof with this method that none of the above is really

equivalent to others, alas. I have

some (empirical) evidence that this may be the case, but obviously better

arguments are needed here.

That is, if I understood the problem at all correctly, which is something =

I

am also starting to doubt.

Regards,

Leonid

On Mon, Dec 20, 2010 at 8:40 AM, Yaroslav Bulatov <yaros...@gmail.com>wro=

te:

Dec 25, 2010, 2:34:43 AM12/25/10

to

Yaroslav,

This wasn't easy, but this was fun. You better have M8 and some C compiler

installed

(if you lack the latter, remove the CompilationTarget->"C" option where it

is present, but

the code will be 2-4 times slower).

(* These first 2 functions courtesy of Norbert Pozar *)

Clear[extractPositionFromSparseArray, positionExtr];

extractPositionFromSparseArray[

HoldPattern[SparseArray[u___]]] := {u}[[4, 2, 2]]

positionExtr[x_List, n_] :=

Flatten@extractPositionFromSparseArray[

SparseArray[Unitize[x - n], Automatic, 1]]

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}];

Clear[permuteNTimes ];

permuteNTimes =

Compile[{{list, _Integer, 2}, {permutation, _Integer,

1}, {n, _Integer}},

Nest[#[[All, permutation]] &, list, n]];

Clear[reducePermutationsC];

reducePermutationsC =

Compile[{{list, _Integer, 2}, {permutation, _Integer, 1}},

Module[{perms = list, n = 1, pos = Table[-1, {Length[list]}],

permuted = list, psame = Table[{-1, -1}, {Length[list]}],

done = False},

While[! done,

permuted = permuteNTimes[perms , permutation, n++];

psame = positionsOfSame[perms , permuted];

If[psame != {{0}},

pos =

positionExtr[UnitStep[(#[[2]] - #[[1]]) &@Transpose[psame]], 0];

If[pos != {},

perms = Delete[perms, List /@ pos],

(* else *)

done = True

],

(* else *)

done = True

]];

perms

], {{positionExtr[_, _], _Integer, 1}},

CompilationOptions -> {"InlineCompiledFunctions" -> True,

"InlineExternalDefinitions" -> True}, CompilationTarget -> "C"];

Clear[getUniquePermutatations];

getUniquePermutatations[permutations_, actions : {__Cycles}] /;

Length[permutations] > 0 :=

With[{rlen = Range[Length[permutations[[1]]]]},

Fold[reducePermutationsC, permutations,

Map[Permute[rlen, #] &, actions]]]

Now, let's see. The warm-up:

In[17]:= getUniquePermutatations[Permutations[Range[8]],

DeleteCases[GroupElements@DihedralGroup[8],Cycles[{}]]]//Short//Timing

Out[17]=

{0.25,{{1,2,3,5,6,8,7,4},{1,2,3,5,7,4,6,8},{1,2,3,5,7,4,8,6},<<3255>>,{8,6,5,4,3,1,2,7},{8,6,5,4,3,2,1,7}}}

The real deal:

In[18]:= (perms = Permutations[{2,2,2,2,2,2,2,1,1,0,0,0,0,0,0,0}])//Length

Out[18]= 411840

In[26]:= (un = getUniquePermutatations[perms,

DeleteCases[GroupElements@DihedralGroup[16],Cycles[{}]]])//Short//Timing

Out[26]=

{5.469,{{0,0,0,2,0,2,2,0,0,2,2,2,2,1,1,0},<<15183>>,{0,0,0,0,0,0,0,1,1,2,2,2,2,2,2,2}}}

So, about 6 seconds on my machine. If need faster, could probably

parallelize, I did not bother. Hope that

the code is correct. Some simple checks like say this one:

In[27]:= pm =

Permute[#,Cycles[{{2,16},{3,15},{4,14},{5,13},{6,12},{7,11},{8,10}}]]&/@un;

In[28]:= positionsOfSame[un,pm]

Out[28]=

{{11886,11886},{9790,9790},{8547,8547},{7889,7889},{5970,5970},{4821,4821},{4142,4142},{3798,3798},{1333,1333},{105,105}}

suggest that it is: the only elements that are the same in the result and

result permuted according to

some group element are those that are invariant under the transform (have

same positions).

I'd say, it is a good showcase for a new Compile.

Regards and Happy Holidays,

Leonid

Dec 26, 2010, 4:01:17 AM12/26/10

to

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 *)

(* 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,

Dec 28, 2010, 8:44:41 AM12/28/10

to

Conceptually, the code that you want is Length@GroupOrbits[group,

Permutations@lst, Permute]. This isn't going to be very efficient

though (I believe it will be improved in future versions of

Mathematica). Here's a brute force approach:

In[1]:= Module[{gens, verts, edges},

gens = PermutationList /@ GroupGenerators@DihedralGroup[16];

verts =

Permutations@{2, 2, 2, 2, 2, 2, 2, 1, 1, 0, 0, 0, 0, 0, 0, 0};

edges = Join @@ (Transpose@{verts, verts[[All, #]]} &) /@ gens;

Length@ConnectedComponents@Graph[Rule @@@ Union@edges]] // Timing

Out[1]= {10.578, 12940}

Maxim Rytin

m...@inbox.ru

Reply all

Reply to author

Forward

0 new messages

Search

Clear search

Close search

Google apps

Main menu