permlib  0.2.6
Library for permutation computations
 All Classes Functions Variables Typedefs Enumerations Friends
include/permlib/bsgs.h
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_