algorithm - por - metodos de factorizacion

¿Cuál es el algoritmo de factorización más rápido?

Question

He escrito un programa que intenta encontrar Parejas Amigas. Esto requiere encontrar las sumas de los divisores adecuados de los números.

Aquí está mi método sumOfDivisors() actual:

int sumOfDivisors(int n)
{  
    int sum = 1;
    int bound = (int) sqrt(n);
    for(int i = 2; i <= 1 + bound; i++)
    {
        if (n % i == 0)
            sum = sum + i + n / i;
    } 
    return sum;
}

Así que tengo que hacer muchas factorizaciones y eso se está convirtiendo en el verdadero cuello de botella en mi aplicación. Escribí un gran número en MAPLE y lo descompuse increíblemente rápido.

¿Cuál es uno de los algoritmos de factorización más rápidos?

Answer 1

Tiró directamente de mi respuesta a esta otra pregunta .

El método funcionará, pero será lento. "¿Qué tan grandes son tus números?" determina el método a usar:

• Menos de 2 ^ 16 o así: tabla de búsqueda.
• Menos de 2 ^ 70 más o menos: la modificación de Richard Brent del algoritmo rho de Pollard .
• Menos de 10 ^ 50: factorización de curva elíptica Lenstra
• Menos de 10 ^ 100: Tamiz cuadrático
• Más de 10 ^ 100: tamiz de campo de número general

Answer 2

Algoritmo de Shor: http://en.wikipedia.org/wiki/Shor%27s_algorithm

Por supuesto, necesitas una computadora cuántica: D

Answer 3

Depende qué tan grandes son tus números. Si está buscando pares amigables, está haciendo muchas factorizaciones, por lo que la clave puede no ser factorizar lo más rápido posible, sino compartir el mayor trabajo posible entre diferentes llamadas. Para acelerar la división de prueba, puede mirar la memorización y / o calcular previamente los números primos hasta la raíz cuadrada del número más grande que le interese. Es más rápido obtener la factorización principal, luego calcular la suma de todos los factores a partir de eso, de lo que es recorrer todo el camino hasta sqrt (n) para cada número.

Si busca pares realmente grandes, digamos más grandes que 2 ^ 64, entonces en un pequeño número de máquinas no puede hacerlo factorizando cada número individual sin importar qué tan rápido sea su factorización. Los atajos que estás usando para encontrar candidatos pueden ayudarte a factorizarlos.

Answer 4

Este es un documento de la factorización de enteros en Maple.

"A partir de algunas instrucciones muy simples:" hacer que la factorización de enteros sea más rápida en Maple ", hemos implementado el algoritmo de factorización de tamiz cuadrático en una combinación de Maple y C ..."

http://www.cecm.sfu.ca/~pborwein/MITACS/papers/percival.pdf

Answer 5

La pregunta en el título (y la última línea) parece tener poco que ver con el cuerpo real de la pregunta. Si está tratando de encontrar pares amistosos, o calcular la suma de divisores para muchos números, factorizar por separado cada número (incluso con el algoritmo más rápido posible) es absolutamente una manera ineficiente de hacerlo.

La función de suma de divisores , σ(n) = (sum of divisors of n) , es una función multiplicativa : para primos relativamente m y n, tenemos σ(mn) = σ(m)σ(n) , entonces

σ (p ₁ ^{k ₁} ... p _r ^{k _r} ) = [(p ₁ ^{k ₁ +1} -1) / (p ₁ -1)] ... [(p _r ^{k _r +1} -1) / (p _r -1 )].

De modo que utilizaría cualquier tamiz simple (por ejemplo, una versión aumentada del Tamiz de Eratóstenes ) para encontrar los números primos hasta n , y, en el proceso, la factorización de todos los números hasta n. (Por ejemplo, al hacer el tamiz, almacene el factor primo más pequeño de cada n. Luego puede factorizar cualquier número n iterando). Esto sería más rápido (en general) que usar cualquier algoritmo de factorización por separado varias veces.

Por cierto: ya existen varias listas conocidas de pares amistosos (ver, por ejemplo, here y los enlaces en MathWorld ). ¿Estás intentando ampliar el registro o hacerlo solo por diversión?

Answer 6

Sugeriría comenzar con el mismo algoritmo utilizado en Maple, el Tamiz Cuadrático .

1. Elija su número impar n para factorizar,
2. Elige un número natural k ,
3. Buscar todos p <= k para que k ^ 2 no sea congruente con (n mod p) para obtener una base de factores B = p1, p2, ..., pt ,
4. Comenzando desde r > floor (n) buscar al menos t + 1 valores para que y ^ 2 = r ^ 2 - n todos tengan solo factores en B,
5. Por cada y1 , y2 , ..., y (t + 1) que acabas de calcular generas un vector v (yi) = (e1, e2, ..., et) donde ei se calcula reduciendo sobre el módulo 2 el exponente pi en yi ,
6. Utilice la eliminación gaussiana para encontrar que algunos de los vectores que se agregaron dieron un vector nulo
7. Establezca x como el producto de ri relacionado con yi encontrado en el paso anterior y establezca y como p1 ^ a * p2 ^ b * p3 ^ c * .. * pt ^ z donde exponentes son la mitad de los exponentes encontrados en la factorización de yi
8. Calcule d = mcd (xy, n) , si 1 <d <n, entonces d es un factor no trivial de n ; de lo contrario, comience desde el paso 2 eligiendo una k más grande.

El problema con estos algoritmos es que realmente implican mucha teoría en cálculo numérico ...

Answer 7

Una implementación de tabla de búsqueda de C ++ versión 2 ²⁷ más 2015 para 1 GB de memoria:

#include <iostream.h> // cerr, cout, and NULL
#include <string.h>   // memcpy()
#define uint unsigned __int32
uint *factors;
const uint MAX_F=134217728; // 2^27

void buildFactors(){
   factors=new (nothrow) uint [(MAX_F+1)*2]; // 4 * 2 * 2^27 = 2^30 = 1GB
   if(factors==NULL)return; // not able to allocate enough free memory
   int i;
   for(i=0;i<(MAX_F+1)*2;i++)factors[i]=0;

   //Sieve of Eratosthenese
   factors[1*2]=1;
   factors[1*2+1]=1;
   for(i=2;i*i<=MAX_F;i++){
      for(;factors[i*2] && i*i<=MAX_F;i++);
      factors[i*2]=1;
      factors[i*2+1]=i;
      for(int j=2;i*j<=MAX_F;j++){
         factors[i*j*2]=i;
         factors[i*j*2+1]=j;
      }
   }
   for(;i<=MAX_F;i++){
      for(;i<=MAX_F && factors[i*2];i++);
      if(i>MAX_F)return;
      factors[i*2]=1;
      factors[i*2+1]=i;
   }
}

uint * factor(uint x, int &factorCount){
   if(x > MAX_F){factorCount=-1;return NULL;}
   uint tmp[70], at=x; int i=0;
   while(factors[at*2]>1){
      tmp[i++]=factors[at*2];
      cout<<"at:"<<at<<" tmp:"<<tmp[i-1]<<endl;
      at=factors[at*2+1];
   }
   if(i==0){
      cout<<"at:"<<x<<" tmp:1"<<endl;
      tmp[i++]=1;
      tmp[i++]=x;
   }else{
      cout<<"at:"<<at<<" tmp:1"<<endl;
      tmp[i++]=at;
   }
   factorCount=i;
   uint *ret=new (nothrow) uint [factorCount];
   if(ret!=NULL)
      memcpy(ret, tmp, sizeof(uint)*factorCount);
   return ret;
}

void main(){
   cout<<"Loading factors lookup table"<<endl;
   buildFactors(); if(factors==NULL){cerr<<"Need 1GB block of free memory"<<endl;return;}
   int size;
   uint x=30030;
   cout<<"/nFactoring: "<<x<<endl;
   uint *f=factor(x,size);
   if(size<0){cerr<<x<<" is too big to factor. Choose a number between 1 and "<<MAX_F<<endl;return;}
   else if(f==NULL){cerr<<"ran out of memory trying to factor "<<x<<endl;return;}

   cout<<"/nThe factors of: "<<x<<" {"<<f[0];
   for(int i=1;i<size;i++)
      cout<<", "<<f[i];
   cout<<"}"<<endl;
   delete [] f;

   x=30637;
   cout<<"/nFactoring: "<<x<<endl;
   f=factor(x,size);
   cout<<"/nThe factors of: "<<x<<" {"<<f[0];
   for(int i=1;i<size;i++)
      cout<<", "<<f[i];
   cout<<"}"<<endl;
   delete [] f;
   delete [] factors;
}

Actualización: o sacrificar algo de simplicidad por un rango un poco más allá de 2 ²⁸

#include <iostream.h> // cerr, cout, and NULL
#include <string.h>   // memcpy(), memset()

//#define dbg(A) A
#ifndef dbg
#define dbg(A)
#endif

#define uint   unsigned __int32
#define uint8  unsigned __int8
#define uint16 unsigned __int16

uint * factors;
uint8  *factors08;
uint16 *factors16;
uint   *factors32;

const uint LIMIT_16   = 514; // First 16-bit factor, 514 = 2*257
const uint LIMIT_32   = 131074;// First 32-bit factor, 131074 = 2*65537
const uint MAX_FACTOR = 268501119;
//const uint64 LIMIT_64 = 8,589,934,594; // First 64-bit factor, 2^33+1

const uint TABLE_SIZE = 268435456; // 2^28 => 4 * 2^28 = 2^30 = 1GB 32-bit table
const uint o08=1, o16=257 ,o32=65665; //o64=4294934465
// TableSize = 2^37 => 8 * 2^37 = 2^40 1TB 64-bit table
//   => MaxFactor = 141,733,953,600

/* Layout of factors[] array
*  Indicies(32-bit)              i                 Value Size  AFactorOf(i)
*  ----------------           ------               ----------  ----------------
*  factors[0..128]            [1..513]             8-bit       factors08[i-o08]
*  factors[129..65408]        [514..131073]        16-bit      factors16[i-o16]
*  factors[65409..268435455]  [131074..268501119]  32-bit      factors32[i-o32]
*
* Note: stopping at i*i causes AFactorOf(i) to not always be LargestFactor(i)
*/
void buildFactors(){
dbg(cout<<"Allocating RAM"<<endl;)
   factors=new (nothrow) uint [TABLE_SIZE]; // 4 * 2^28 = 2^30 = 1GB
   if(factors==NULL)return; // not able to allocate enough free memory
   uint i,j;
   factors08 = (uint8 *)factors;
   factors16 = (uint16 *)factors;
   factors32 = factors;
dbg(cout<<"Zeroing RAM"<<endl;)
   memset(factors,0,sizeof(uint)*TABLE_SIZE);
   //for(i=0;i<TABLE_SIZE;i++)factors[i]=0;

//Sieve of Eratosthenese
     //8-bit values
dbg(cout<<"Setting: 8-Bit Values"<<endl;)
   factors08[1-o08]=1;
   for(i=2;i*i<LIMIT_16;i++){
      for(;factors08[i-o08] && i*i<LIMIT_16;i++);
dbg(cout<<"Filtering: "<<i<<endl;)
      factors08[i-o08]=1;
      for(j=2;i*j<LIMIT_16;j++)factors08[i*j-o08]=i;
      for(;i*j<LIMIT_32;j++)factors16[i*j-o16]=i;
      for(;i*j<=MAX_FACTOR;j++)factors32[i*j-o32]=i;
   }
   for(;i<LIMIT_16;i++){
      for(;i<LIMIT_16 && factors08[i-o08];i++);
dbg(cout<<"Filtering: "<<i<<endl;)
      if(i<LIMIT_16){
         factors08[i-o08]=1;
         j=LIMIT_16/i+(LIMIT_16%i>0);
         for(;i*j<LIMIT_32;j++)factors16[i*j-o16]=i;
         for(;i*j<=MAX_FACTOR;j++)factors32[i*j-o32]=i;
      }
   }i--;

dbg(cout<<"Setting: 16-Bit Values"<<endl;)
     //16-bit values
   for(;i*i<LIMIT_32;i++){
      for(;factors16[i-o16] && i*i<LIMIT_32;i++);
      factors16[i-o16]=1;
      for(j=2;i*j<LIMIT_32;j++)factors16[i*j-o16]=i;
      for(;i*j<=MAX_FACTOR;j++)factors32[i*j-o32]=i;
   }
   for(;i<LIMIT_32;i++){
      for(;i<LIMIT_32 && factors16[i-o16];i++);
      if(i<LIMIT_32){
         factors16[i-o16]=1;
         j=LIMIT_32/i+(LIMIT_32%i>0);
         for(;i*j<=MAX_FACTOR;j++)factors32[i*j-o32]=i;
      }
   }i--;

dbg(cout<<"Setting: 32-Bit Values"<<endl;)
     //32-bit values
   for(;i*i<=MAX_FACTOR;i++){
      for(;factors32[i-o32] && i*i<=MAX_FACTOR;i++);
      factors32[i-o32]=1;
      for(j=2;i*j<=MAX_FACTOR;j++)factors32[i*j-o32]=i;
   }
   for(;i<=MAX_FACTOR;i++){
      for(;i<=MAX_FACTOR && factors32[i-o32];i++);
      if(i>MAX_FACTOR)return;
      factors32[i-o32]=1;
   }
}

uint * factor(uint x, int &factorCount){
   if(x > MAX_FACTOR){factorCount=-1;return NULL;}
   uint tmp[70], at=x; int i=0;
   while(at>=LIMIT_32 && factors32[at-o32]>1){
      tmp[i++]=factors32[at-o32];
dbg(cout<<"at32:"<<at<<" tmp:"<<tmp[i-1]<<endl;)
      at/=tmp[i-1];
   }
   if(at<LIMIT_32){
      while(at>=LIMIT_16 && factors16[at-o16]>1){
         tmp[i++]=factors16[at-o16];
dbg(cout<<"at16:"<<at<<" tmp:"<<tmp[i-1]<<endl;)
         at/=tmp[i-1];
      }
      if(at<LIMIT_16){
         while(factors08[at-o08]>1){
            tmp[i++]=factors08[at-o08];
dbg(cout<<"at08:"<<at<<" tmp:"<<tmp[i-1]<<endl;)
            at/=tmp[i-1];
         }
      }
   }
   if(i==0){
dbg(cout<<"at:"<<x<<" tmp:1"<<endl;)
      tmp[i++]=1;
      tmp[i++]=x;
   }else{
dbg(cout<<"at:"<<at<<" tmp:1"<<endl;)
      tmp[i++]=at;
   }
   factorCount=i;
   uint *ret=new (nothrow) uint [factorCount];
   if(ret!=NULL)
      memcpy(ret, tmp, sizeof(uint)*factorCount);
   return ret;
}
uint AFactorOf(uint x){
   if(x > MAX_FACTOR)return -1;
   if(x < LIMIT_16) return factors08[x-o08];
   if(x < LIMIT_32) return factors16[x-o16];
                    return factors32[x-o32];
}

void main(){
   cout<<"Loading factors lookup table"<<endl;
   buildFactors(); if(factors==NULL){cerr<<"Need 1GB block of free memory"<<endl;return;}
   int size;
   uint x=13855127;//25255230;//30030;
   cout<<"/nFactoring: "<<x<<endl;
   uint *f=factor(x,size);
   if(size<0){cerr<<x<<" is too big to factor. Choose a number between 1 and "<<MAX_FACTOR<<endl;return;}
   else if(f==NULL){cerr<<"ran out of memory trying to factor "<<x<<endl;return;}

   cout<<"/nThe factors of: "<<x<<" {"<<f[0];
   for(int i=1;i<size;i++)
      cout<<", "<<f[i];
   cout<<"}"<<endl;
   delete [] f;

   x=30637;
   cout<<"/nFactoring: "<<x<<endl;
   f=factor(x,size);
   cout<<"/nThe factors of: "<<x<<" {"<<f[0];
   for(int i=1;i<size;i++)
      cout<<", "<<f[i];
   cout<<"}"<<endl;
   delete [] f;
   delete [] factors;
}