permlib  0.2.9
Library for permutation computations
type_recognition.h
1 // ---------------------------------------------------------------------------
2 //
3 // This file is part of PermLib.
4 //
5 // Copyright (c) 2009-2011 Thomas Rehn <thomas@carmen76.de>
6 // All rights reserved.
7 //
8 // Redistribution and use in source and binary forms, with or without
9 // modification, are permitted provided that the following conditions
10 // are met:
11 // 1. Redistributions of source code must retain the above copyright
12 // notice, this list of conditions and the following disclaimer.
13 // 2. Redistributions in binary form must reproduce the above copyright
14 // notice, this list of conditions and the following disclaimer in the
15 // documentation and/or other materials provided with the distribution.
16 // 3. The name of the author may not be used to endorse or promote products
17 // derived from this software without specific prior written permission.
18 //
19 // THIS SOFTWARE IS PROVIDED BY THE AUTHOR ``AS IS'' AND ANY EXPRESS OR
20 // IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES
21 // OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED.
22 // IN NO EVENT SHALL THE AUTHOR BE LIABLE FOR ANY DIRECT, INDIRECT,
23 // INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT
24 // NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
25 // DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
26 // THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
27 // (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF
28 // THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
29 //
30 // ---------------------------------------------------------------------------
31 
32 #ifndef TYPE_RECOGNITION_H_
33 #define TYPE_RECOGNITION_H_
34 
35 #include <permlib/prime_helper.h>
36 #include <permlib/transversal/schreier_tree_transversal.h>
37 #include <permlib/generator/bsgs_generator.h>
38 #include <permlib/construct/schreier_sims_construction.h>
39 #include <permlib/transversal/orbit_set.h>
40 #include <permlib/test/giant_test.h>
41 #include <permlib/test/group_type.h>
42 #include <permlib/test/primitivity_sgs_test.h>
43 #include <permlib/test/primitivity_test.h>
44 #include <permlib/permlib_api.h>
45 
46 #include <boost/shared_ptr.hpp>
47 #include <boost/math/common_factor_rt.hpp>
48 #include <iostream>
49 
50 
51 namespace permlib {
52 
54 
58 template<class PERM, class TRANSVERSAL = SchreierTreeTransversal<PERM> >
60 public:
62  typedef boost::shared_ptr<BSGS<PERM, TRANSVERSAL> > PermutationGroupPtr;
63 
69  template<class InputIterator>
70  TypeRecognition(unsigned int n, InputIterator genBegin, InputIterator genEnd);
71 
76 
78  GroupType* analyze() { return analyze(true); }
80  PermutationGroupPtr bsgs() const { return m_bsgs; }
81 
83 
89  GroupType* largeSymmetricDiagonalSubgroup(std::vector<dom_int>& orbitCharacteristic) const;
90 private:
91  const unsigned int m_n;
92  std::list<typename PERM::ptr> m_generators;
93  const PermutationGroupPtr m_externalBSGS;
94  PermutationGroupPtr m_bsgs;
95 
99  GroupType* analyze(bool allowRecursiveCalls);
100 
102 
105  GroupType* checkOrbitActions();
106 
108  static const unsigned int m_bsgsComputationDegreeLimit = 50;
109 };
110 
111 template<class PERM, class TRANSVERSAL>
112 template<class InputIterator>
113 TypeRecognition<PERM,TRANSVERSAL>::TypeRecognition(unsigned int n, InputIterator genBegin, InputIterator genEnd)
114  : m_n(n), m_generators(genBegin, genEnd)
115 { }
116 
117 template<class PERM, class TRANSVERSAL>
119  : m_n(bsgs_->n), m_generators(bsgs_->S.begin(), bsgs_->S.end()), m_externalBSGS(bsgs_), m_bsgs(m_externalBSGS)
120 { }
121 
122 
123 template<class PERM, class TRANSVERSAL>
124 GroupType* TypeRecognition<PERM,TRANSVERSAL>::analyze(bool allowRecursiveCalls) {
125  if (m_n < m_bsgsComputationDegreeLimit && ! m_bsgs) {
127  m_bsgs = PermutationGroupPtr(new BSGS<PERM, TRANSVERSAL>(ssc.construct(m_generators.begin(), m_generators.end())));
128  }
129 
130  // error probability for randomized algorithms
131  const double eps = 1e-3;
132 
133  if (m_bsgs) {
134  const unsigned int groupIntOrder = m_bsgs->template order<unsigned int>();
135 
136  if (groupIntOrder == 1)
137  return new TrivialGroupType(m_n);
138 
139  if (groupIntOrder == 2)
140  return new SymmetricGroupType(groupIntOrder, m_n);
141 
142  //
143  // check for cyclic groups
144  //
145  if (PrimeHelper::isPrimeNumber(groupIntOrder, false)) {
146  BOOST_ASSERT( m_n % groupIntOrder == 0 );
147  return new CyclicGroupType(groupIntOrder, m_n);
148  }
149 
150  if (groupIntOrder == 4) {
151  // distinguish between C_4 and the Klein four group
152  BSGSGenerator<TRANSVERSAL> bsgsGen(m_bsgs->U);
153  while (bsgsGen.hasNext()) {
154  PERM el = bsgsGen.next();
155  if (el.order() == 4)
156  return new CyclicGroupType(4, m_n);
157  }
158  return new DirectProductGroupType(new SymmetricGroupType(2,2), new SymmetricGroupType(2,2), m_n);
159  }
160 
161  //
162  // check for symmetric groups, natural action
163  //
164  unsigned int symmetricGroupCandidateDegree = 0;
165 
166  std::vector<unsigned int> transversalSizes(m_bsgs->U.size());
167  for (unsigned int i = 0; i < m_bsgs->U.size(); ++i) {
168  transversalSizes[i] = m_bsgs->U[i].size();
169  }
170  std::sort(transversalSizes.begin(), transversalSizes.end());
171 
172  // test whether group order is a factorial of some natural number
173 
174  unsigned int expectedFactor = 2;
175  for (std::vector<unsigned int>::const_iterator it = transversalSizes.begin(); it != transversalSizes.end(); ++it) {
176  if (*it == 1)
177  continue;
178  if (*it == expectedFactor)
179  ++expectedFactor;
180  else {
181  expectedFactor = 0;
182  break;
183  }
184  }
185 
186  if (expectedFactor > 0)
187  symmetricGroupCandidateDegree = expectedFactor - 1;
188 
189  if (symmetricGroupCandidateDegree == m_n)
190  return new SymmetricGroupType(symmetricGroupCandidateDegree, m_n);
191 
192  // check for wreath products of symmetric groups if group is transitive
193  if (m_bsgs->U[0].size() == m_n) {
194  std::list<typename PERM::ptr> S1;
195  PointwiseStabilizerPredicate<PERM> stabPred(m_bsgs->B[0]);
196  BOOST_FOREACH(const typename PERM::ptr& p, m_bsgs->S) {
197  if (stabPred(p))
198  S1.push_back(p);
199  }
200  PrimitivitySGSTest<TRANSVERSAL> primitivityTest(m_bsgs->S.begin(), m_bsgs->S.end(), m_bsgs->B[0], S1.begin(), S1.end(), m_bsgs->U[0]);
201  std::vector<dom_int> minimalBlock;
202  bool groupIsPrimitive = primitivityTest.blockOfImprimitivity(&minimalBlock);
203  if ( ! groupIsPrimitive ) {
204  // TODO: avoid integer overflow in order computation
205  unsigned int degreeG = minimalBlock.size();
206  unsigned int degreeH = m_n / degreeG;
207  if (m_n % degreeG == 0) {
208  boost::uint64_t wreathSize = 1;
209  // group order must equal factorial(deg G)^(deg H) * factorial(deg H)
210  for (unsigned int i = 1; i <= degreeH; ++i) {
211  for (unsigned int j = 2; j <= degreeG; ++j) {
212  wreathSize *= j;
213  }
214  wreathSize *= i;
215  }
216  if (wreathSize == m_bsgs->order())
217  return new WreathSymmetricGroupType(degreeG, degreeH, m_n);
218  }
219  } else {
220  GiantTest<PERM> giantTest;
221  GiantTestBase::GiantGroupType giant = giantTest.determineGiantType(eps,
222  m_n, m_bsgs->S.begin(), m_bsgs->S.end(), true);
223  if (giant == GiantTestBase::Symmetric)
224  return new SymmetricGroupType(m_n, m_n);
225  else if (giant == GiantTestBase::Alternating)
226  return new AlternatingGroupType(m_n, m_n);
227  }
228  }
229 
230 
231  if (allowRecursiveCalls) {
232  // check for multiple linked copies of S_k
233  // TODO: check for inner direct products of S_k \times S_l
234  GroupType* t = checkOrbitActions();
235  if (t) {
236  return t;
237  }
238  }
239  } // end if m_bsgs
240  else // else if ! m_bsgs
241  {
242  GiantTest<PERM> giantTest;
243  GiantTestBase::GiantGroupType giant = giantTest.determineGiantType(eps,
244  m_n, m_generators.begin(), m_generators.end());
245  if (giant == GiantTestBase::Symmetric)
246  return new SymmetricGroupType(m_n, m_n);
247  else if (giant == GiantTestBase::Alternating)
248  return new AlternatingGroupType(m_n, m_n);
249 
250  if (allowRecursiveCalls) {
251  GroupType* t = checkOrbitActions();
252  if (t)
253  return t;
254  }
255  }
256 
257  if (m_bsgs)
258  return new AnonymousGroupType<>(m_n, m_bsgs->order());
259  return new AnonymousGroupType<>(m_n);
260 }
261 
262 template<class PERM, class TRANSVERSAL>
263 GroupType* TypeRecognition<PERM,TRANSVERSAL>::checkOrbitActions() {
264  if (m_generators.empty())
265  return NULL;
266 
267  typedef typename PERM::ptr PERMptr;
268  boost::dynamic_bitset<> worked(m_n);
269  GroupType* candidateType = 0;
270  unsigned int groupSupport = m_n;
271 
272  for (unsigned int i = 0; i < m_n; ++i) {
273  if (worked[i])
274  continue;
275  worked.set(i, true);
276  OrbitSet<PERM, dom_int> orbit;
277  orbit.orbit(i, m_generators, typename Transversal<PERM>::TrivialAction());
278  for (typename OrbitSet<PERM, dom_int>::const_iterator it = orbit.begin(); it != orbit.end(); ++it) {
279  worked.set(*it, true);
280  }
281  // ignore fix points
282  if (orbit.size() == 1) {
283  --groupSupport;
284  continue;
285  }
286 
287  // restrict group to this orbit
288  std::list<PERMptr> orbitGenerators;
289  BOOST_FOREACH(const PERMptr& p, m_generators) {
290  orbitGenerators.push_back(PERMptr(p->project(orbit.size(), orbit.begin(), orbit.end())));
291  }
292  TypeRecognition<PERM, TRANSVERSAL> orbitTypeRecognition(orbit.size(), orbitGenerators.begin(), orbitGenerators.end());
293  GroupType* orbitType = orbitTypeRecognition.analyze(false);
294  if ( ! candidateType )
295  candidateType = orbitType;
296  else {
297  const bool isEqualType = candidateType->equals(orbitType);
298  delete orbitType;
299  if ( ! isEqualType )
300  return NULL;
301  }
302  }
303  if (candidateType)
304  candidateType->setNonNaturalAction(groupSupport);
305  return candidateType;
306 }
307 
308 template<typename PERM>
310  std::vector<dom_int> operator()(const PERM& p, const std::vector<dom_int>& v) const {
311  std::vector<dom_int> ret(v.size());
312  for (unsigned int i = 0; i < v.size(); ++i)
313  ret[i] = p / v[i];
314  std::sort(ret.begin(), ret.end());
315  return ret;
316  }
317 };
318 
319 template<class PERM, class TRANSVERSAL>
320 GroupType* TypeRecognition<PERM,TRANSVERSAL>::largeSymmetricDiagonalSubgroup(std::vector<dom_int>& orbitCharacteristic) const {
321  if (m_generators.empty())
322  return NULL;
323 
324  typedef boost::shared_ptr<OrbitSet<PERM, dom_int> > OrbitPtr;
325  typedef typename PERM::ptr PERMptr;
326 
327  orbitCharacteristic.clear();
328  orbitCharacteristic.resize(m_n);
329  unsigned int orbitCharacteristicCounter = 0;
330 
331  std::list<OrbitPtr> orbits;
332  boost::dynamic_bitset<> worked(m_n);
333  for (unsigned int i = 0; i < m_n; ++i) {
334  if (worked[i])
335  continue;
336  worked.set(i, true);
337  OrbitPtr orbit(new OrbitSet<PERM, dom_int>());
338  orbit->orbit(i, m_generators, typename Transversal<PERM>::TrivialAction());
339  for (typename OrbitSet<PERM, dom_int>::const_iterator it = orbit->begin(); it != orbit->end(); ++it) {
340  worked.set(*it, true);
341  }
342  orbits.push_back(orbit);
343  }
344 
345  size_t orbitGCD = orbits.front()->size();
346  BOOST_FOREACH(const OrbitPtr& orbit, orbits) {
347  orbitGCD = boost::math::gcd(orbitGCD, orbit->size());
348  }
349 
350  GroupType* lastType = 0;
351  BOOST_FOREACH(const OrbitPtr& orbit, orbits) {
352  GroupType* orbitType = 0;
353  // restrict group to this orbit
354  std::list<PERMptr> orbitGenerators;
355  // copy orbit elements into vector because we need a fixed order to recover from 'project'ion
356  std::vector<dom_int> orbitV(orbit->begin(), orbit->end());
357  BOOST_FOREACH(const PERMptr& p, m_generators) {
358  orbitGenerators.push_back(PERMptr(p->project(orbit->size(), orbitV.begin(), orbitV.end())));
359  //TODO: filter identities
360  }
361  PermutationGroupPtr orbitGroup = construct(orbit->size(), orbitGenerators.begin(), orbitGenerators.end());
362 
363  // compute all minimal blocks
364  PrimitivityTest<PERM> primitivityTest(orbit->size(), orbitGenerators.begin(), orbitGenerators.end());
365  std::list<std::vector<dom_int> > blocks;
366  bool groupIsPrimitive = primitivityTest.allBlocks(&blocks);
367 
368  if (!groupIsPrimitive) {
369  BOOST_FOREACH(const std::vector<dom_int>& b, blocks) {
370  // if group is transitive (#orbits == 1), then we do not need the gcd condition
371  if (orbits.size() != 1 && b.size() % orbitGCD != 0)
372  continue;
373 
374  // look at group action within block b
375  // 1. compute stabilizer of b
376  // 2. project action of stabilizer onto b
377  // 3. check whether induced action is a symmetric group
378  PermutationGroupPtr stab = setStabilizer(*orbitGroup, b.begin(), b.end());
379  std::list<PERMptr> suborbitGenerators;
380  BOOST_FOREACH(const PERMptr& p, stab->S) {
381  suborbitGenerators.push_back(PERMptr(p->project(b.size(), b.begin(), b.end())));
382  }
383  TypeRecognition<PERM, TRANSVERSAL> suborbitRecognition(b.size(), suborbitGenerators.begin(), suborbitGenerators.end());
384  GroupType* subType = suborbitRecognition.analyze(false);
385  if (subType->type() == GroupType::Named && subType->isNaturalAction()) {
386  NamedGroupType* namedType = reinterpret_cast<NamedGroupType*>(subType);
387  if (strcmp(namedType->name(), "S") == 0) {
388  orbitType = subType;
389  }
390  }
391  if (orbitType) {
393  blockOrbit.orbit(b, orbitGenerators, BlockVectorAction<PERM>());
394  BOOST_FOREACH(const std::vector<dom_int>& block, std::make_pair(blockOrbit.begin(), blockOrbit.end())) {
395  BOOST_FOREACH(const dom_int j, block) {
396  orbitCharacteristic[orbitV[j]] = orbitCharacteristicCounter;
397  }
398  ++orbitCharacteristicCounter;
399  }
400  break;
401  } else {
402  delete subType;
403  }
404  }
405  } else {
406  TypeRecognition<PERM, TRANSVERSAL> suborbitRecognition(orbit->size(), orbitGenerators.begin(), orbitGenerators.end());
407  orbitType = suborbitRecognition.analyze(false);
408  BOOST_FOREACH(const dom_int& o, orbitV) {
409  orbitCharacteristic[o] = orbitCharacteristicCounter;
410  }
411  ++orbitCharacteristicCounter;
412  }
413  if (!orbitType) {
414  delete lastType;
415  return NULL;
416  }
417 
418  if (!lastType) {
419  lastType = orbitType;
420  } else {
421  if (!lastType->equals(orbitType)) {
422  delete orbitType;
423  delete lastType;
424  return NULL;
425  }
426  delete orbitType;
427  }
428  }
429 
430  if (lastType)
431  lastType->setNonNaturalAction(m_n);
432 
433  return lastType;
434 }
435 
436 }
437 
438 #endif
439 
action of a PERM on unsigned long element
Definition: transversal.h:112
const_iterator begin() const
begin-iterator to orbit elements
Definition: orbit_set.h:66
GroupType * analyze()
analyzes the given group and attempts to determine the group type
Definition: type_recognition.h:78
stores an orbit in a set for fast contains() operation
Definition: orbit_set.h:42
abstract base class for named groups (such as cyclic and symmetric groups)
Definition: group_type.h:157
const_iterator end() const
end-iterator to orbit elements
Definition: orbit_set.h:68
GroupType * largeSymmetricDiagonalSubgroup(std::vector< dom_int > &orbitCharacteristic) const
attempts to find a large diagonal subgroup of symmetric groups
Definition: type_recognition.h:320
const char * name() const
the name of the group
Definition: group_type.h:164
bool equals(const GroupType *type_) const
checks if two group types represent the same permutation group
Definition: group_type.h:72
static bool isPrimeNumber(unsigned int x, bool checkBounds)
checks if a given number is prime
Definition: prime_helper.h:52
bool allBlocks(std::list< std::vector< dom_int > > *blocks, bool findOnlyOneBlock=false) const
Definition: primitivity_test.h:132
boost::shared_ptr< BSGS< PERM, TRANSVERSAL > > PermutationGroupPtr
abbreviation for a pointer to a BSGS structure
Definition: type_recognition.h:62
Definition: type_recognition.h:309
GiantGroupType
Enumeration of "giant" groups, i.e. Alternating and Symmetric group.
Definition: giant_test.h:52
PermutationGroupPtr bsgs() const
returns a BSGS if one was constructed during the analysis
Definition: type_recognition.h:80
void setNonNaturalAction(unsigned int realDegree_)
stores the information that this group acts non-naturally on realDegree many elements ...
Definition: group_type.h:83
abstract base class for permutation group types
Definition: group_type.h:40
Class for a basic type recognition of permutation groups.
Definition: type_recognition.h:59
Represents a base and strong generating set (BSGS)
Definition: bsgs.h:58
void orbit(const PDOMAIN &beta, const std::list< typename PERM::ptr > &generators, Action a)
computes orbit of beta under generators
Definition: orbit_set.h:92
Tests a transitive group is availble for primitivity.
Definition: primitivity_test.h:55
TypeRecognition(unsigned int n, InputIterator genBegin, InputIterator genEnd)
Definition: type_recognition.h:113
BSGS construction with classic Schreier-Sims algorithm.
Definition: schreier_sims_construction.h:46
Definition: abstract_bsgs.h:49