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