permlib  0.2.6
Library for permutation computations
 All Classes Functions Variables Typedefs Enumerations Friends
include/permlib/search/base_search.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 BASE_SEARCH_H_
00034 #define BASE_SEARCH_H_
00035 
00036 #include <permlib/change/conjugating_base_change.h>
00037 #include <permlib/change/random_base_transpose.h>
00038 
00039 #include <boost/dynamic_bitset.hpp>
00040 
00041 namespace permlib {
00042 
00044 template<class BSGSIN, class TRANSRET>
00045 class BaseSearch {
00046 public:
00047         typedef typename BSGSIN::PERMtype PERM;
00048         typedef typename BSGSIN::TRANStype TRANS;
00049         typedef std::list<typename PERM::ptr> PERMlistType;
00050         
00052 
00057         BaseSearch(const BSGSIN& bsgs, unsigned int pruningLevelDCM, bool stopAfterFirstElement);
00058         
00060         virtual ~BaseSearch(){}
00061         
00063 
00066         bool minOrbit(unsigned long alpha, BSGS<PERM,TRANSRET> &groupK, unsigned int i, unsigned long beta_i) const;
00067         
00069         virtual typename PERM::ptr searchCosetRepresentative();
00070         
00072 
00078         virtual typename PERM::ptr searchCosetRepresentative(BSGS<PERM,TRANSRET> &groupK, BSGS<PERM,TRANSRET> &groupL) = 0;
00079         
00081         unsigned long m_statNodesVisited;
00083         unsigned long m_statNodesPrunedCosetMinimality;
00085         unsigned long m_statNodesPrunedCosetMinimality2;
00087         unsigned long m_statNodesPrunedChildRestriction;
00088 
00089 protected:
00091         BSGSIN m_bsgs;
00093         BSGSIN* m_bsgs2;
00095         boost::scoped_ptr<SubgroupPredicate<PERM> > m_pred;
00096         
00098         std::vector<unsigned long> m_order;
00100         boost::scoped_ptr<BaseSorterByReference> m_sorter;
00101         
00103         ConjugatingBaseChange<PERM,TRANS,RandomBaseTranspose<PERM,TRANS> > m_baseChange;
00104 
00106         const unsigned int m_pruningLevelDCM;
00108         bool m_limitInitialized;
00110         unsigned int m_limitBase;
00112         unsigned int m_limitLevel;
00113         
00115         bool pruneDCM(const PERM& t, unsigned int backtrackLevel, BSGS<PERM,TRANSRET> &groupK, BSGS<PERM,TRANSRET> &groupL);
00117         bool checkLeaf(unsigned int level);
00119         unsigned int processLeaf(const PERM& t, unsigned int level, unsigned int backtrackLevel, unsigned int completed, BSGS<PERM,TRANSRET> &groupK, BSGS<PERM,TRANSRET> &groupL);
00121         virtual const std::vector<dom_int>& subgroupBase() const = 0;
00122         
00124         void setupEmptySubgroup(BSGS<PERM,TRANSRET>& group) const;
00125         
00127         const bool m_stopAfterFirstElement;
00129         typename PERM::ptr m_lastElement;
00130 private:
00131         static PERMlistType ms_emptyList;
00132 };
00133 
00134 //
00135 // IMPLEMENTATION
00136 //
00137 
00138 template<class BSGSIN,class TRANSRET>
00139 typename BaseSearch<BSGSIN,TRANSRET>::PERMlistType BaseSearch<BSGSIN,TRANSRET>::ms_emptyList;
00140 
00141 
00142 template<class BSGSIN,class TRANSRET>
00143 BaseSearch<BSGSIN,TRANSRET>::BaseSearch(const BSGSIN& bsgs, unsigned int pruningLevelDCM, bool stopAfterFirstElement) 
00144         : m_statNodesVisited(0), m_statNodesPrunedCosetMinimality(0), m_statNodesPrunedCosetMinimality2(0),
00145           m_statNodesPrunedChildRestriction(0),
00146           m_bsgs(bsgs), m_bsgs2(0), m_pred(0), m_baseChange(m_bsgs),
00147           m_pruningLevelDCM(pruningLevelDCM),
00148           m_limitInitialized(false), m_limitBase(0), m_limitLevel(0),
00149           m_stopAfterFirstElement(stopAfterFirstElement),
00150           m_lastElement()
00151 { 
00152 }
00153 
00154 
00155 template<class BSGSIN,class TRANSRET>
00156 bool BaseSearch<BSGSIN,TRANSRET>::minOrbit(unsigned long alpha, BSGS<PERM,TRANSRET> &groupK, unsigned int i, unsigned long beta_i) const {
00157         PERMlistType S_i;
00158         std::copy_if(groupK.S.begin(), groupK.S.end(), std::back_inserter(S_i), PointwiseStabilizerPredicate<PERM>(groupK.B.begin(), groupK.B.begin() + i));
00159         if (S_i.empty()) {
00160                 if (alpha == beta_i)
00161                         return true;
00162                 return (*m_sorter)(beta_i, alpha);
00163         }
00164         
00165         //TODO: avoid multiple allocation?
00166         boost::dynamic_bitset<> orbitCharacteristic(m_bsgs.n);
00167         orbitCharacteristic.set(alpha, 1);
00168         std::list<unsigned long> orbit;
00169         orbit.push_back(alpha);
00170         for (std::list<unsigned long>::const_iterator it = orbit.begin(); it != orbit.end(); ++it) {
00171                 unsigned long beta = *it;
00172                 BOOST_FOREACH(const typename PERM::ptr& p, S_i) {
00173                         unsigned long beta_p = *p / beta;
00174                         if (!orbitCharacteristic[beta_p]) {
00175                                 orbitCharacteristic.set(beta_p, 1);
00176                                 orbit.push_back(beta_p);
00177                                 if ((*m_sorter)(beta_p, beta_i)) {
00178                                         PERMLIB_DEBUG(std::cout << "DCM2 beta_p = " << beta_p+1 << " , beta_i = " << beta_i+1 << std::endl;)
00179                                         return false;
00180                                 }
00181                         }
00182                 }
00183         }
00184         return true;
00185 }
00186 
00187 template<class BSGSIN,class TRANSRET>
00188 bool BaseSearch<BSGSIN,TRANSRET>::pruneDCM(const PERM& t, unsigned int backtrackLevel, BSGS<PERM,TRANSRET> &groupK, BSGS<PERM,TRANSRET> &groupL) {
00189         // change base only for the lower nodes in the tree
00190         if (backtrackLevel < m_pruningLevelDCM) {
00191                 //TODO: avoid multiple allocation?
00192                 std::vector<unsigned long> newBaseImage(subgroupBase().begin(), subgroupBase().end());
00193                 for (unsigned int j=0; j<=backtrackLevel; ++j)
00194                         newBaseImage[j] = t / newBaseImage[j];
00195                 //print_iterable(newBaseImage.begin(), newBaseImage.begin() + (backtrackLevel+1), 1, "new base image");
00196                 ConjugatingBaseChange<PERM,TRANSRET,RandomBaseTranspose<PERM,TRANSRET> > cbc(groupL);
00197                 cbc.change(groupL, newBaseImage.begin(), newBaseImage.begin() + (backtrackLevel+1), false);
00198                 //print_iterable(groupL.B.begin(), groupL.B.end(), 1, "new base");
00199         }
00200 
00201         const unsigned long alpha = groupK.B[backtrackLevel];
00202         
00203         for (unsigned int i = 0; i <= backtrackLevel; ++i) {
00204                 if (i == backtrackLevel || groupK.U[i].contains(alpha)) {
00205                         PERMLIB_DEBUG(std::cout << "DCM2 found " << (alpha+1) << " in U_" << i << "  btLevel " << backtrackLevel << std::endl;)
00206                         PERMLIB_DEBUG(std::cout << "     t = " << t << std::endl;)
00207                         
00208                         if (!minOrbit(t / alpha, groupL, i, t / groupK.B[i])) {
00209                                 PERMLIB_DEBUG(std::cout << "DCM2 : " << ((t / groupK.B[i]) + 1) << " // " << ((t / alpha) + 1) << std::endl;)
00210                                 PERMLIB_DEBUG(std::cout << " K = " << groupK << std::endl;)
00211                                 PERMLIB_DEBUG(std::cout << " L = " << groupL << std::endl;)
00212                                 return true;
00213                         }
00214                 }
00215                 if (t / groupK.B[i] != groupL.B[i])
00216                         return false;
00217         }
00218         return false;
00219 }
00220 
00221 template<class BSGSIN,class TRANSRET>
00222 bool BaseSearch<BSGSIN,TRANSRET>::checkLeaf(unsigned int level) {
00223         return m_limitInitialized && level >= m_limitLevel;
00224 }
00225 
00226 template<class BSGSIN,class TRANSRET>
00227 unsigned int BaseSearch<BSGSIN,TRANSRET>::processLeaf(const PERM& t, unsigned int level, unsigned int backtrackLevel, unsigned int completed, BSGS<PERM,TRANSRET> &groupK, BSGS<PERM,TRANSRET> &groupL) {
00228         PERMLIB_DEBUG(std::cout << "XXX level " << level << "  bLevel " << backtrackLevel << std::endl;)
00229         PERMLIB_DEBUG(std::cout << "XXX limitLevel " << m_limitLevel << "  limitBase " << m_limitBase << std::endl;)
00230         typedef typename PERM::ptr PERMptr;
00231         if ((*m_pred)(t)) {
00232                 if (m_stopAfterFirstElement) {
00233                         m_lastElement = typename PERM::ptr(new PERM(t));
00234                         return 0;
00235                 }
00236                 const bool isIdentity = t.isIdentity();
00237                 int posK = 0, posL = 0;
00238                 if (m_limitInitialized && level == m_limitLevel && isIdentity) {
00239                         PointwiseStabilizerPredicate<PERM> stabPred(m_bsgs.B.begin(), m_bsgs.B.begin() + m_limitBase);
00240                         BOOST_FOREACH(const PERMptr &s, m_bsgs.S) {
00241                                 if (stabPred(s)) {
00242                                         PERMLIB_DEBUG(std::cout << *s << " extended gen\n";)
00243                                         BOOST_ASSERT((*m_pred)(*s));
00244                                         PERMptr sK(new PERM(*s));
00245                                         PERMptr sL(new PERM(*s));
00246                                         posK = std::max(posK, groupK.insertGenerator(sK, false));
00247                                         posL = std::max(posL, groupL.insertGenerator(sL, false));
00248                                 }
00249                         }
00250                         //return completed;
00251                 }
00252                 if (!isIdentity) {
00253                         PERMptr genK(new PERM(t));
00254                         posK = std::max(posK, groupK.insertGenerator(genK, false));
00255                         PERMptr genL(new PERM(t));
00256                         posL = std::max(posL, groupL.insertGenerator(genL, false));
00257                         PERMLIB_DEBUG(std::cout << "-- accepted" << std::endl;)
00258                 }
00259                 groupK.updateOrbits(posK);
00260                 groupL.updateOrbits(posL);
00261                 return completed;
00262         }
00263         return level;
00264 }
00265 
00266 template<class BSGSIN,class TRANSRET>
00267 void BaseSearch<BSGSIN,TRANSRET>:: setupEmptySubgroup(BSGS<PERM,TRANSRET>& group) const {
00268         group.B = subgroupBase();
00269         group.U.resize(subgroupBase().size(), TRANSRET(this->m_bsgs.n));
00270         for (unsigned int i=0; i<subgroupBase().size(); ++i)
00271                 group.orbit(i, ms_emptyList);
00272 }
00273 
00274 template<class BSGSIN,class TRANSRET>
00275 typename BaseSearch<BSGSIN,TRANSRET>::PERM::ptr BaseSearch<BSGSIN,TRANSRET>::searchCosetRepresentative() {
00276         BSGS<PERM,TRANSRET> groupK(this->m_bsgs.n);
00277         BSGS<PERM,TRANSRET> groupL(this->m_bsgs.n);
00278         
00279         setupEmptySubgroup(groupK);
00280         setupEmptySubgroup(groupL);
00281         
00282         return this->searchCosetRepresentative(groupK, groupL);
00283 }
00284 
00285 }
00286 
00287 #endif // -- BASE_SEARCH_H_