Sieving.H

Go to the documentation of this file.
00001 
00006 #ifndef SIEVING_HEADER_
00007 #define SIEVING_HEADER_
00008 
00009 #include "qsieve-fwd.H"
00010 #include "StaticFactorbase.H"
00011 #include <cmath>
00012 
00013 
00014 
00015 extern int LogicalSieveSize; // sieving interval will be [-LogicalSieveSize,LogicalSieveSize] for each MPQS polynomial
00016 
00017 
00018 
00019 // !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00020 // !!! IMPORTANT !!! #define and typedef MUST match !!!
00021 // !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00022 #ifndef TSIEVEELEMENTSIZE
00023 #warning "TSIEVEELEMENTSIZE not defined. Using default value: 1"
00024 #define TSIEVEELEMENTSIZE 1 /* 1 = signed char, 2 = signed short int, 4 = signed int, 0 = no assembler inline code (=generic code) */
00025 #endif
00026 
00027 // if ASM_386 is defined, then TSIEVEELEMENTSIZE defines the type which is used in inline-assembler subroutines
00028 
00029 #if TSIEVEELEMENTSIZE == 0
00030  typedef signed char TSieveElement; // must be signed!
00031  #define log_SieveEntryMultiplier 1 /* sieve entries are casted to TSiebelement using this multiplier */
00032  #ifdef ASM_386
00033   #warning "assembly code is disabled for sieving!"
00034  #endif
00035 #elif TSIEVEELEMENTSIZE == 1
00036  typedef signed char TSieveElement; // must be signed!
00037  #define log_SieveEntryMultiplier 1 /* sieve entries are casted to TSiebelement using this multiplier */
00038 #elif TSIEVEELEMENTSIZE == 2
00039  typedef signed short int TSieveElement; // must be signed!
00040  #define log_SieveEntryMultiplier 16 /* sieve entries are casted to TSiebelement using this multiplier */
00041 #elif TSIEVEELEMENTSIZE == 4
00042  typedef signed int TSieveElement; // must be signed!
00043  #define log_SieveEntryMultiplier 256 /* sieve entries are casted to TSiebelement using this multiplier */
00044 #else
00045  #error "invalid value for TSIEVEELEMENTSIZE"
00046 #endif
00047 
00048 /*
00049   The sieve needs less memory (and less cache), if TSieveElement is chosen
00050   of the smallest possible type. Sieving is faster for short types.
00051 
00052   On the other hand, the granularity to detect relations is worse for
00053   shorter datatypes than for longer ones.
00054 
00055   You have also to consider processor specific effects (eg. alignment penalties)!
00056 
00057   Following datatypes are possible:
00058     - signed char (1 byte)
00059     - signed short int (2 bytes)
00060     - signed int (4 bytes)
00061     - float, double (floating point types, several (>3) bytes)
00062 
00063   Size of types may vary between systems. I haven't checked, whether
00064   the program works for other constraints than
00065     int >= 4 bytes, short int >= 2 bytes.
00066 
00067   Important remark:
00068    If you use floating point types for TSieveElement (which I don't
00069    recommend anymore), then calling "ceil" is counterproductive (and such
00070    calls should be removed in the source code).  But there is really no
00071    reason, why one should do sieving using floating point types.
00072 */
00073 
00074 /*
00075   about: log_SieveEntryMultiplier <VALUE>
00076   (-> refer TSIEVEELEMENTSIZE)
00077 
00078   Sieve entries are casted to TSieveElement using this multiplier to gain a
00079   better accuracy.  If the TSieveElement is of floating point type, then
00080   log_SieveEntryMultiplier can be set to 1. In all other cases you should
00081   set it to the highest possible value, such that the sieving threshold (for
00082   numbers usually to factorize) fits into TSieveElement.
00083 
00084   In order not to miss any detection of a relation when using short datatypes,
00085   you must use "ceil" while subtracting hits in the sieve. (Thereby the
00086   sieve gets close meshed.)
00087 */
00088 
00089 
00090 
00091 #if defined(IS_SERVER) || defined(IS_VALIDATOR)
00092 
00093  // The server does no sieving, but the server source code needs some
00094  // information about how sieving is done. Additionally the server sources
00095  // contain some code that is depending on data structures defined here. We
00096  // could eliminate these lines by conditional compiling (#ifdef IS_SERVER),
00097  // but this would be to bloated. Instead of this, we provide a stubby class
00098  // for the SERVER to satisfy the linker.
00099  
00100 class SieveControl : protected StaticFactorbaseSettings
00101 {
00102  private:
00103   static const int FBLowestStartIndex = 5;
00104   /* This is the first index of the factorbase to sieve with.
00105      This value needs to be greater than 1, because the first two factors of
00106      the factorbase are -1 and (probably) 2. Fore these two numbers as well
00107      as for multiplier of kN we are unable to sieve. (Well, at least it would
00108      be more difficult.) For a fast, but dirty sieving you can choose a
00109      greater value as well.  But the higher the value is, the more relations
00110      remain undetected...
00111    */
00112 
00113   static void set_logVal_for_Primenumber(const int pos, const TSieveElement val, const TSieveElement dirtyval);
00114   static void set_logVal_for_Primepower(const int pos, const TSieveElement val);
00115 
00116  public:
00117   friend void StaticFactorbase::compute_StaticFactorbase();
00118 };
00119 
00120 
00121 // two dummy static functions, not inlined because they must be linked...
00122 void SieveControl::set_logVal_for_Primenumber(const int, const TSieveElement, const TSieveElement)
00123 {
00124   //cout << "set_logVal_for_Primenumber: do nothing..." << endl;
00125   // do nothing, as the server doesn't need this information!
00126 }
00127 
00128 void SieveControl::set_logVal_for_Primepower(const int, const TSieveElement)
00129 {
00130   //cout << "set_logVal_for_Primepower: do nothing..." << endl;
00131   // do nothing, as the server doesn't need this information!
00132 }
00133 
00134 #else /* IS_SERVER not defined */
00135 
00136 // these are the "full" declarations
00137 #include <iosfwd>
00138 #include <gmp.h>
00139 #include "utils.H"
00140 #include "mpz_wrapper.H"
00141 using namespace my_mpz_wrapper;
00142 
00143 using std::ostream;
00144 using std::istream;
00145 
00146 using std::cout;
00147 using std::cerr;
00148 
00149 
00150 extern int SieveOffset;
00151 /* In general the logical sieve interval for each mpqs polynomial will exceed
00152    the size of the physical sieve. The logical sieve has to be sieved in
00153    several parts/steps. SieveOffset is the offset of the logical sieve at
00154    the starting point of a physical sieve.
00155  */
00156 
00157 
00158 extern TSieveElement SieveArray_[PhysicalSieveSize+64+8] __attribute__ ((aligned (64))); // memory for the sieve
00159 TSieveElement* const SieveArray = &SieveArray_[64];
00160 // SieveArray[-64..-1] are -1 (for fast value-based end-of-loop range detections)
00161 // SieveArray[0..PhysicalSieveSize-1] contain actually sieve values
00162 // SieveArray[PhysicalSieveSize..PhysicalSieveSize+7] contain dummy values
00163 
00164 
00165 extern TSieveElement log_Primzahl_of_PrimePowers[StaticFactorbase::max_additional_Powers];
00166 // converted logarithms for further powers of prime numbers
00167 
00168 
00169 extern mpz_t n, kN, r;
00170 extern double Factor_Threshold;
00171 
00172 inline TSieveElement log_DynamicLargePrime (const int p)
00173 {
00174 #if defined(ASM_386) && TSIEVEELEMENTSIZE == 1 && log_SieveEntryMultiplier == 1
00175  #ifdef DEBUG
00176   #warning "assembler i386 TSIEVEELEMENT log lookup enabled"
00177  #endif
00178 
00179   register TSieveElement logval;
00180 
00181 #if 0
00182   // "old" crafted i386 assembler code, very fast,
00183   // we assume 8103<p<=1318815734, otherwise result is either 10 or 21.
00184   asm ( \
00185        "movb $21,%[L] \n\t" \
00186        "cmpl $485165195,%[P] # exp(20) \n\t" \
00187        "sbbb $0,%[L] \n\t" \
00188        "cmpl $178482300,%[P] # exp(19) \n\t" \
00189        "sbbb $0,%[L] \n\t" \
00190        "cmpl $65659969,%[P] # exp(18) \n\t" \
00191        "sbbb $0,%[L] \n\t" \
00192        "cmpl $24154952,%[P] # exp(17) \n\t" \
00193        "sbbb $0,%[L] \n\t" \
00194        "cmpl $8886110,%[P] # exp(16) \n\t" \
00195        "sbbb $0,%[L] \n\t" \
00196        "cmpl $3269017,%[P] # exp(15) \n\t" \
00197        "sbbb $0,%[L] \n\t" \
00198        "cmpl $1202604,%[P] # exp(14) \n\t" \
00199        "sbbb $0,%[L] \n\t" \
00200        "cmpl $442413,%[P] # exp(13) \n\t" \
00201        "sbbb $0,%[L] \n\t" \
00202        "cmpl $162754,%[P] # exp(12) \n\t" \
00203        "sbbb $0,%[L] \n\t" \
00204        "cmpl $59874,%[P] # exp(11) \n\t" \
00205        "sbbb $0,%[L] \n\t" \
00206        "cmpl $22026,%[P] # exp(10) \n\t" \
00207        "sbbb $0,%[L]" \
00208        : [L] "=&q" (logval) : [P] "r" (p) : "cc");
00209 #endif
00210 
00211 
00212    static const int TREE[14] __attribute__ ((aligned (64))) =
00213     // heap of floor(exp(x+6)) values
00214      {
00215        // TREE[0..1], x: 4, 12, 
00216        22026, 65659969,
00217        // TREE[2..5], x: 2, 6, 10, 14,
00218        2980, 162754, 8886110, 485165195,
00219        // TREE[6..13], x: 1, 3, 5, 7, 9, 11, 13, 15
00220        1096, 8103, 59874, 442413, 3269017, 24154952, 178482300, 1318815734
00221      };
00222 
00223   // compute ceil(ln(p)) for 403<p<=2147483647
00224   // via binary tree search (without branch instructions).
00225   // (out of range values: 7 for 0<=p<404, 22 for (unsigned int)(p)>2147483647)
00226 
00227 #if 0
00228   // i386 standard code,
00229   // ugly slow, because of memory access latencies
00230   asm ( \
00231        "xorl %%ecx,%%ecx \n\t" \
00232        "cmpl $1202604,%[P] \n\t" \
00233        "setnc %%cl \n\t" \
00234        "cmpl %[P],%[F](,%%ecx,4) \n\t" \
00235        "rclb $1,%%cl \n\t" \
00236        "cmpl %[P],%[F]+4*2(,%%ecx,4) \n\t" \
00237        "rclb $1,%%cl \n\t" \
00238        "cmpl %[P],%[F]+4*6(,%%ecx,4) \n\t" \
00239        "rclb $1,%%cl \n\t" \
00240        "addb $7,%%cl \n" \
00241        : [L] "=&c" (logval) : [F] "o" (TREE[0]), [P] "r" (p) : "cc");
00242 #endif
00243 
00244 #if 1
00245   // crafted Athlon code, i386 compatible,
00246   // fast, because most latencies are resolved
00247   asm ( \
00248        "xorl %%ecx,%%ecx \n\t" \
00249        "cmpl $65659969,%[P] \n\t" \
00250        "sbb $-1,%%cl \n\t" \
00251        "cmp $22026,%[P] \n\t" \
00252        "sbb $-1,%%cl \n\t" \
00253        "cmpl $1202604,%[P] \n\t" \
00254        "sbb $-1,%%cl \n\t" \
00255        "cmpl %[P],%[F]+4*2(,%%ecx,4) \n\t" \
00256        "lea (,%%ecx,4),%%eax \n\t" \
00257        "adcb $0,%%al \n\t" \
00258        "cmpl %[P],%[F]+4*6(,%%ecx,8) \n\t" \
00259        "adcb $0,%%al \n\t" \
00260        "cmpl %[P],%[F]+4*7(,%%ecx,8) \n\t" \
00261        "adcb $7,%%al \n" \
00262        : [L] "=&a" (logval) : [F] "o" (TREE[0]), [P] "r" (p) : "cc", "ecx");
00263 #endif
00264 
00265 
00266 
00267 #if 0 || defined(DEBUG)
00268   // remark: do not use std::log(), as it might be ambiguous for some
00269   // versions of g++ (-> see utils.H)
00270   if (logval!=ceil(log_SieveEntryMultiplier*log(p)))
00271    {
00272      MARK;
00273      cout << "P: " << p << " log:" << static_cast<int>(logval) << "  " << ceil(log_SieveEntryMultiplier*log(p)) << endl;
00274      exit(1);
00275    }
00276 #endif
00277   return logval;
00278 #else
00279   // generic version
00280   return static_cast<TSieveElement>(ceil(log_SieveEntryMultiplier*log(p)));
00281 #endif
00282 }
00283 
00284 
00286 class SieveControl : protected StaticFactorbaseSettings
00287 {
00288  private:
00289   static const int FBLowestStartIndex = 5;
00290   /* factorbase elements below this index will be elided while sieving.
00291      Index must be greater than 1, because the first two elements of the
00292      static factorbase are -1 and probably 2. For these values (and the
00293      prefactor of kN) sieving is not (easily) possible.
00294      For fast (very) dirty-sieving a greater value can be chosen. But too
00295      big indices decrease rate of detection of relations.
00296    */
00297 
00298   inline static double DirtyFactor(const int i)
00299    {
00300      if (i<FBLowestStartIndex+4) return 1.0;
00301      return 1.0/sqrt(i); // a function slowly going towards 0
00302    }
00303   /*
00304      DirtyFactor denotes the weight, with which the i-th prime number of the
00305      static factor base will hin into the sieve, if the hit is "dirty".
00306 
00307      Example:
00308       If you sieve the normal way using a prime number p, then you have a hit
00309       every p-th position in the sieve. This hit will be weighted by log(p).
00310 
00311       If  p is "dirty" then a hit will be weighted using DirtyFactor*log(p) for
00312       *every* sieve entry! Only if a relation is detected later, then this
00313       error (dirty-sieving effect) will be corrected to check, whether the
00314       hit would be a hit under normal circumstances. -- Doing this speeds up
00315       the sieving because less primes (especially the lower ones) need to be
00316       sieved at the cost of correcting the effect later. As normally only a
00317       small fraction of the sieve entries lead to relations, this shouldn't
00318       be a problem if the dirty prime numbers are dynamically adapted to the
00319       rate of hits.
00320 
00321       A DirtyFactor of 1.0 behaves neutral compared with normal sieving except
00322       for (potential) rounding errors, but it speeds up the sieving.
00323 
00324       The smaller the DirtyFactor, the more prime numbers can be dirty.
00325       (However, the rate of detection per sieve interval decreases in favour
00326       of sieving speed.)
00327 
00328       A DirtyFactor of 1/n behaves roughly as if each n-th dirty prime
00329       number would be a member in every relation.
00330    */
00331 
00332 protected:
00333   static TSieveElement log_PrimeNumbers[StaticFactorbase::MaxSize] __attribute__ ((aligned (16)));
00334 private:
00335    static unsigned int PhysicalSieveIntervalIndex; // number of the current physical interval to sieve with
00336    // starts with 0 for the start of each new logical sieve interval
00337    // and is increased by 1 for each new physical sieve interval
00338 
00339    static unsigned int BestPSI[2];
00340    // indices of the two best PhysicalSieveIntervals
00341 
00342    static TSieveElement PhysicalSieveIntervalThresholdCorrector[1024];
00343    // idea: the logical sieve interval differs slightly in quality for different
00344    // physical sieve intervals. So let's roughly correct the threshold for
00345    // each physical sieve interval from time to time. This may increase the quality of
00346    // detected relations.
00347    // SieveControl::GetLogThreshold() will always add this correction value to its
00348    // computed value. So it behaves totally neutral. You may change these values on the fly...
00349    // initial default values: 0, therefore neutral and no correction.
00350 
00351   static TSieveElement log_PrimeNumbers_dirty[StaticFactorbase::MaxSize] __attribute__ ((aligned (16)));
00352     // these are the converted logarithms of primenumbers from the static factorbase
00353 
00354   static TSieveElement raw_log_Threshold, log_Threshold, log_ThresholdDeltaHint; // will be initialized later!
00355   static int myFBStartIndex, myFBStartIndexHint;
00356   static signed int DirtinessPegel;
00357 
00358 public:
00359   inline static TSieveElement GetRawLogThreshold() { return raw_log_Threshold; }
00360   inline static TSieveElement GetLogThreshold()
00361    {
00362      return log_Threshold+PhysicalSieveIntervalThresholdCorrector[PhysicalSieveIntervalIndex];
00363    }
00364   inline static TSieveElement GetLogDirtyCorrection() { return raw_log_Threshold-log_Threshold; }
00365   inline static void RenewThresholds()
00366    {
00367      log_Threshold+=log_ThresholdDeltaHint;
00368      log_ThresholdDeltaHint=0; // reset
00369      myFBStartIndex=myFBStartIndexHint;
00370    }
00371   static void RecalibrateThresholds();
00372   static void create_sieve_image(int* &Image, int mySieveOffset);
00373 
00374   inline static int FBStartIndex() { return myFBStartIndex; }
00375   static void compute_SieveThreshold();
00376   static bool Hit_after_DirtyCorrection(const signed int SievePos, TSieveElement Value);
00377 
00378   inline static TSieveElement PrimzahlSieveWeight(const int i)
00379    {
00380      // returns the weight of the prime concerning the sieve;
00381      // usually this is the (scaled) logarithm which has been put
00382      // in the table before.
00383      // This function encapsulates access to this value to make it read-only.
00384      return log_PrimeNumbers[i];
00385    }
00386   inline static TSieveElement PrimzahlSieveWeight_dirty(const int i)
00387    {
00388      // returns the weight of the prime concerning the sieve, if the
00389      // prime has been marked as "dirty" before, that is: for performance
00390      // reasons we omit using this prime for sieving.
00391      // Usually this value is used to lower the threshold beforehand.
00392      // This function encapsulates access to this value to make it read-only.
00393      return log_PrimeNumbers_dirty[i];
00394    }
00395 
00396 private:
00397 
00398   static void set_logVal_for_Primenumber(const int pos, const TSieveElement val, const TSieveElement dirtyval);
00399   static void set_logVal_for_Primepower(const int pos, const TSieveElement val);
00400 
00401   static void IncreaseDirtiness(void)
00402    {
00403      // decrease threshold and increment index for first element to sieve with
00404      if (myFBStartIndexHint>(StaticFactorbase::Size()>>2)) return; // this would be far too much dirtiness!
00405      log_ThresholdDeltaHint-=PrimzahlSieveWeight_dirty(myFBStartIndexHint);
00406      ++myFBStartIndexHint;
00407 #ifdef VERBOSE_INFO
00408      cout << "suggested Dirtiness now " << myFBStartIndexHint << endl;
00409 #endif
00410    }
00411 
00412   static void DecreaseDirtiness(void)
00413    {
00414      // increment threshold and decrement index for first element to sieve with
00415      if (myFBStartIndexHint<=FBLowestStartIndex) return; // impossible to decrease
00416      log_ThresholdDeltaHint+=PrimzahlSieveWeight_dirty(--myFBStartIndexHint);
00417 #ifdef VERBOSE_INFO
00418      cout << "suggested Dirtiness now " << myFBStartIndexHint << endl;
00419 #endif
00420    }
00421 
00422   static void RequestDirtiness(const signed int Level)
00423    {
00424      DirtinessPegel+=Level;
00425      if (DirtinessPegel>1000000)
00426       {
00427         DirtinessPegel=0;
00428         IncreaseDirtiness();
00429       } else
00430      if (DirtinessPegel<-1000000)
00431       {
00432         DirtinessPegel=0;
00433         DecreaseDirtiness();
00434       }
00435    }
00436 
00437  public:
00438   static inline void RequestMoreDirtiness(const signed int Level = 1000)
00439    {
00440      RequestDirtiness(Level);
00441    }
00442   static inline void RequestLessDirtiness(const signed int Level = 50)
00443    {
00444      RequestDirtiness(-Level);
00445    }
00446 
00447   friend void StaticFactorbase::compute_StaticFactorbase();
00448   friend void initialize_Sieve();
00449   friend void initialize_Sieve(const int Image[]);
00450   friend void do_sieving_full_intervals();
00451   friend void do_sieving_partial_intervals();
00452 };
00453 
00454 #include <qsieve.H>
00455 #include <vector>
00456 
00457 
00458 // just a declaration:
00459 namespace DynamicFactorArrays
00460 {
00461   extern unsigned int DynamicFactorsInUse;
00462   extern double DYNFB_threshold;
00463   void compute_Deltas_for_DynamicFactors(const int offset);
00464 }
00465 
00466 
00467 // helper class for computing Deltas
00468 class CDeltaComputations : protected StaticFactorbase
00469 {
00470  private:
00471              // a bit tricky:
00472              // the first two values are undefined,
00473              // then follow pairs of (FBPrime[i]*FBPrime[i+1],D_mod(FBPrime[i]*FBPrime[i+1])
00474              // then follow the unpaired "normal" values...
00475              // (see the update()-function to understand this!)
00476   static int D_mod_FBPrime[StaticFactorbase::MaxSize];
00477   static signed int D_difference_to_last; // 0 for undefined; all other values represent the real difference
00478   class CmpqsDmodPUpdater
00479    {
00480     public:
00481      mpz_t act_D;
00482      CmpqsDmodPUpdater()
00483       {
00484         mpz_init_set_d(act_D,-10e40); // to guarantee that the first MPQS interval triggers an update!
00485         for (int i=0; i<StaticFactorbase::MaxSize; ++i) D_mod_FBPrime[i]=0;
00486       }
00487      ~CmpqsDmodPUpdater()
00488       {
00489         mpz_clear(act_D);
00490       }
00491      void update();
00492    };
00493 
00494   static CmpqsDmodPUpdater mpqsDmodPUpdater;
00495   static unsigned int get_A2_mod_FBPrime(const int Nr)
00496    {
00497 #ifdef DEBUG
00498      if (Nr<3 || Nr>StaticFactorbase::Size()) { MARK; exit(777); }
00499 #endif
00500      // A2 := A*2 = 2*D^2; see mpqsPolynom.H: method get_A2_mod(const unsigned int m)
00501      unsigned int h = numtheory::squaremod(D_mod_FBPrime[Nr],PrimeNumbers[Nr])<<1;
00502      if (h>=static_cast<unsigned int>(PrimeNumbers[Nr])) h-=PrimeNumbers[Nr];
00503      return h;
00504    }
00505   static unsigned int get_A2_mod_FBPrimePair(const int Nr)
00506    {
00507 #ifdef DEBUG
00508      if (Nr&1 || PrimeNumbers[Nr+1]>=46340) { MARK; exit(777); } // see comments above!
00509 #endif
00510      // A2 := A*2 = 2*D^2; see mpqsPolynom.H: method get_A2_mod(const unsigned int m)
00511      unsigned int h = numtheory::squaremod(D_mod_FBPrime[Nr+1],D_mod_FBPrime[Nr])<<1;
00512      if (h>=static_cast<unsigned int>(D_mod_FBPrime[Nr])) h-=D_mod_FBPrime[Nr];
00513      return h;
00514    }
00515 
00516 
00517  protected:
00518   static int Delta_of_PrimeNumbers[StaticFactorbase::MaxSize][2] __attribute__ ((aligned (16)));
00519   static int Delta_of_PrimePowers[StaticFactorbase::max_additional_Powers][2] __attribute__ ((aligned (16)));
00520   static int LocalPhysInterval_Delta_of_PrimeNumbers[StaticFactorbase::MaxSize][2] __attribute__ ((aligned (16)));
00521   static int LocalPhysInterval_Delta_of_PrimePowers[StaticFactorbase::max_additional_Powers][2] __attribute__ ((aligned (16)));
00522 
00523  private:
00524   static void recompute_all_Deltas(const bool DoReallyAll=true);
00525   static void GetStaticFactorHits(std::vector<unsigned int> &Factors, const signed int SievePos);
00526 
00527  friend void do_sieving_full_intervals();
00528  friend void do_sieving_partial_intervals();
00529  friend CRelation::CRelation(const signed int SievePos, short int HitCount /* =0 */);
00530  friend bool SieveControl::Hit_after_DirtyCorrection(const signed int SievePos, TSieveElement Value);
00531  friend void DynamicFactorArrays::compute_Deltas_for_DynamicFactors(const int offset);
00532 };
00533 
00534 
00535 // helper class to encapsulate private data
00536 class CSieveStaticFactors : private CDeltaComputations, protected SieveControl
00537 {
00538  private:
00539 
00540 #if (TSIEVEELEMENTSIZE == 1)
00541   // helper routines for speeding up asm sieving greatly!
00542   // idea: since LogValues for prime numbers don't differ for most
00543   // consecutive static factors (and these values are monotone increasing),
00544   // we can count runs of LogValues.
00545   static unsigned short int LogValRuns[128];
00546   static unsigned short int PrimeNumberDiffs[MaxSize];
00547   static void init_LogValRuns()
00548    {
00549      for (int i=StaticFactorbase::Size()-1; i>1; --i)
00550       {
00551         unsigned int h=PrimeNumbers[i]-PrimeNumbers[i-1];
00552         //if (h>100) cout << "i: " << i << " h: " << h << endl;
00553         PrimeNumberDiffs[i]=h;
00554       }
00555 
00556      register TSieveElement LogVal = PrimzahlSieveWeight(StaticFactorbase::Size()-1);
00557      register unsigned int c = 0;
00558      int pos = 0;
00559      for (int i = StaticFactorbase::Size()-1; PrimeNumbers[i]>=PhysicalSieveSize; --i)
00560       {
00561         if (SieveControl::PrimzahlSieveWeight(i)==LogVal) ++c;
00562         else
00563          {
00564            if (c>65534)
00565             {
00566               MARK;
00567               cerr << "overflow!" << endl;
00568             }
00569            // store length of run of same LogVals
00570            LogValRuns[pos]=c; c=1; ++pos;
00571            if (--LogVal != PrimzahlSieveWeight(i))
00572            {
00573              cerr << "unexpected LogValue!" << endl;
00574              MARK;
00575              exit(1);
00576            }
00577          }
00578       }
00579      LogValRuns[pos++]=c;
00580      LogValRuns[pos++]=0; // stop marker
00581      LogValRuns[pos]=0; // additional stop marker
00582      //cout << "Runs: " << pos << ": ";
00583      //for (int i=0; i<pos; ++i) cout << LogValRuns[i] << "  ";
00584      //cout << endl; 
00585    }
00586 #endif
00587 
00588   static void create_Tables_for_init_LocalDeltas(int* &Table1, int* &Table2, const int SieveOffset);
00589   static void init_LocalDeltas() __attribute__ ((deprecated));
00590   static void init_LocalDeltas(const int Table1[], const int Table2[]);
00591   static void do_sieving_StaticFactors();
00592 
00593   static int* DefaultHelperTable[2];
00594 
00595 public:
00596  static void initialize()
00597   {
00598 #if (TSIEVEELEMENTSIZE == 1)
00599     init_LogValRuns();
00600 #endif
00601     create_Tables_for_init_LocalDeltas(DefaultHelperTable[0],DefaultHelperTable[1], -LogicalSieveSize);
00602      // this should be the default for sieving data: starting at -LogicalSieveSize
00603   }
00604 
00605  static void destruct()
00606   {
00607     delete [] DefaultHelperTable[1]; DefaultHelperTable[1]=NULL;
00608     delete [] DefaultHelperTable[0]; DefaultHelperTable[0]=NULL;
00609   }
00610 
00611  friend void do_sieving_full_intervals();
00612  friend void do_sieving_partial_intervals();
00613  friend int main(const int argc, const char* const argv[]);
00614 };
00615 
00616 
00617 void do_scanning_Sieve(void);
00618 
00619 void do_sieving_Squares();
00620 void do_sieving_DynamicFactors();
00621 void Werteverteilung_im_Sieb_ausgeben(void);
00622 
00623 extern void (*do_sieving)(); // function pointer to sieve method that is used
00624 void do_sieving_full_intervals();
00625 
00626 #include "DynamicFactorRelations.H"
00627 extern TDynamicFactorRelations DynamicFactorRelations;
00628 
00629 #endif /* IS_SERVER not defined */
00630 
00631 
00632 // this is used to store the sieved dynamic factors.
00633 // Remark: To collect the dynamic factors efficiently in the sieving phase, I decided to make this array global.
00634 //         (The constructor of "CRelation" is called with the number of found dynamic factors for each relation.)
00635 class CHits
00636 {
00637  public: static const short int MaxHits=64; // "64 dynamic factors should be enough for everyone"
00638  public: int Faktor;
00639  private: int Exponent;
00640  public:
00641   inline int GetExponent() const { return Exponent; }
00642  friend CRelation::CRelation(const signed int SievePos, short int HitCount /* =0 */);
00643 };
00644 
00645 extern CHits Hits[CHits::MaxHits]; // collected dynamic factors
00646 
00647 #endif /* SIEVING_HEADER_ */

Generated on Wed Nov 7 23:29:26 2007 for Qsieve by  doxygen 1.5.4