0votos

Números primos hasta N en C

por josejuan hace 4 años

Usando la misma bruta que antes pero considerando sólo los impares (el bitarray enumera los impares) toma 3 segundos.

Obtener todos los números primos menores o iguales a N en el menor tiempo posible.

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
#include <stdio.h>  
#include <stdlib.h>  
#include <string.h>  
#include <math.h>  
#include <omp.h>  
  
#define WORD unsigned char  
#define BITS 8  
#define POW 3  
  
int main(int argc, char **argv) {  
    int n = 1000;  
    if (argc > 1)  
        n = atoi(argv[1]); 
 
    // ajustamos al tamaño de palabra indicado para tener los indicados o más 
    n = (n >> POW) << POW; 
    printf("Checking to %i\n", n); 
 
    // el mayor impar estará en (i << 1) + 1 
    int i = n >> 1; 
 
    int t = 1 + (int) sqrt(n); 
 
    WORD *M = (WORD *) malloc(sizeof(WORD) * (i + 1));  
    memset(M, 0, i + 1);  
  
#pragma omp parallel  
    {  
        int th  = omp_get_thread_num();  
        int nth = omp_get_num_threads();  
        printf("Thread %i / %i\n", th, nth);  
        WORD *W = M;  
        int j;  
        for (j = 1 + th; j < i; j += nth) {  
            if ((W[j >> POW] & (1 << (j & (BITS - 1)))) == 0) {  
                int p = (j << 1) + 1; // es primo en 1 sólo hilo, puede ser re-descartado con varios 
                if(p <= t) { 
                  int h = (p * p) >> 1; // su cuadrado indexado (los múltiplos anteriores han sido descartados por primos menores) 
                  for (; h < i; h += p)  
                      W[h >> POW] |= 1 << (h & (BITS - 1)); // atómica?  
            }  
        }  
    }  
  
    {  
        int z;  
        i--; 
        for (z = 0; z < 15; z++, i--) {  
            while (M[i >> POW] & (1 << (i & (BITS - 1)))) i--;  
            printf("%i\n", (i << 1) + 1);  
        }  
    }  
    return 0;  
 
/* 
 
   [josejuan@centella primes]$ export OMP_NUM_THREADS=14 
   [josejuan@centella primes]$ gcc -lm -fopenmp -m64 -O3 primes.c 
   [josejuan@centella primes]$ time -f "%E, %M" ./a.out 1000000000 
   Checking to 1000000000 
   Thread 3 / 14 
   Thread 5 / 14 
   Thread 1 / 14 
   Thread 6 / 14 
   Thread 0 / 14 
   Thread 7 / 14 
   Thread 2 / 14 
   Thread 8 / 14 
   Thread 4 / 14 
   Thread 9 / 14 
   Thread 10 / 14 
   Thread 11 / 14 
   Thread 12 / 14 
   Thread 13 / 14 
   999999937 
   999999929 
   999999893 
   999999883 
   999999797 
   999999761 
   999999757 
   999999751 
   999999739 
   999999733 
   999999677 
   999999667 
   999999613 
   999999607 
   999999599 
   0:03.08, 490288 
   [josejuan@centella primes]$ cat /proc/cpuinfo | grep -i '\(model name\|mhz\)' | head -n 2 
   model name      : AMD Phenom(tm) II X6 1045T Processor 
   cpu MHz         : 2700.000 
 
 */ 
1 comentario
0votos

Escrito por josejuan hace 4 años

Podemos bajar a 2,7 segundos preprocesando múltiplos de 3
    ...

    n = (n >> 6) << 6; // forzar múltiplo 64

    ...

    //Cambiado el `memset(M, 0, i + 1);` por
    {
      uint64_t *m = (uint64_t *) M;
      int L = i >> 6;
      int j = 0;
      while(j < L) {
        m[j++] = 0x2492492492492492;
        m[j++] = 0x9249249249249249;
        m[j++] = 0x4924924924924924;
      }
      if(j < L) m[j++] = 0x2492492492492492;
      if(j < L) m[j++] = 0x9249249249249249;
      if(j < L) m[j++] = 0x4924924924924924;
    }

    ...

    // y empezando en 5
    for (j = 2 + th; ...

Comenta la solución

Tienes que identificarte para poder publicar tu comentario.