permlib  0.2.9
Library for permutation computations
bsgs.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 
33 #ifndef BSGS_H_
34 #define BSGS_H_
35 
36 #include <map>
37 #include <list>
38 #include <vector>
39 
40 #include <boost/cstdint.hpp>
41 #include <boost/foreach.hpp>
42 #include <boost/scoped_ptr.hpp>
43 #include <boost/shared_ptr.hpp>
44 #include <boost/utility.hpp>
45 
46 #include <permlib/bsgs_core.h>
47 
48 #include <permlib/transversal/orbit_set.h>
49 #include <permlib/transversal/transversal.h>
50 #include <permlib/predicate/pointwise_stabilizer_predicate.h>
51 #include <permlib/predicate/stabilizes_point_predicate.h>
52 
53 #include <permlib/redundant_base_point_insertion_strategy.h>
54 
55 namespace permlib {
56 
57 template <class PERM, class TRANS>
58 struct BSGS;
59 
60 template <class PERM, class TRANS>
61 std::ostream &operator<< (std::ostream &out, const BSGS<PERM,TRANS> &bsgs) {
62  out << "BASE[" << bsgs.B.size() << "]" << std::endl;
63  BOOST_FOREACH(unsigned long beta, bsgs.B) {
64  out << static_cast<unsigned int>(beta+1) << ",";
65  }
66  out << std::endl;
67  out << "SGS[" << bsgs.S.size() << "]" << std::endl;
68  BOOST_FOREACH(const typename PERM::ptr &g, bsgs.S) {
69  out << *g << ",";
70  }
71  out << std::endl;
72  out << "U" << std::endl;
73  BOOST_FOREACH(const TRANS &U, bsgs.U) {
74  for (unsigned int i=0; i<bsgs.n; ++i)
75  // trigger transversal depth reload
76  boost::scoped_ptr<PERM> dummy(U.at(i));
77  out << U.size() << "{" << U.m_statMaxDepth << "}" << ",";
78  }
79  out << " = " << bsgs.order() << std::endl;
80  BOOST_FOREACH(const TRANS &U, bsgs.U) {
81  out << U << std::endl;
82  }
83  return out;
84 }
85 
87 template <class PERM, class TRANS>
88 struct BSGS : public BSGSCore<PERM,TRANS> {
89  typedef typename BSGSCore<PERM,TRANS>::PERMlist PERMlist;
90 
92  explicit BSGS(dom_int n);
94 
98  BSGS(const BSGS<PERM,TRANS>& bsgs);
99 
101 
106 
108 
111  template<typename Integer>
112  Integer order() const;
113 
115 
118  boost::uint64_t order() const;
119 
121 
126  unsigned int sift(const PERM& g, PERM& siftee, unsigned int j = 0) const;
128 
134  unsigned int sift(const PERM& g, PERM& siftee, unsigned int j, unsigned int k) const;
136  bool sifts(const PERM& g) const;
137 
139 
145  bool chooseBaseElement(const PERM &h, dom_int &beta) const;
147 
152  unsigned int insertRedundantBasePoint(unsigned int beta, unsigned int minPos = 0);
154  void stripRedundantBasePoints(int minPos = 0);
155 
157 
165 
167 
170  void orbit(unsigned int j, const PERMlist &generators);
172 
175  void orbitUpdate(unsigned int j, const PERMlist &generators, const typename PERM::ptr &g);
176 
178 
183  int insertGenerator(const typename PERM::ptr& g, bool updateOrbit);
185 
189  void updateOrbits(int pos);
190 
192 
195  PERM random(const int i = 0) const;
196 
198 
202  void conjugate(const PERM& g);
203 
205  friend std::ostream &operator<< <> (std::ostream &out, const BSGS<PERM,TRANS> &bsgs);
206 private:
208  template <class BaseIterator, class TransversalIterator>
209  unsigned int sift(const PERM& g, PERM& siftee, BaseIterator begin, BaseIterator end, TransversalIterator beginT, TransversalIterator endT) const;
210 
212  static int ms_bsgsId;
213 
215  void copyTransversals(const BSGS<PERM,TRANS>& bsgs);
216 };
217 
218 //
219 // ---- IMPLEMENTATION
220 //
221 
222 template <class PERM, class TRANS>
223 int BSGS<PERM,TRANS>::ms_bsgsId = 0;
224 
225 template <class PERM, class TRANS>
227  : BSGSCore<PERM,TRANS>(++ms_bsgsId, n_, 0)
228 {}
229 
230 template <class PERM, class TRANS>
232  : BSGSCore<PERM,TRANS>(bsgs.m_id, bsgs.B, bsgs.U, bsgs.n)
233 {
234  copyTransversals(bsgs);
235 }
236 
237 template <class PERM, class TRANS>
239  if (this == &bsgs)
240  return *this;
241 
242  this->B = bsgs.B;
243  this->n = bsgs.n;
244  this->m_id = bsgs.m_id;
245 
246  copyTransversals(bsgs);
247  return *this;
248 }
249 
250 template <class PERM, class TRANS>
251 template <class BaseIterator, class TransversalIterator>
252 unsigned int BSGS<PERM, TRANS>::sift(const PERM& g, PERM& siftee, BaseIterator begin, BaseIterator end, TransversalIterator beginT, TransversalIterator endT) const{
253  unsigned int k = 0;
254  siftee = g;
255  BaseIterator baseIt;
256  TransversalIterator transIt;
257  for (baseIt = begin, transIt = beginT; baseIt != end && transIt != endT; ++baseIt, ++transIt) {
258  unsigned long b = *baseIt;
259  const TRANS& U_i = *transIt;
260  //std::cout << " ~~~ sift " << siftee << " b" << b << std::endl;
261  boost::scoped_ptr<PERM> u_b(U_i.at(siftee / b));
262  if (u_b == 0)
263  return k;
264  u_b->invertInplace();
265  siftee *= *u_b;
266  ++k;
267  }
268  return k;
269 }
270 
271 template <class PERM, class TRANS>
272 unsigned int BSGS<PERM, TRANS>::sift(const PERM& g, PERM& siftee, unsigned int j) const {
273  return sift(g, siftee, this->B.begin() + j, this->B.end(), this->U.begin() + j, this->U.end());
274 }
275 
276 template <class PERM, class TRANS>
277 unsigned int BSGS<PERM, TRANS>::sift(const PERM& g, PERM& siftee, unsigned int j, unsigned int k) const {
278  return sift(g, siftee, this->B.begin() + j, this->B.begin() + k, this->U.begin() + j, this->U.begin() + k);
279 }
280 
281 template <class PERM, class TRANS>
282 bool BSGS<PERM, TRANS>::sifts(const PERM& g) const {
283  PERM siftee(this->n);
284  unsigned int m = sift(g, siftee);
285  return this->B.size() == m && siftee.isIdentity();
286 }
287 
288 template <class PERM, class TRANS>
289 bool BSGS<PERM, TRANS>::chooseBaseElement(const PERM &h, dom_int &beta) const {
290  for (beta = 0; beta < this->n; ++beta) {
291  if (std::find(this->B.begin(), this->B.end(), beta) != this->B.end())
292  continue;
293  if (h / beta != beta)
294  return true;
295  }
296  return false;
297 }
298 
299 template <class PERM, class TRANS>
300 void BSGS<PERM, TRANS>::orbit(unsigned int j, const PERMlist &generators) {
301  this->U[j].orbit(this->B[j], generators);
302 }
303 
304 template <class PERM, class TRANS>
305 void BSGS<PERM, TRANS>::orbitUpdate(unsigned int j, const PERMlist &generators, const typename PERM::ptr &g) {
306  this->U[j].orbitUpdate(this->B[j], generators, g);
307 }
308 
309 template <class PERM, class TRANS>
310 PERM BSGS<PERM, TRANS>::random(const int i) const {
311  BOOST_ASSERT( i >= 0 );
312  PERM g(this->n);
313  for (int l = this->U.size()-1; l>=i ; --l) {
314  //std::cout << l << " : " << U[l] << " : " << U[l].size() << std::endl;
315  unsigned long beta = *(boost::next(this->U[l].begin(), randomInt(this->U[l].size())));
316  boost::scoped_ptr<PERM> u_beta(this->U[l].at(beta));
317  g *= *u_beta;
318  }
319  return g;
320 }
321 
322 template <class PERM, class TRANS>
323 void BSGS<PERM, TRANS>::conjugate(const PERM& g) {
324  PERM gInv(g);
325  gInv.invertInplace();
326 
327  //
328  // to conjugate a BSGS, all three components (B,S,U) have to be adjusted
329  //
330 
331  // STEP 1: conjugate generating set S
332  BOOST_FOREACH(typename PERM::ptr& p, this->S) {
333  *p ^= gInv;
334  *p *= g;
335  }
336 
337  std::vector<dom_int> oldB(this->B);
338  for (unsigned int i = 0; i < this->U.size(); ++i) {
339  // STEP 2: adapt base B
340  this->B[i] = g / oldB[i];
341  // STEP 3: conjugate transversal U
342  this->U[i].permute(g, gInv);
343  }
344 }
345 
346 template <class PERM, class TRANS>
347 int BSGS<PERM, TRANS>::insertGenerator(const typename PERM::ptr& g, bool updateOrbit) {
348  int pos = 0;
349  for (; static_cast<unsigned int>(pos) < this->B.size(); ++pos) {
350  if (*g / this->B[pos] != this->B[pos])
351  break;
352  }
353 
354  if (static_cast<unsigned int>(pos) == this->B.size()) {
355  dom_int beta;
356  bool newBaseElement __attribute__((unused)) = chooseBaseElement(*g, beta);
357  BOOST_ASSERT( newBaseElement );
358  this->B.push_back(beta);
359  this->U.push_back(TRANS(this->n));
360  }
361 
362  const int insertionPosition = pos;
363  this->S.push_back(g);
364 
365  if (updateOrbit) {
366  bool groupOrderChanged = false;
367  for (; pos >= 0; --pos) {
368  PERMlist orbitGenerators;
369  const unsigned int oldTransversalSize = this->U[pos].size();
370  //std::cout << "INSERT orbit @ " << pos << " : " << g << std::endl;
371  std::copy_if(this->S.begin(), this->S.end(), std::back_inserter(orbitGenerators),
372  PointwiseStabilizerPredicate<PERM>(this->B.begin(), this->B.begin() + pos));
373  if (orbitGenerators.size() > 0) {
374  orbitUpdate(pos, orbitGenerators, g);
375 
376  // group order can only increase by adding generators
377  if (this->U[pos].size() > oldTransversalSize)
378  groupOrderChanged = true;
379  }
380  }
381 
382  if (!groupOrderChanged) {
383  this->S.pop_back();
384  return -1;
385  }
386  }
387 
388  return insertionPosition;
389 }
390 
391 template <class PERM, class TRANS>
393  if (pos < 0)
394  return;
395  for (; pos >= 0; --pos) {
396  PERMlist orbitGenerators;
397  std::copy_if(this->S.begin(), this->S.end(), std::back_inserter(orbitGenerators),
398  PointwiseStabilizerPredicate<PERM>(this->B.begin(), this->B.begin() + pos));
399  if (orbitGenerators.size() > 0)
400  orbit(pos, orbitGenerators);
401  }
402 }
403 
404 template <class PERM, class TRANS>
405 template <typename Integer>
406 Integer BSGS<PERM, TRANS>::order() const {
407  Integer orderValue(1);
408  BOOST_FOREACH(const TRANS &Ui, this->U) {
409  orderValue *= Ui.size();
410  }
411  return orderValue;
412 }
413 
414 template <class PERM, class TRANS>
415 boost::uint64_t BSGS<PERM, TRANS>::order() const {
416  return order<boost::uint64_t>();
417 }
418 
419 
420 template <class PERM, class TRANS>
421 unsigned int BSGS<PERM, TRANS>::insertRedundantBasePoint(unsigned int beta, unsigned int minPos) {
422  PERMlist S_i;
424  int pos = is.findInsertionPoint(beta, S_i);
425  if (pos < 0)
426  return -(pos+1);
427  pos = std::max(static_cast<unsigned int>(pos), minPos);
428 
429  this->B.insert(this->B.begin() + pos, beta);
430  this->U.insert(this->U.begin() + pos, TRANS(this->n));
431  this->U[pos].orbit(beta, S_i);
432  return pos;
433 }
434 
435 template <class PERM, class TRANS>
437  for (int i = this->B.size()-1; i >= minPos; --i) {
438  if (this->U[i].size() <= 1) {
439  if (i == static_cast<int>(this->B.size()-1)) {
440  this->B.pop_back();
441  this->U.pop_back();
442  } else {
443  this->B.erase(this->B.begin() + i);
444  this->U.erase(this->U.begin() + i);
445  }
446  }
447  }
448 }
449 
450 
452 
456 template <class PERM>
457 class StrongGeneratingSetSorter : public std::binary_function<typename PERM::ptr, typename PERM::ptr, bool> {
458 public:
463  template<class InputIterator>
464  StrongGeneratingSetSorter(InputIterator baseBegin, InputIterator baseEnd) : m_base(baseBegin, baseEnd) { }
465 
467  bool operator()(const typename PERM::ptr& p1, const typename PERM::ptr& p2) const {
468  BOOST_FOREACH(const dom_int b, m_base) {
469  if ( p1->at(b) == b && p2->at(b) != b )
470  return true;
471  if ( p1->at(b) != b )
472  return false;
473  }
474  return false;
475  }
476 private:
477  std::vector<dom_int> m_base;
478 };
479 
480 template <class PERM, class TRANS>
482  PERMlist sortedSGS(this->S);
483  sortedSGS.sort(StrongGeneratingSetSorter<PERM>(this->B.begin(), this->B.end()));
484 
485  PERMlist filteredSGS;
487  dom_int oldBaseElement = static_cast<dom_int>(-1);
488  BOOST_FOREACH(const typename PERM::ptr& gen, sortedSGS) {
489  if (gen->isIdentity())
490  continue;
491  filteredSGS.push_back(gen);
492 
493  // Compute to which base element this strong generator belongs.
494  // That is, gen stabilizes all base points up to baseElement.
495  // No generator can stabilize all base elements because we excluded
496  // identity permutations before.
497  dom_int baseElement = this->B.front();
498  BOOST_FOREACH(const dom_int b, this->B) {
499  baseElement = b;
500  if (*gen / b != b)
501  break;
502  }
503  PERMLIB_DEBUG(std::cout << "gen " << *gen << " @ " << baseElement << std::endl;)
504 
506  newOrbit->orbit(baseElement, filteredSGS, typename Transversal<PERM>::TrivialAction());
507  if (oldBaseElement == baseElement && newOrbit->size() == oldOrbit->size()) {
508  delete newOrbit;
509  PERMLIB_DEBUG(std::cout << " removed\n";)
510  filteredSGS.pop_back();
511  } else {
512  delete oldOrbit;
513  oldOrbit = newOrbit;
514  }
515  oldBaseElement = baseElement;
516  }
517  delete oldOrbit;
518 
519  this->S = filteredSGS;
520 }
521 
522 template <class PERM, class TRANS>
524  std::map<PERM*,typename PERM::ptr> genMap;
525  BOOST_FOREACH(const typename PERM::ptr& p, bsgs.S) {
526  typename PERM::ptr deepcopy(new PERM(*p));
527  //std::cout << "found " << p.get() << " = " << *p << std::endl;
528  genMap.insert(std::make_pair(p.get(), deepcopy));
529  this->S.push_back(deepcopy);
530  }
531 
532  BOOST_ASSERT(this->B.size() == bsgs.B.size());
533  BOOST_ASSERT(bsgs.B.size() == bsgs.U.size());
534  this->U.clear();
535  this->U.resize(bsgs.U.size(), TRANS(bsgs.n));
536  BOOST_ASSERT(this->U.size() == bsgs.U.size());
537 
538  for (unsigned int i=0; i<this->U.size(); ++i) {
539  this->U[i] = bsgs.U[i].clone(genMap);
540  }
541 }
542 
543 }
544 
545 #endif // -- BSGS_H_
dom_int n
degree of group
Definition: bsgs_core.h:61
action of a PERM on unsigned long element
Definition: transversal.h:112
unsigned int insertRedundantBasePoint(unsigned int beta, unsigned int minPos=0)
inserts a redundant base beta
Definition: bsgs.h:421
BSGS< PERM, TRANS > & operator=(const BSGS< PERM, TRANS > &)
assignment operator
Definition: bsgs.h:238
insertion position after first non-trivial transversal
Definition: redundant_base_point_insertion_strategy.h:68
Class that can be used to sort a strong generating set.
Definition: bsgs.h:457
bool sifts(const PERM &g) const
true iff g sifts through transversal system
Definition: bsgs.h:282
void orbitUpdate(unsigned int j, const PERMlist &generators, const typename PERM::ptr &g)
updates the j-th fundamental orbit with the given orbit generators and a new generator g ...
Definition: bsgs.h:305
stores an orbit in a set for fast contains() operation
Definition: orbit_set.h:42
BSGS(dom_int n)
constructs an empty group of degree n
Definition: bsgs.h:226
bool chooseBaseElement(const PERM &h, dom_int &beta) const
tries to find a new base element
Definition: bsgs.h:289
core data of a base and strong generating set (BSGS)
Definition: bsgs_core.h:42
predicate matching a permutation if it stabilizes a given list of points pointwise ...
Definition: pointwise_stabilizer_predicate.h:42
Integer order() const
order of the group
Definition: bsgs.h:406
StrongGeneratingSetSorter(InputIterator baseBegin, InputIterator baseEnd)
Definition: bsgs.h:464
unsigned int sift(const PERM &g, PERM &siftee, unsigned int j=0) const
sifts an element through the specified transversal range
Definition: bsgs.h:272
int m_id
id of this BSGS instance
Definition: bsgs_core.h:81
void conjugate(const PERM &g)
conjugate group with a permutation
Definition: bsgs.h:323
void stripRedundantBasePoints(int minPos=0)
strips redundant base points from the end to the minPos-th base element
Definition: bsgs.h:436
std::vector< dom_int > B
base
Definition: bsgs_core.h:55
size_t size() const
number of orbit elements
Definition: orbit_set.h:60
PERM random(const int i=0) const
generates a uniformly distributed random element of
Definition: bsgs.h:310
int insertGenerator(const typename PERM::ptr &g, bool updateOrbit)
adds a new group generator
Definition: bsgs.h:347
bool operator()(const typename PERM::ptr &p1, const typename PERM::ptr &p2) const
true iff p1 stabilizes more base points (in increasing order) than p2
Definition: bsgs.h:467
Represents a base and strong generating set (BSGS)
Definition: bsgs.h:58
virtual int findInsertionPoint(dom_int beta, std::list< typename PERM::ptr > &S_i) const
finds possible insertion point for base point
Definition: redundant_base_point_insertion_strategy.h:73
std::vector< TRANS > U
transversals along the stabilizer chain
Definition: bsgs_core.h:59
void orbit(unsigned int j, const PERMlist &generators)
re-computes the j-th fundamental orbit with the given orbit generators
Definition: bsgs.h:300
void stripRedundantStrongGenerators()
removes redundant generators
Definition: bsgs.h:481
void updateOrbits(int pos)
updates transversals/orbits
Definition: bsgs.h:392
PERMlist S
strong generating set
Definition: bsgs_core.h:57
Definition: abstract_bsgs.h:49