permlib
0.2.6
Library for permutation computations
|
00001 // --------------------------------------------------------------------------- 00002 // 00003 // This file is part of PermLib. 00004 // 00005 // Copyright (c) 2009-2011 Thomas Rehn <thomas@carmen76.de> 00006 // All rights reserved. 00007 // 00008 // Redistribution and use in source and binary forms, with or without 00009 // modification, are permitted provided that the following conditions 00010 // are met: 00011 // 1. Redistributions of source code must retain the above copyright 00012 // notice, this list of conditions and the following disclaimer. 00013 // 2. Redistributions in binary form must reproduce the above copyright 00014 // notice, this list of conditions and the following disclaimer in the 00015 // documentation and/or other materials provided with the distribution. 00016 // 3. The name of the author may not be used to endorse or promote products 00017 // derived from this software without specific prior written permission. 00018 // 00019 // THIS SOFTWARE IS PROVIDED BY THE AUTHOR ``AS IS'' AND ANY EXPRESS OR 00020 // IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES 00021 // OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. 00022 // IN NO EVENT SHALL THE AUTHOR BE LIABLE FOR ANY DIRECT, INDIRECT, 00023 // INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT 00024 // NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, 00025 // DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY 00026 // THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT 00027 // (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF 00028 // THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. 00029 // 00030 // --------------------------------------------------------------------------- 00031 00032 00033 #ifndef BSGS_H_ 00034 #define BSGS_H_ 00035 00036 #include <map> 00037 #include <list> 00038 #include <vector> 00039 00040 #include <boost/cstdint.hpp> 00041 #include <boost/foreach.hpp> 00042 #include <boost/scoped_ptr.hpp> 00043 #include <boost/shared_ptr.hpp> 00044 #include <boost/utility.hpp> 00045 00046 #include <permlib/bsgs_core.h> 00047 00048 #include <permlib/predicate/pointwise_stabilizer_predicate.h> 00049 #include <permlib/predicate/stabilizes_point_predicate.h> 00050 00051 #include <permlib/redundant_base_point_insertion_strategy.h> 00052 00053 namespace permlib { 00054 00055 template <class PERM, class TRANS> 00056 struct BSGS; 00057 00058 template <class PERM, class TRANS> 00059 std::ostream &operator<< (std::ostream &out, const BSGS<PERM,TRANS> &bsgs) { 00060 out << "BASE[" << bsgs.B.size() << "]" << std::endl; 00061 BOOST_FOREACH(unsigned long beta, bsgs.B) { 00062 out << static_cast<unsigned int>(beta+1) << ","; 00063 } 00064 out << std::endl; 00065 out << "SGS[" << bsgs.S.size() << "]" << std::endl; 00066 BOOST_FOREACH(const typename PERM::ptr &g, bsgs.S) { 00067 out << *g << ","; 00068 } 00069 out << std::endl; 00070 out << "U" << std::endl; 00071 BOOST_FOREACH(const TRANS &U, bsgs.U) { 00072 for (unsigned int i=0; i<bsgs.n; ++i) 00073 // trigger transversal depth reload 00074 boost::scoped_ptr<PERM> dummy(U.at(i)); 00075 out << U.size() << "{" << U.m_statMaxDepth << "}" << ","; 00076 } 00077 out << " = " << bsgs.order() << std::endl; 00078 BOOST_FOREACH(const TRANS &U, bsgs.U) { 00079 out << U << std::endl; 00080 } 00081 return out; 00082 } 00083 00085 template <class PERM, class TRANS> 00086 struct BSGS : public BSGSCore<PERM,TRANS> { 00087 typedef typename BSGSCore<PERM,TRANS>::PERMlist PERMlist; 00088 00090 explicit BSGS(dom_int n); 00092 00096 BSGS(const BSGS<PERM,TRANS>& bsgs); 00097 00099 00103 BSGS<PERM,TRANS>& operator=(const BSGS<PERM,TRANS>&); 00104 00106 00109 template<typename Integer> 00110 Integer order() const; 00111 00113 00116 boost::uint64_t order() const; 00117 00119 00124 unsigned int sift(const PERM& g, PERM& siftee, unsigned int j = 0) const; 00126 00132 unsigned int sift(const PERM& g, PERM& siftee, unsigned int j, unsigned int k) const; 00134 bool sifts(const PERM& g) const; 00135 00137 00143 bool chooseBaseElement(const PERM &h, dom_int &beta) const; 00145 00150 unsigned int insertRedundantBasePoint(unsigned int beta, unsigned int minPos = 0); 00152 void stripRedundantBasePoints(int minPos = 0); 00153 00155 00158 void orbit(unsigned int j, const PERMlist &generators); 00160 00163 void orbitUpdate(unsigned int j, const PERMlist &generators, const typename PERM::ptr &g); 00164 00166 00171 int insertGenerator(const typename PERM::ptr& g, bool updateOrbit); 00173 00177 void updateOrbits(int pos); 00178 00180 00183 PERM random(const int i = 0) const; 00184 00186 00190 void conjugate(const PERM& g); 00191 00193 friend std::ostream &operator<< <> (std::ostream &out, const BSGS<PERM,TRANS> &bsgs); 00194 private: 00196 template <class BaseIterator, class TransversalIterator> 00197 unsigned int sift(const PERM& g, PERM& siftee, BaseIterator begin, BaseIterator end, TransversalIterator beginT, TransversalIterator endT) const; 00198 00200 static int ms_bsgsId; 00201 00203 void copyTransversals(const BSGS<PERM,TRANS>& bsgs); 00204 }; 00205 00206 // 00207 // ---- IMPLEMENTATION 00208 // 00209 00210 template <class PERM, class TRANS> 00211 int BSGS<PERM,TRANS>::ms_bsgsId = 0; 00212 00213 template <class PERM, class TRANS> 00214 BSGS<PERM,TRANS>::BSGS(dom_int n_) 00215 : BSGSCore<PERM,TRANS>(++ms_bsgsId, n_, 0) 00216 {} 00217 00218 template <class PERM, class TRANS> 00219 BSGS<PERM,TRANS>::BSGS(const BSGS<PERM,TRANS>& bsgs) 00220 : BSGSCore<PERM,TRANS>(bsgs.m_id, bsgs.B, bsgs.U, bsgs.n) 00221 { 00222 copyTransversals(bsgs); 00223 } 00224 00225 template <class PERM, class TRANS> 00226 BSGS<PERM,TRANS>& BSGS<PERM,TRANS>::operator=(const BSGS<PERM,TRANS>& bsgs) { 00227 if (this == &bsgs) 00228 return *this; 00229 00230 this->B = bsgs.B; 00231 this->n = bsgs.n; 00232 this->m_id = bsgs.m_id; 00233 00234 copyTransversals(bsgs); 00235 return *this; 00236 } 00237 00238 template <class PERM, class TRANS> 00239 template <class BaseIterator, class TransversalIterator> 00240 unsigned int BSGS<PERM, TRANS>::sift(const PERM& g, PERM& siftee, BaseIterator begin, BaseIterator end, TransversalIterator beginT, TransversalIterator endT) const{ 00241 unsigned int k = 0; 00242 siftee = g; 00243 BaseIterator baseIt; 00244 TransversalIterator transIt; 00245 for (baseIt = begin, transIt = beginT; baseIt != end && transIt != endT; ++baseIt, ++transIt) { 00246 unsigned long b = *baseIt; 00247 const TRANS& U_i = *transIt; 00248 //std::cout << " ~~~ sift " << siftee << " b" << b << std::endl; 00249 boost::scoped_ptr<PERM> u_b(U_i.at(siftee / b)); 00250 if (u_b == 0) 00251 return k; 00252 u_b->invertInplace(); 00253 siftee *= *u_b; 00254 ++k; 00255 } 00256 return k; 00257 } 00258 00259 template <class PERM, class TRANS> 00260 unsigned int BSGS<PERM, TRANS>::sift(const PERM& g, PERM& siftee, unsigned int j) const { 00261 return sift(g, siftee, this->B.begin() + j, this->B.end(), this->U.begin() + j, this->U.end()); 00262 } 00263 00264 template <class PERM, class TRANS> 00265 unsigned int BSGS<PERM, TRANS>::sift(const PERM& g, PERM& siftee, unsigned int j, unsigned int k) const { 00266 return sift(g, siftee, this->B.begin() + j, this->B.begin() + k, this->U.begin() + j, this->U.begin() + k); 00267 } 00268 00269 template <class PERM, class TRANS> 00270 bool BSGS<PERM, TRANS>::sifts(const PERM& g) const { 00271 PERM siftee(this->n); 00272 unsigned int m = sift(g, siftee); 00273 return this->B.size() == m && siftee.isIdentity(); 00274 } 00275 00276 template <class PERM, class TRANS> 00277 bool BSGS<PERM, TRANS>::chooseBaseElement(const PERM &h, dom_int &beta) const { 00278 for (beta = 0; beta < this->n; ++beta) { 00279 if (std::find(this->B.begin(), this->B.end(), beta) != this->B.end()) 00280 continue; 00281 if (h / beta != beta) 00282 return true; 00283 } 00284 return false; 00285 } 00286 00287 template <class PERM, class TRANS> 00288 void BSGS<PERM, TRANS>::orbit(unsigned int j, const PERMlist &generators) { 00289 this->U[j].orbit(this->B[j], generators); 00290 } 00291 00292 template <class PERM, class TRANS> 00293 void BSGS<PERM, TRANS>::orbitUpdate(unsigned int j, const PERMlist &generators, const typename PERM::ptr &g) { 00294 this->U[j].orbitUpdate(this->B[j], generators, g); 00295 } 00296 00297 template <class PERM, class TRANS> 00298 PERM BSGS<PERM, TRANS>::random(const int i) const { 00299 BOOST_ASSERT( i >= 0 ); 00300 PERM g(this->n); 00301 for (int l = this->U.size()-1; l>=i ; --l) { 00302 //std::cout << l << " : " << U[l] << " : " << U[l].size() << std::endl; 00303 unsigned long beta = *(boost::next(this->U[l].begin(), randomInt(this->U[l].size()))); 00304 boost::scoped_ptr<PERM> u_beta(this->U[l].at(beta)); 00305 g *= *u_beta; 00306 } 00307 return g; 00308 } 00309 00310 template <class PERM, class TRANS> 00311 void BSGS<PERM, TRANS>::conjugate(const PERM& g) { 00312 PERM gInv(g); 00313 gInv.invertInplace(); 00314 00315 // 00316 // to conjugate a BSGS, all three components (B,S,U) have to be adjusted 00317 // 00318 00319 // STEP 1: conjugate generating set S 00320 BOOST_FOREACH(typename PERM::ptr& p, this->S) { 00321 *p ^= gInv; 00322 *p *= g; 00323 } 00324 00325 std::vector<dom_int> oldB(this->B); 00326 for (unsigned int i = 0; i < this->U.size(); ++i) { 00327 // STEP 2: adapt base B 00328 this->B[i] = g / oldB[i]; 00329 // STEP 3: conjugate transversal U 00330 this->U[i].permute(g, gInv); 00331 } 00332 } 00333 00334 template <class PERM, class TRANS> 00335 int BSGS<PERM, TRANS>::insertGenerator(const typename PERM::ptr& g, bool updateOrbit) { 00336 int pos = 0; 00337 for (; static_cast<unsigned int>(pos) < this->B.size(); ++pos) { 00338 if (*g / this->B[pos] != this->B[pos]) 00339 break; 00340 } 00341 00342 if (static_cast<unsigned int>(pos) == this->B.size()) { 00343 dom_int beta; 00344 bool newBaseElement __attribute__((unused)) = chooseBaseElement(*g, beta); 00345 BOOST_ASSERT( newBaseElement ); 00346 this->B.push_back(beta); 00347 this->U.push_back(TRANS(this->n)); 00348 } 00349 00350 const int insertionPosition = pos; 00351 this->S.push_back(g); 00352 00353 if (updateOrbit) { 00354 bool groupOrderChanged = false; 00355 for (; pos >= 0; --pos) { 00356 PERMlist orbitGenerators; 00357 const unsigned int oldTransversalSize = this->U[pos].size(); 00358 //std::cout << "INSERT orbit @ " << pos << " : " << g << std::endl; 00359 std::copy_if(this->S.begin(), this->S.end(), std::back_inserter(orbitGenerators), 00360 PointwiseStabilizerPredicate<PERM>(this->B.begin(), this->B.begin() + pos)); 00361 if (orbitGenerators.size() > 0) { 00362 orbitUpdate(pos, orbitGenerators, g); 00363 00364 // group order can only increase by adding generators 00365 if (this->U[pos].size() > oldTransversalSize) 00366 groupOrderChanged = true; 00367 } 00368 } 00369 00370 if (!groupOrderChanged) { 00371 this->S.pop_back(); 00372 return -1; 00373 } 00374 } 00375 00376 return insertionPosition; 00377 } 00378 00379 template <class PERM, class TRANS> 00380 void BSGS<PERM, TRANS>::updateOrbits(int pos) { 00381 if (pos < 0) 00382 return; 00383 for (; pos >= 0; --pos) { 00384 PERMlist orbitGenerators; 00385 std::copy_if(this->S.begin(), this->S.end(), std::back_inserter(orbitGenerators), 00386 PointwiseStabilizerPredicate<PERM>(this->B.begin(), this->B.begin() + pos)); 00387 if (orbitGenerators.size() > 0) 00388 orbit(pos, orbitGenerators); 00389 } 00390 } 00391 00392 template <class PERM, class TRANS> 00393 template <typename Integer> 00394 Integer BSGS<PERM, TRANS>::order() const { 00395 Integer orderValue(1); 00396 BOOST_FOREACH(const TRANS &Ui, this->U) { 00397 orderValue *= Ui.size(); 00398 } 00399 return orderValue; 00400 } 00401 00402 template <class PERM, class TRANS> 00403 boost::uint64_t BSGS<PERM, TRANS>::order() const { 00404 return order<boost::uint64_t>(); 00405 } 00406 00407 00408 template <class PERM, class TRANS> 00409 unsigned int BSGS<PERM, TRANS>::insertRedundantBasePoint(unsigned int beta, unsigned int minPos) { 00410 PERMlist S_i; 00411 TrivialRedundantBasePointInsertionStrategy<PERM,TRANS> is(*this); 00412 int pos = is.findInsertionPoint(beta, S_i); 00413 if (pos < 0) 00414 return -(pos+1); 00415 pos = std::max(static_cast<unsigned int>(pos), minPos); 00416 00417 this->B.insert(this->B.begin() + pos, beta); 00418 this->U.insert(this->U.begin() + pos, TRANS(this->n)); 00419 this->U[pos].orbit(beta, S_i); 00420 return pos; 00421 } 00422 00423 template <class PERM, class TRANS> 00424 void BSGS<PERM, TRANS>::stripRedundantBasePoints(int minPos) { 00425 for (int i = this->B.size()-1; i >= minPos; --i) { 00426 if (this->U[i].size() <= 1) { 00427 if (i == static_cast<int>(this->B.size()-1)) { 00428 this->B.pop_back(); 00429 this->U.pop_back(); 00430 } else { 00431 this->B.erase(this->B.begin() + i); 00432 this->U.erase(this->U.begin() + i); 00433 } 00434 } 00435 } 00436 } 00437 00438 template <class PERM, class TRANS> 00439 void BSGS<PERM,TRANS>::copyTransversals(const BSGS<PERM,TRANS>& bsgs) { 00440 std::map<PERM*,typename PERM::ptr> genMap; 00441 BOOST_FOREACH(const typename PERM::ptr& p, bsgs.S) { 00442 typename PERM::ptr deepcopy(new PERM(*p)); 00443 //std::cout << "found " << p.get() << " = " << *p << std::endl; 00444 genMap.insert(std::make_pair(p.get(), deepcopy)); 00445 this->S.push_back(deepcopy); 00446 } 00447 00448 BOOST_ASSERT(this->B.size() == bsgs.B.size()); 00449 BOOST_ASSERT(bsgs.B.size() == bsgs.U.size()); 00450 this->U.clear(); 00451 this->U.resize(bsgs.U.size(), TRANS(bsgs.n)); 00452 BOOST_ASSERT(this->U.size() == bsgs.U.size()); 00453 00454 for (unsigned int i=0; i<this->U.size(); ++i) { 00455 this->U[i] = bsgs.U[i].clone(genMap); 00456 } 00457 } 00458 00459 } 00460 00461 #endif // -- BSGS_H_