permlib  0.2.6
Library for permutation computations
 All Classes Functions Variables Typedefs Enumerations Friends
include/permlib/test/primitivity_test.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 PRIMITIVITY_TEST_H_
00034 #define PRIMITIVITY_TEST_H_
00035 
00036 #include <permlib/prime_helper.h>
00037 
00038 #include <boost/foreach.hpp>
00039 #include <boost/utility.hpp>
00040 #include <vector>
00041 #include <list>
00042 
00043 namespace permlib {
00044 
00046 
00053 template<typename PERM>
00054 class PrimitivityTest {
00055 public:
00063         template<typename InputIterator>
00064         PrimitivityTest(const unsigned int n, InputIterator genBegin, InputIterator genEnd);
00065         
00071         bool blockOfImprimitivity(std::vector<dom_int>* minimalBlock) const;
00072         
00076         bool isPrimitive() const { return blockOfImprimitivity(NULL); }
00077         
00078 private:
00079         const unsigned int m_n;
00080         unsigned int m_primeLimit;
00081         const std::list<typename PERM::ptr> m_generators;
00082         
00083         bool fillTrivialBlock(std::vector<dom_int>* minimalBlock) const;
00084         
00085         static dom_int rep(dom_int kappa, std::vector<dom_int>& p);
00086         
00087         bool merge(dom_int kappa, dom_int lambda, std::vector<dom_int>& c, std::vector<dom_int>& p, std::vector<dom_int>& q, unsigned int& l) const;
00088 };
00089 
00090 
00091 
00092 template<typename PERM>
00093 template<typename InputIterator>
00094 PrimitivityTest<PERM>::PrimitivityTest(const unsigned int n, InputIterator genBegin, InputIterator genEnd)
00095         : m_n(n), m_primeLimit(m_n), m_generators(genBegin, genEnd)
00096 {
00097         for (const unsigned int* p = PrimeHelper::firstPrime(); p != PrimeHelper::lastPrime(); ++p) {
00098                 if (m_n % (*p) == 0) {
00099                         m_primeLimit = m_n / (*p);
00100                         break;
00101                 }
00102         }
00103 }
00104 
00105 
00106 template<typename PERM>
00107 bool PrimitivityTest<PERM>::blockOfImprimitivity(std::vector<dom_int>* minimalBlock) const {
00108         std::vector<dom_int> alphas(2);
00109         alphas[0] = 0;
00110         
00111         for (dom_int a = 1; a < m_n; ++a) {
00112                 alphas[1] = a;
00113                 
00114                 const unsigned int k = alphas.size();
00115                 unsigned int l = k - 1;
00116                 std::vector<dom_int> p(m_n);
00117                 std::vector<dom_int> q(m_n);
00118                 std::vector<dom_int> c(m_n);
00119                 
00120                 for (unsigned int i = 0; i < m_n; ++i) {
00121                         c[i] = 1;
00122                         p[i] = i;
00123                 }
00124                 
00125                 for (unsigned int i = 0; i < k - 1; ++i) {
00126                         p[alphas[i+1]] = alphas[0];
00127                         q[i] = alphas[i+1];
00128                 }
00129                 
00130                 bool tryNextAlpha = false;
00131                 c[alphas[0]] = k;
00132                 for (unsigned int i = 0; i < l; ++i) {
00133                         const dom_int gamma = q[i];
00134                         BOOST_FOREACH(const typename PERM::ptr& x, m_generators) {
00135                                 const dom_int delta = rep(gamma, p);
00136                                 if (merge(x->at(gamma), x->at(delta), c, p, q, l)) {
00137                                         tryNextAlpha = true;
00138                                         goto TRY_NEXT_ALPHA;
00139                                 }
00140                         }
00141                 }
00142 TRY_NEXT_ALPHA:
00143                 if (tryNextAlpha)
00144                         continue;
00145                 
00146                 for (unsigned int i = 0; i < m_n; ++i)
00147                         rep(i, p);
00148                 
00149                 const unsigned int minBlockSize = c[rep(alphas[0], p)];
00150                 if (minBlockSize < m_n) {
00151                         if (minimalBlock) {
00152                                 minimalBlock->clear();
00153                                 minimalBlock->reserve(minBlockSize);
00154                                 for (unsigned int i = 0; i < m_n; ++i)
00155                                         if (p[i] == p[alphas[0]])
00156                                                 minimalBlock->push_back(i);
00157                         }
00158                         return false;
00159                 }
00160         }
00161         
00162         return fillTrivialBlock(minimalBlock);
00163 }
00164 
00165 template<typename PERM>
00166 bool PrimitivityTest<PERM>::fillTrivialBlock(std::vector<dom_int>* minimalBlock) const {
00167         if (minimalBlock) {
00168                 minimalBlock->clear();
00169                 minimalBlock->resize(m_n);
00170                 for (unsigned int i = 0; i < m_n; ++i)
00171                         minimalBlock->at(i) = i;
00172         }
00173         return true;
00174 }
00175 
00176 template<typename PERM>
00177 dom_int PrimitivityTest<PERM>::rep(dom_int kappa, std::vector<dom_int>& p) {
00178         dom_int lambda = kappa;
00179         dom_int rho = p[lambda];
00180         while (rho != lambda) {
00181                 lambda = rho;
00182                 rho = p[lambda];
00183         }
00184         
00185         dom_int mu = kappa;
00186         rho = p[mu];
00187         while (rho != lambda) {
00188                 p[mu] = lambda;
00189                 mu = rho;
00190                 rho = p[mu];
00191         }
00192         
00193         return lambda;
00194 }
00195 
00196 template<typename PERM>
00197 bool PrimitivityTest<PERM>::merge(dom_int kappa, dom_int lambda, std::vector<dom_int>& c, std::vector<dom_int>& p, std::vector<dom_int>& q, unsigned int& l) const {
00198         dom_int phi = rep(kappa, p);
00199         dom_int psi = rep(lambda, p);
00200         
00201         if (phi != psi) {
00202                 dom_int mu, nu;
00203                 if (c[phi] >= c[psi]) {
00204                         mu = phi;
00205                         nu = psi;
00206                 } else {
00207                         mu = psi;
00208                         nu = phi;
00209                 }
00210                 p[nu] = mu;
00211                 c[mu] += c[nu];
00212                 if (c[mu] > m_primeLimit)
00213                         return true;
00214                 
00215                 q[l] = nu;
00216                 ++l;
00217         }
00218         
00219         return false;
00220 }
00221 
00222 }
00223 
00224 #endif