permlib  0.2.9
Library for permutation computations
r_base.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 RBASE_H_
34 #define RBASE_H_
35 
36 #include <permlib/predicate/subgroup_predicate.h>
37 
38 #include <permlib/search/base_search.h>
39 
40 #include <permlib/search/partition/partition.h>
41 #include <permlib/search/partition/refinement_family.h>
42 #include <permlib/search/partition/backtrack_refinement.h>
43 
44 #include <permlib/change/conjugating_base_change.h>
45 #include <permlib/change/random_base_transpose.h>
46 
47 #include <permlib/sorter/base_sorter.h>
48 
49 #include <utility>
50 
51 namespace permlib {
52 namespace partition {
53 
55 template<class BSGSIN,class TRANSRET>
56 class RBase : public BaseSearch<BSGSIN,TRANSRET> {
57 public:
58  typedef typename BaseSearch<BSGSIN,TRANSRET>::PERM PERM;
59  typedef typename BaseSearch<BSGSIN,TRANSRET>::TRANS TRANS;
60 
62 
67  RBase(const BSGSIN& bsgs, unsigned int pruningLevelDCM, bool stopAfterFirstElement = false);
68 
69  typedef typename Refinement<PERM>::RefinementPtr RefinementPtr;
70  typedef typename RefinementFamily<PERM>::PartitionPtr PartitionPtr;
71  typedef typename std::list<std::pair<PartitionPtr,RefinementPtr> >::const_iterator PartitionIt;
72 
74  void search(BSGS<PERM,TRANSRET> &groupK);
75 
78 protected:
81  Partition m_partition2;
82 
84 
89  void construct(SubgroupPredicate<PERM>* pred, RefinementFamily<PERM>* predRefinement);
90 
92  virtual unsigned int processNewFixPoints(const Partition& pi, unsigned int level);
93 
94  virtual const std::vector<dom_int>& subgroupBase() const;
95 private:
97  std::vector<dom_int> m_subgroupBase;
99  std::list<std::pair<PartitionPtr,RefinementPtr> > partitions;
100 
102  unsigned int search(PartitionIt pIt, Partition &pi, const PERM& t, const PERM* t2, unsigned int level, unsigned int backtrackLevel, unsigned int& completed, BSGS<PERM,TRANSRET> &groupK, BSGS<PERM,TRANSRET> &groupL);
103 
105  bool updateMappingPermutation(const BSGSIN& bsgs, const Partition& sigma, const Partition& pi2, PERM& t2) const;
106 };
107 
108 template<class BSGSIN,class TRANSRET>
109 RBase<BSGSIN,TRANSRET>::RBase(const BSGSIN& bsgs, unsigned int pruningLevelDCM, bool stopAfterFirstElement)
110  : BaseSearch<BSGSIN,TRANSRET>(bsgs, pruningLevelDCM, stopAfterFirstElement),
111  m_partition(bsgs.n), m_partition2(bsgs.n)
112 { }
113 
114 template<class BSGSIN,class TRANSRET>
116  this->m_pred.reset(pred);
117  typedef typename boost::shared_ptr<RefinementFamily<PERM> > RefinementFamilyPtr;
118  std::list<RefinementFamilyPtr> refinements;
119 
120  if (!this->m_bsgs.isSymmetricGroup()) {
121  RefinementFamilyPtr gr( new GroupRefinementFamily<PERM,TRANS>(this->m_bsgs) );
122  refinements.push_back(gr);
123  }
124 
125  if (predRefinement) {
126  RefinementFamilyPtr predR( predRefinement );
127  refinements.push_back(predR);
128  }
129 
130  PERMLIB_DEBUG(print_iterable(this->m_bsgs.B.begin(), this->m_bsgs.B.end(), 1, "orig BSGS");)
131 
132  Partition pi(m_partition);
133  while (pi.cells() < this->m_bsgs.n) {
134  PERMLIB_DEBUG(std::cout << std::endl << "PI0 = " << pi << std::endl;)
135  bool found = false;
136  do {
137  found = false;
138  unsigned int foo = 0;
139  BOOST_FOREACH(RefinementFamilyPtr ref, refinements) {
140  const unsigned int oldFixPointsSize = pi.fixPointsSize();
141  std::pair<PartitionPtr,RefinementPtr> newRef = ref->apply(pi);
142  if (newRef.first) {
143  partitions.push_back(newRef);
144  if (oldFixPointsSize < pi.fixPointsSize()) {
145  processNewFixPoints(pi, partitions.size());
146  }
147  //std::cout << "BSGS " << this->m_bsgs;
148  found = true;
149  }
150  ++foo;
151  }
152  } while(found);
153 
154  PERMLIB_DEBUG(std::cout << std::endl << "PI1 = " << pi << std::endl;)
155 
156  if (pi.cells() < this->m_bsgs.n) {
157  unsigned int alpha = -1;
158  //print_iterable(pi.fixPointsBegin(), pi.fixPointsEnd(), 1, " fix0");
159  //print_iterable(this->m_bsgs.B.begin(), this->m_bsgs.B.end(), 1, "bsgs0");
160  if (pi.fixPointsSize() < this->m_bsgs.B.size())
161  alpha = this->m_bsgs.B[pi.fixPointsSize()];
162  if (alpha >= this->m_bsgs.n) {
163  for (unsigned int i = 0; i < this->m_bsgs.n; ++i) {
164  if (std::find(pi.fixPointsBegin(), pi.fixPointsEnd(), i) == pi.fixPointsEnd()) {
165  alpha = i;
166  break;
167  }
168  }
169  }
170  BOOST_ASSERT( alpha < this->m_bsgs.n );
171  PERMLIB_DEBUG(std::cout << "choose alpha = " << alpha << std::endl;)
172  RefinementPtr br(new BacktrackRefinement<PERM>(this->m_bsgs.n, alpha));
173  BacktrackRefinement<PERM>* ref = dynamic_cast<BacktrackRefinement<PERM> *>(br.get());
174  ref->initializeAndApply(pi);
175  PartitionPtr newPi(new Partition(pi));
176  PERMLIB_DEBUG(std::cout << "BACKTRACK " << (ref->alpha()+1) << " in " << pi << " --> " << *newPi << std::endl;)
177  partitions.push_back(std::make_pair(newPi, br));
178 
179  processNewFixPoints(pi, partitions.size());
180 
181  //std::cout << "BSGS " << this->m_bsgs;
182  m_subgroupBase.push_back(ref->alpha());
183  }
184  }
185 
186  this->m_order = BaseSorterByReference::createOrder(this->m_bsgs.n, pi.fixPointsBegin(), pi.fixPointsEnd());
187  this->m_sorter.reset(new BaseSorterByReference(this->m_order));
188  for (typename std::list<std::pair<PartitionPtr,RefinementPtr> >::iterator pIt = partitions.begin(); pIt != partitions.end(); ++pIt) {
189  (*pIt).second->sort(*this->m_sorter, 0);
190  PERMLIB_DEBUG(std::cout << "SIGMA = " << *(*pIt).first << std::endl;)
191  }
192 
193  PERMLIB_DEBUG(print_iterable(this->m_order.begin(), this->m_order.end(), 0, "ORDER");)
194 }
195 
196 template<class BSGSIN,class TRANSRET>
197 unsigned int RBase<BSGSIN,TRANSRET>::processNewFixPoints(const Partition& pi, unsigned int level) {
198  const unsigned int basePos = this->m_baseChange.change(this->m_bsgs, pi.fixPointsBegin(), pi.fixPointsEnd(), true);
199  if (this->m_bsgs2)
200  this->m_baseChange.change(*this->m_bsgs2, pi.fixPointsBegin(), pi.fixPointsEnd(), true);
201  //print_iterable(pi.fixPointsBegin(), pi.fixPointsEnd(), 1, " fix");
202  PERMLIB_DEBUG(print_iterable(this->m_bsgs.B.begin(), this->m_bsgs.B.end(), 1, "change base");)
203  return basePos;
204 }
205 
206 template<class BSGSIN,class TRANSRET>
208  BOOST_ASSERT( this->m_pred != 0 );
209 
210  this->setupEmptySubgroup(groupK);
211 
212  unsigned int completed = partitions.size();
213  BSGS<PERM,TRANSRET> groupL(groupK);
214  PERM identH(this->m_bsgs.n);
215  search(partitions.begin(), m_partition2, PERM(this->m_bsgs.n), &identH, 0, 0, completed, groupK, groupL);
216 }
217 
218 template<class BSGSIN,class TRANSRET>
220  BOOST_ASSERT( this->m_pred != 0 );
221 
222  // !!!
223  //
224  // TODO: check that groupK and groupL have the right base (starting with subgroupBase)
225  //
226  // !!!
227 
228  unsigned int completed = partitions.size();
229  //BSGS<PERM,TRANS> groupL(groupK);
230  PERM t(this->m_bsgs.n);
231  PERM t2(this->m_bsgs.n);
232  BOOST_ASSERT(partitions.begin() != partitions.end());
233  const Partition& sigma = *((*(partitions.begin())).first);
234  if (sigma.fixPointsSize()) {
235  updateMappingPermutation(this->m_bsgs, sigma, this->m_partition2, t);
236  if (this->m_bsgs2)
237  updateMappingPermutation(*this->m_bsgs2, sigma, this->m_partition2, t2);
238  }
239  search(partitions.begin(), m_partition2, t, &t2, 0, 0, completed, groupK, groupL);
240 
242 }
243 
244 
245 
246 template<class BSGSIN,class TRANSRET>
247 unsigned int RBase<BSGSIN,TRANSRET>::search(PartitionIt pIt, Partition &pi, const PERM& t, const PERM* t2, unsigned int level, unsigned int backtrackLevel, unsigned int& completed, BSGS<PERM,TRANSRET> &groupK, BSGS<PERM,TRANSRET> &groupL) {
248  ++this->m_statNodesVisited;
249 
250  if (pIt == partitions.end() || this->checkLeaf(level)) {
251  PERMLIB_DEBUG(std::cout << "LEAF: " << pi << " with t = " << t << std::endl;)
252  return this->processLeaf(t, level, backtrackLevel, completed, groupK, groupL);
253  }
254 
255  const Partition& sigma = *((*pIt).first);
256  const RefinementPtr& ref = (*pIt).second;
257  ++pIt;
258 
259  unsigned int s = ref->alternatives();
260  const bool isBacktrack = ref->type() == Backtrack;
261  const bool isGroup = ref->type() == Group;
262  const PERM* tForRefinement = &t;
263 
264  if (isGroup) {
265  GroupRefinement<PERM,TRANS>* gref = static_cast<GroupRefinement<PERM,TRANS>*>(ref.get());
266  if (this->m_bsgs2 && gref->bsgs() == *this->m_bsgs2) {
267  tForRefinement = t2;
268  }
269  }
270 
271  ref->sort(*this->m_sorter, &pi);
272  typedef typename Refinement<PERM>::RefinementPtrIterator RefIt;
273  for (RefIt rIt = ref->backtrackBegin(); rIt != ref->backtrackEnd(); ++rIt) {
274  if (isBacktrack && s < groupK.U[backtrackLevel].size()) {
275  PERMLIB_DEBUG(std::cout << "PRUNE the rest: s=" << s << " < " << groupK.U[backtrackLevel].size() << std::endl;)
276  this->m_statNodesPrunedCosetMinimality += s;
277  break;
278  }
279 
280  --s;
281  RefinementPtr ref2 = *rIt;
282 
283  const unsigned int oldFixPointsSize = pi.fixPointsSize();
284  PERMLIB_DEBUG(std::cout << " refinement from " << pi << std::endl;)
285  const unsigned int strictRefinement = ref2->apply2(pi, *tForRefinement);
286  PERMLIB_DEBUG(std::cout << " to " << pi << " with " << strictRefinement << std::endl;)
287  PERMLIB_DEBUG(for(unsigned int jj=0; jj<level; ++jj) std::cout << " ";)
288  PERMLIB_DEBUG(std::cout << "NODE " << sigma << " ~~~> " << pi << std::endl;)
289  /*
290  for (unsigned int q = 0; q < level; ++q) std::cout << " ";
291  std::cout << " " << level << ": " << sigma << " <-> " << pi2 << " from " << pi << std::endl;
292  for (unsigned int q = 0; q < level; ++q) std::cout << " ";
293  std::cout << " t = " << t << std::endl;
294  */
295  if (!strictRefinement) {
296  PERMLIB_DEBUG(std::cout << "no strict refinement " << sigma << " -- " << pi << std::endl;)
297  ++this->m_statNodesPrunedChildRestriction;
298  continue;
299  }
300  if (pi.cells() != sigma.cells()) {
301  PERMLIB_DEBUG(std::cout << "cell number mismatch " << sigma << " -- " << pi << std::endl;)
302  ref2->undo(pi, strictRefinement);
303  ++this->m_statNodesPrunedChildRestriction;
304  continue;
305  }
306  if (pi.fixPointsSize() != sigma.fixPointsSize()) {
307  PERMLIB_DEBUG(std::cout << "fix point number mismatch " << sigma << " -- " << pi << std::endl;)
308  ref2->undo(pi, strictRefinement);
309  ++this->m_statNodesPrunedChildRestriction;
310  continue;
311  }
312  PERM tG(t);
313  PERM* tH = 0;
314  if (pi.fixPointsSize() != oldFixPointsSize) {
315  if (!updateMappingPermutation(this->m_bsgs, sigma, pi, tG)) {
316  PERMLIB_DEBUG(std::cout << "no t found " << sigma << " -- " << pi << "; tG = " << tG << std::endl;)
317  ref2->undo(pi, strictRefinement);
318  ++this->m_statNodesPrunedChildRestriction;
319  continue;
320  }
321  if (this->m_bsgs2) {
322  tH = new PERM(*t2);
323  if (!updateMappingPermutation(*this->m_bsgs2, sigma, pi, *tH)) {
324  PERMLIB_DEBUG(std::cout << "no t found " << sigma << " -- " << pi << "; tH = " << tH << std::endl;)
325  ref2->undo(pi, strictRefinement);
326  ++this->m_statNodesPrunedChildRestriction;
327  continue;
328  }
329  }
330  }
331  if (this->m_pruningLevelDCM && isBacktrack) {
332  if (this->pruneDCM(tG, backtrackLevel, groupK, groupL)) {
333  ++this->m_statNodesPrunedCosetMinimality2;
334  ref2->undo(pi, strictRefinement);
335  continue;
336  }
337  }
338  unsigned int ret = search(pIt, pi, tG, tH ? tH : t2, level+1, isBacktrack ? (backtrackLevel + 1) : backtrackLevel, completed, groupK, groupL);
339  delete tH;
340  PERMLIB_DEBUG(std::cout << "retract " << strictRefinement << " from " << pi << " to ";)
341  ref2->undo(pi, strictRefinement);
342  PERMLIB_DEBUG(std::cout << pi << std::endl;)
344  return 0;
345  if (ret < level)
346  return ret;
347  }
348 
349  completed = std::min(completed, level);
350  return level;
351 }
352 
353 template<class BSGSIN,class TRANSRET>
354 bool RBase<BSGSIN,TRANSRET>::updateMappingPermutation(const BSGSIN& bsgs, const Partition& sigma, const Partition& pi, PERM& t2) const {
355  typedef std::vector<unsigned int>::const_iterator FixIt;
356  std::vector<dom_int>::const_iterator bIt;
357  unsigned int i = 0;
358  FixIt fixSigmaIt = sigma.fixPointsBegin();
359  const FixIt fixSigmaEndIt = sigma.fixPointsEnd();
360  FixIt fixPiIt = pi.fixPointsBegin();
361  PERMLIB_DEBUG(print_iterable(bsgs.B.begin(), bsgs.B.end(), 1, "B ");)
362  PERMLIB_DEBUG(print_iterable(fixSigmaIt, fixSigmaEndIt, 1, "Sigma");)
363  for (bIt = bsgs.B.begin(); bIt != bsgs.B.end(); ++bIt, ++i) {
364  PERMLIB_DEBUG(std::cout << " base: " << (*bIt)+1 << std::endl;)
365  while (fixSigmaIt != fixSigmaEndIt && *fixSigmaIt != *bIt) {
366  PERMLIB_DEBUG(std::cout << " skipping " << (*fixSigmaIt)+1 << " for " << (*bIt)+1 << std::endl;)
367  ++fixSigmaIt;
368  ++fixPiIt;
369  }
370  if (fixSigmaIt == fixSigmaEndIt) {
371  PERMLIB_DEBUG(std::cout << " no more fix point found for " << (*bIt)+1 << std::endl;)
372  return true;
373  }
374  const unsigned int alpha = *fixSigmaIt;
375  const unsigned int beta = *fixPiIt;
376  if (t2 / alpha != beta) {
377  boost::scoped_ptr<PERM> u_beta(bsgs.U[i].at(t2 % beta));
378  if (u_beta) {
379  //std::cout << " multiply with " << *u_beta << " for " << alpha+1 << "," << beta+1 << " // base " << bsgs.B[i] + 1<< std::endl;
380  t2 ^= *u_beta;
381  } else {
382  //std::cout << "could not find a u_b with " << (t2 % beta) << " at " << i << "--" << bsgs.B[i] << " -- " << &bsgs << std::endl;
383  return false;
384  }
385  }
386 
387  ++fixSigmaIt;
388  ++fixPiIt;
389  }
390  return true;
391 }
392 
393 template<class BSGSIN,class TRANSRET>
394 const std::vector<dom_int>& RBase<BSGSIN,TRANSRET>::subgroupBase() const {
395  return m_subgroupBase;
396 }
397 
398 }
399 }
400 
401 #endif // -- RBASE_H_
bool initializeAndApply(Partition &pi)
applies (left-)refinement to partition and initializes refinement for future use in R-base ...
Definition: refinement.h:136
const bool m_stopAfterFirstElement
true iff the search can be stopped after the first element found with the desired property ...
Definition: base_search.h:127
virtual PERM::ptr searchCosetRepresentative()
searches for a coset representative if one exists
Definition: base_search.h:271
vector_t::const_iterator fixPointsEnd() const
iterator to the end of fix points
Definition: partition.h:167
base class for searching in a group
Definition: base_search.h:45
STL namespace.
partition
Definition: partition.h:48
A sorter that sorts a sequence (e.g. ) with respect to a given input ordering (e.g. a base)
Definition: base_sorter.h:113
-refinements for group membership
Definition: refinement_family.h:67
unsigned long alpha() const
alpha point chosen for backtracking
Definition: backtrack_refinement.h:110
static std::vector< unsigned long > createOrder(unsigned int size, InputIterator begin, InputIterator end)
constructs an ordering array with the same parameters as BaseSorter for use with BaseSorterByReferenc...
Definition: base_sorter.h:121
virtual const std::vector< dom_int > & subgroupBase() const
base of the sought subgroup
Definition: r_base.h:394
unsigned long cells() const
number of cells in this partition
Definition: partition.h:157
RBase(const BSGSIN &bsgs, unsigned int pruningLevelDCM, bool stopAfterFirstElement=false)
constructor
Definition: r_base.h:109
R-base for partition backtracking.
Definition: r_base.h:56
unsigned int fixPointsSize() const
number of fix points in this partition
Definition: partition.h:161
void construct(SubgroupPredicate< PERM > *pred, RefinementFamily< PERM > *predRefinement)
constructs an R-base for given predicate and refinement family
Definition: r_base.h:115
Partition m_partition
partition to base the backtrack tree on
Definition: r_base.h:80
backtrack refinement
Definition: backtrack_refinement.h:45
virtual unsigned int processNewFixPoints(const Partition &pi, unsigned int level)
callback when a new fix point appears during R-base construction
Definition: r_base.h:197
void search(BSGS< PERM, TRANSRET > &groupK)
perform search and store result in groupK
Definition: r_base.h:207
abstract base class for subgroup (and coset) predicates
Definition: subgroup_predicate.h:45
Represents a base and strong generating set (BSGS)
Definition: bsgs.h:58
std::vector< TRANS > U
transversals along the stabilizer chain
Definition: bsgs_core.h:59
represents a class of -refinements for a given problem
Definition: refinement_family.h:49
vector_t::const_iterator fixPointsBegin() const
iterator to the begin of fix points
Definition: partition.h:164
Definition: abstract_bsgs.h:49