h1. Ciurul lui Erathostene
h1. Ciurul lui Eratostene
(Creat de '_Cosmin_':user/Cosmin la data de _2004-11-24_ categoria _Teoria numerelor_, autor(i) _Cosmin_)
(Categoria _Matematica_, Autor _Cosmin Negruseri_)
Articolul de fata incearca o implementarea mai eficienta a acestui algoritm clasic. Se poate optimiza pentru a folosi doar O(sqrt(n)) memorie, varianta prezentata aici folosind O(n / log n) memorie, unde log n e numarul de biti al unui cuvant.
Articolul de fata incearca o implementarea mai eficienta a acestui algoritm clasic. Se poate optimiza pentru a folosi doar $O(sqrt(n))$ memorie, varianta prezentata aici folosind $O(n / log n)$ memorie, unde log n e numarul de biti al unui cuvant.
==Include(page="template/raw")==
Ciurul lui Erathostene e un algoritm clasic care se invata la scoala impreuna cu conceptul de numere prime inca din clasa a 6. Acest algoritm determina toate numerele prime mai mici decat un numar dat ca parametru. Ideea lui de abordare a acestei probleme poate fi modificata pentru a rezolva si alte probleme precum problema Fractii din arhiva infoarena, problema "Riemann vs Mertens":http://acm.uva.es/p/v107/10738.html din arhiva "uva":http://acm.uva.es/, problema "Divizibilitate":http://algoritmus.org/probleme/Probleme_Runda04.php a rundei 4 a concursului algoritmus, sau problema "Square Free":http://www.topcoder.com/stat?c=problem_statement&pm=2342&rd=4770 a SRM-ului 190 de pe TOPCODER, rezolvarea ei o gasiti "aici":http://www.topcoder.com/index?t=statistics&c=srm190_prob.
Ciurul lui Eratostene e un algoritm clasic care se invata la scoala impreuna cu conceptul de numere prime inca din clasa a 6. Acest algoritm determina toate numerele prime mai mici decat un numar dat ca parametru. Ideea lui de abordare a acestei probleme poate fi modificata pentru a rezolva si alte probleme precum problema "Fractii":problema/Fractii din arhiva infoarena, problema "Riemann vs Mertens":http://acm.uva.es/p/v107/10738.html din arhiva "uva":http://acm.uva.es/, problema "Divizibilitate":http://algoritmus.org/probleme/Probleme_Runda04.php a rundei 4 a concursului algoritmus, sau problema "Square Free":http://www.topcoder.com/stat?c=problem_statement&pm=2342&rd=4770 a SRM-ului 190 de pe TOPCODER, rezolvarea ei o gasiti "aici":http://www.topcoder.com/index?t=statistics&c=srm190_prob.
In articolul acesta nu ne vom concentra asupra acestor probleme ci asupra implementarii optimizate ale algoritmului original.
Am mai scris despre aceste optimizari intr-un articol din ginfo, dar codul de acolo nu era testat si nu merge :).
Ideea la acest algoritm e ca marcam intr-un sir fiecare multiplu al unui numar prim si numerele ramase nemarcate sunt numere prime, o descriere mai grafica gasiti la "adresa":http://mathworld.wolfram.com/SieveofEratosthenes.html.
Sa incercam o prima implementare a acestui algoritm (implementarile vor folosi limbajul java, dar sunt foarte usor transformabile in C/C++).
Sa incercam o prima implementare a acestui algoritm (implementarile vor folosi limbajul Java, dar sunt foarte usor transformabile in C/C++).
p(pre). {@ //class PrimeNumbersSieve1@}
{@ final int MAXSIZE = 1000001;@}
{@ char[] p = new char[MAXSIZE];@}
{@ @}
{@ //p[i] == 0 if i is prime@}
{@ @}
{@ public int getTheNumber(int n) {@}
{@ int i, j, nr = 0;@}
{@ for (i = 2; i <= n; ++i) {@}
{@ if (p[i] @}=={@ 0) {@}
{@ nr++;@}
{@ for (j = i + i; j <= n; j += i) {@}
{@ p[j] = 1;@}
{@ }@}
{@ }@}
{@ }@}
{@ return nr;@}
{@ }@}
== code(java) |
//class PrimeNumbersSieve1
final int MAXSIZE = 1000001;
char[] p = new char[MAXSIZE];
//p[i] == 0 if i is prime
public int getTheNumber(int n) {
int i, j, nr = 0;
for (i = 2; i <= n; ++i) {
if (p[i] == 0) {
nr++;
for (j = i + i; j <= n; j += i) {
p[j] = 1;
}
}
}
return nr;
}
==
O prima idee de optimizare ar fi sa nu mai luam in calcul numerele pare pentru ca stim ca singurul numar prim par e {$2$}. Deci sa vedem noua varianta a programului:
p(pre). {@ //class PrimeNumbersSieve2@}
{@ final int MAXSIZE = 1000001;@}
{@ char[] p = new char[MAXSIZE];@}
{@ @}
{@ //p[i] == 0 if i is prime@}
{@ @}
{@ public int getTheNumber(int n) {@}
{@ int i, j, nr = 1;@}
{@ for (i = 3; i <= n; i += 2) {@}
{@ if (p[i] @}=={@ 0) {@}
{@ nr++;@}
{@ for (j = i + i + i; j <= n; j += i << 1) {@}
{@ p[j] = 1;@}
{@ }@}
{@ }@}
{@ }@}
{@ return nr;@}
{@ }@}
== code(java) |
//class PrimeNumbersSieve2
final int MAXSIZE = 1000001;
char[] p = new char[MAXSIZE];
//p[i] == 0 if i is prime
public int getTheNumber(int n) {
int i, j, nr = 1;
for (i = 3; i <= n; i += 2) {
if (p[i] == 0) {
nr++;
for (j = i + i + i; j <= n; j += i << 1) {
p[j] = 1;
}
}
}
return nr;
}
==
Putem incerca o optimizare de memorie, pentru ca nu mai avem nevoie de elementele cu index par din sirul {$p$}. Acum semnificatia lui $p{~i~}$ s-a schimbat $p{~i~}$ fiind 0 daca $2*i+1$ e numar prim si $1$ daca nu.
p(pre). {@ // class PrimeNumbersSieve3 @}
{@ final int MAXSIZE = 1000000/2+1;@}
{@ char[] p = new char[MAXSIZE];@}
{@ @}
{@ //p[i] @}=={@ 0 if 2*i + 1 is prime@}
{@ @}
{@ public int getTheNumber(int n) {@}
{@ int i, j, nr = 1;@}
{@ for (i = 1; (i << 1) + 1 <= n; i += 1) {@}
{@ if (p[i] @}=={@ 0) {@}
{@ nr++;@}
{@ for (j = i + i + i + 1; (j << 1) + 1 <= n; j += (i << 1) + 1) {@}
{@ p[j] = 1;@}
{@ }@}
{@ }@}
{@ }@}
{@ return nr;@}
{@ }@}
Urmatoarea optimizare va fi marcarea multiplilor numarului prim $i$ de la $i*i$ nu de la $2*i$ cum am facut in prima varianta sau de la $3*i$ cum am facut in a {$2$}-a. Aceasta optimizare este evidenta: orice numar prim compus multiplu de $i$ mai mic decat $i*i$ are un factor prim mai mic decat {$i$}, si acel factor l-a marcat mai devreme, deci nu are rost sa il marcam si la pasul {$i$}. Ideea aceasta este exact ideea ce se foloseste la rezolvarea problemei "Numere Prime":http://www.infoarena.ro/task/prim din arhiva infoarena. Sa vedem acum codul sursa:
p(pre). {@ //class PrimeNumbersSieve4 @}
{@ final int MAXSIZE = 1000000/2+1;@}
{@ char[] p = new char[MAXSIZE];@}
{@ @}
{@ //p[i] @}=={@ 0 if 2*i + 1 is prime@}
{@ @}
{@ public int getTheNumber(int n) {@}
{@ int i, j, nr = 1;@}
{@ for (i = 1; ((i * i) << 1) + (i << 1) <= n; i += 1) {@}
{@ if (p[i] @}=={@ 0) {@}
{@ for (j = ((i * i) << 1) + (i << 1); (j << 1) + 1 <= n; j += (i << 1) + 1) {@}
{@ p[j] = 1;@}
{@ }@}
{@ }@}
{@ }@}
{@ for (i=1; 2 * i + 1 <= n; ++i) @}
{@ if (p[i] @}=={@ 0) nr++;@}
{@ return nr;@}
{@ }@}
== code(java) |
// class PrimeNumbersSieve3
final int MAXSIZE = 1000000/2+1;
char[] p = new char[MAXSIZE];
//p[i] == 0 if 2*i + 1 is prime
public int getTheNumber(int n) {
int i, j, nr = 1;
for (i = 1; (i << 1) + 1 <= n; i += 1) {
if (p[i] == 0) {
nr++;
for (j = i + i + i + 1; (j << 1) + 1 <= n; j += (i << 1) + 1) {
p[j] = 1;
}
}
}
return nr;
}
==
Urmatoarea optimizare va fi marcarea multiplilor numarului prim $i$ de la $i*i$ nu de la $2*i$ cum am facut in prima varianta sau de la $3*i$ cum am facut in a {$2$}-a. Aceasta optimizare este evidenta: orice numar prim compus multiplu de $i$ mai mic decat $i*i$ are un factor prim mai mic decat {$i$}, si acel factor l-a marcat mai devreme, deci nu are rost sa il marcam si la pasul {$i$}. Ideea aceasta este exact ideea ce se foloseste la rezolvarea problemei "Numere Prime":problema/prim din arhiva infoarena. Sa vedem acum codul sursa:
== code(java) |
//class PrimeNumbersSieve4
final int MAXSIZE = 1000000/2+1;
char[] p = new char[MAXSIZE];
//p[i] == 0 if 2*i + 1 is prime
public int getTheNumber(int n) {
int i, j, nr = 1;
for (i = 1; ((i * i) << 1) + (i << 1) <= n; i += 1) {
if (p[i] == 0) {
for (j = ((i * i) << 1) + (i << 1); (j << 1) + 1 <= n; j += (i << 1) + 1) {
p[j] = 1;
}
}
}
for (i=1; 2 * i + 1 <= n; ++i)
if (p[i] == 0) nr++;
return nr;
}
==
Codul sursa arata putin urat pentru ca nu lucram direct cu $i$ ci cu {$2*i+1$}, am mai facut optimizarea ce apare si in "mathworld":http://mathworld.wolfram.com/, nu parcurgem numerele pana la $n$ pentru marcarea multiplilor ci pana la $sqrt(n)$ lucru care e evident dupa cele explicate mai sus.
Ultima imbunatatire care o vom aduce este aceea de a folosi mai putina memorie. Cum pentru fiecare numar e necesara doar o informatie booleana, aceasta o putem tine intr-un bit, nu este necesar un char intreg. Sa vedem cum arata codul:
p(pre). {@/class PrimeNumbersSieve5 @}
{@ final int MAXSIZE = 100000000/2/8+1;@}
{@ char[] p = new char[MAXSIZE];@}
{@ @}
{@ //p[i] @}=={@ 0 if 2*i + 1 is prime@}
{@ @}
{@ public int getTheNumber(int n) {@}
{@ int i, j, nr = 1;@}
{@ for (i = 1; ((i * i) << 1) + (i << 1) <= n; i += 1) { @}
{@ if ((p[i >> 3] & (1 << (i & 7))) == 0) {@}
{@ for (j = ((i * i) << 1) + (i << 1); (j << 1) + 1 <= n; j += (i << 1) + 1) {@}
{@ p[j >> 3] |= (1 << (j & 7));@}
{@ }@}
{@ }@}
{@ }@}
{@ for (i = 1; 2 * i + 1 <= n; ++i) @}
{@ if ((p[i >> 3] & (1 << (i & 7))) @}=={@ 0) @}
{@ nr++;@}
{@ return nr;@}
{@ }@}
== code(java) |
//class PrimeNumbersSieve5
final int MAXSIZE = 100000000/2/8+1;
char[] p = new char[MAXSIZE];
//(p[i>>3]&(1<<(i&7))) == 0 if 2*i + 1 is prime
public int getTheNumber(int n) {
int i, j, nr = 1;
for (i = 1; ((i * i) << 1) + (i << 1) <= n; i += 1) {
if ((p[i >> 3] & (1 << (i & 7))) == 0) {
for (j = ((i * i) << 1) + (i << 1); (j << 1) + 1 <= n; j += (i << 1) + 1) {
p[j >> 3] |= (1 << (j & 7));
}
}
}
for (i = 1; 2 * i + 1 <= n; ++i)
if ((p[i >> 3] & (1 << (i & 7))) == 0)
nr++;
return nr;
}
==
Codul {@p[i >> 3] & (1 << (i & 7)) @}=={@ 0@} testeaza daca al {@i@}-lea bit din sirul de biti e $0$ (deci daca {@2*i+1@} e prim). {@i >> 3@} e echivalent cu {@i / 8@} deci se gaseste pt bitul al {@i@}-lea in ce char e el stocat din sirul nostru de charuri. {@i & 7@} e echivalent cu {@i % 8@} ca sa aflam pe ce bit al charului e stocat numarul prim {@i@}.
Codul folosit pentru testare e:
p(pre). {@public static void main(String[] arg) {@}
{@ double tick = System.currentTimeMillis();@}
{@ for (int i = 0; i < 10; ++i) {@}
{@ System.out.println("Run " + i @}
{@ + " result: "+@}
{@ new PrimeNumbersSieve1().@}
{@ getTheNumber(1000000));@}
{@ }@}
{@ System.out.println(" 10 runs of "@}
{@ + "PrimeNumbersSieve1()." @}
{@ + "getTheNumber(1000000) " @}
{@ + "lasted: " @}
{@ + (System.currentTimeMillis() @}
{@ - tick) @}
{@ * 1e-3 + " seconds");@}
{@ }@}
== code(java) |
public static void main(String[] arg) {
double tick = System.currentTimeMillis();
for (int i = 0; i < 10; ++i) {
System.out.println("Run " + i
+ " result: "+
new PrimeNumbersSieve1().
getTheNumber(1000000));
}
System.out.println(" 10 runs of "
+ "PrimeNumbersSieve1()."
+ "getTheNumber(1000000) "
+ "lasted: "
+ (System.currentTimeMillis()
- tick)
* 1e-3 + " seconds");
}
==
S-a rulat fiecare cod de $10$ ori pentru a da ocazia optimizarilor Just In Time din java sa isi faca treaba.
Se pare ca optimizarile codului initial au dus la o imbunatatire a vitezei cu un factor de {$10$}.
Problemele sugerate la inceputul articolului sunt destul de frumoase si ar trebui rezolvate. Pentru fiecare dintre ele puteti sa va testati rezolvarea: la algoritmus sunt disponibile fisierele de intrare si iesire, infoarena si uva au evaluatoare automate si pentru problema de pe TOPCODER puteti intra in applet-ul lor si sa incercati sa o rezolvati in practice roomul asociat SRM-ului 190 dupa ce o rezolvati folositi optiunea Run System Test din Practice Room Options.
Problemele sugerate la inceputul articolului sunt destul de frumoase si ar trebui rezolvate. Pentru fiecare dintre ele puteti sa va testati rezolvarea: la algoritmus sunt disponibile fisierele de intrare si iesire, infoarena si uva au evaluatoare automate si pentru problema de pe "TOPCODER":http://www.topcoder.com/tc puteti intra in applet-ul lor si sa incercati sa o rezolvati in practice roomul asociat SRM-ului 190 dupa ce o rezolvati folositi optiunea Run System Test din Practice Room Options.
