Crivello di Eratostene in Java

[marakid]

Non conosco bene il Crivello, ma la prima idea che mi viene in mente è di suddividere il calcolo su più thread per parallelozzare il parallelizzabile

[sottovento]

Originariamente inviato da DoctorT

intanto posto la mia implementazione in java :

import java.util.BitSet;
import java.util.Vector;
import java.lang.Math;
import java.io.*;

class Eratostene{
    
	static int segmented_sieve(long limit, int segment_size) {
		int sqrt = (int) Math.sqrt((double)limit);
		int count = 1;
		long s = 2;
		long n = 3;
		
		BitSet sieve = new BitSet(segment_size);
		// togli il commento se vuoi i numeri stampati a video
		// System.out.println(2);
		System.out.println("sqrt="+ sqrt);
		
		// crea un array di numeri primi < di sqrt(limit)
		
		boolean[] is_prime = new boolean[sqrt+1];
		for (int i=0;i<=sqrt;i++) is_prime[i] = true;
		for (int i = 2; i * i <= sqrt; i++)
		  {
			if (is_prime[i])
				for (int j = i * i; j <= sqrt; j += i)
					is_prime[j] = false;
		}
	
		Vector<Integer> primes = new Vector<Integer> ();
		Vector<Integer> next = new Vector<Integer> ();
		
		for (long low = 0; low <= limit; low += segment_size)
		{
			sieve.set(0, segment_size-1);		
			// current segment = interval [low, high]
			long high = Math.min(low + segment_size - 1,  limit) ;
			// memorizza i primi occorrenti per incrociare i multiple
			for (; s * s <= high; s++)
			{
				if (is_prime[(int)s])
				{
					primes.add((int)s);
					next.add((int)(s * s - low));
				}
			}
		
			// crivello per il segmento corrente
			for (int i = 1; i < primes.size(); i++)
			{
				int j = next.get(i);
				for (int k = primes.get(i) * 2; j < segment_size; j += k)
					sieve.clear(j);
				next.set(i, j - segment_size);
			}

			for (; n <= high; n += 2) // ricerca i bit con indice dispari
				if (sieve.get((int)(n - low))) // n is a prime
				{
					count++;
					// togli il commento se vuoi i numeri stampati a video
					// System.out.println(n); 
				}
		};
		return count;
	}
	
	public static void main(String [] args) throws IOException {	

		int sieve_size = 65536;
		long n = 8000000000L;
		long inizio = System.nanoTime();

		System.out.println("Attendere");	
		System.out.println("generazione numeri primi ..."); 
		System.out.println("I primi generati sono: " + segmented_sieve(n,sieve_size));

		long fine = System.nanoTime();
		float time = (fine-inizio)/1000000;
		System.out.println(time + " ms");

	} // end main
} // end class

il tempo di esecuzione sul mio sistema è di circa 130 secondi ... mi sa che c'è un bel po' da ottimizzare

Prova a cambiare alla fine. Alla stampa del numero primo (la stampa commentata, prova a scrivere:

                        if (n == 14680063L)
                            System.out.println ("Ho trovato 14680063");
                        else if (n > 14680063L)
                            System.out.println ("Peccato, 14680063 mi e' scappato");

Cosa ti aspetti che stampi?

[ingframin]

Ho capito dove sbaglio! Perdonami sottovento
Qui è spiegato con una chiarezza incredibile:
http://www.cs.hmc.edu/~oneill/papers/Sieve-JFP.pdf

[DoctorT]

Originariamente inviato da sottovento

Prova a cambiare alla fine. Alla stampa del numero primo (la stampa commentata, prova a scrivere:

                        if (n == 14680063L)
                            System.out.println ("Ho trovato 14680063");
                        else if (n > 14680063L)
                            System.out.println ("Peccato, 14680063 mi e' scappato");

Cosa ti aspetti che stampi?

oops mi è scappato.
dovevo scrivere: sieve.set(0, segment_size);
invece di: sieve.set(0, segment_size -1);

errore subdolo ... mi faceva perdere 30 primi sul primo milione.

[wingman87]

Originariamente inviato da ingframin

Ho capito dove sbaglio! Perdonami sottovento
Qui è spiegato con una chiarezza incredibile:
http://www.cs.hmc.edu/~oneill/papers/Sieve-JFP.pdf

Scusa, non ho più trovato il tempo di risponderti... grazie per il link.

[sottovento]

Originariamente inviato da ingframin

Ho capito dove sbaglio! Perdonami sottovento
Qui è spiegato con una chiarezza incredibile:
http://www.cs.hmc.edu/~oneill/papers/Sieve-JFP.pdf

Ci mancherebbe! Non hai offeso nessuno! Anzi, dai tuoi interventi c'e' sempre da imparare. Grazie per il link

[sottovento]

Originariamente inviato da DoctorT

oops mi è scappato.
dovevo scrivere: sieve.set(0, segment_size);
invece di: sieve.set(0, segment_size -1);

errore subdolo ... mi faceva perdere 30 primi sul primo milione.

Ottimizzazione: prova a sostituire i Vector con gli ArrayList
Secondo ma vai sotto il minuto

[ingframin]

Questo gia'somiglia di più al crivello di Eratostene.
Prendetelo con le pinze, non sono ancora sicuro che tutt i risultati siano corretti al 100%.
(Ci sto mettendo mano a lavoro tra una cosa e l'altra :P)
Il problema è che non ho abbastanza ram nel cassone che ho a lavoro...
Se aggiungo 4 miliardi di elementi al vettore il pc mi muore e se do più di 2GB di memoria all'heap della JVM il sistema comincia a fare swapping.
Se aggiungo solo i candidati ad essere primi (tutti i numeri che non sono divisibili per 2, per 3 e che soddisfano la condizione (n+1)%6 ==0 o (n-1)%6==0, da casa vi posto la fonte) e metto lheap a 3 GB (si swappa ma che devo fare...) perdo il riferimento degli indici e il meccanismo del crivello (che è quello di eliminare per indice e non facendo il test di divisibilità) non funziona più.
Diciamo che questo è un finto crivello.
Appena arrivo stasera (zita permettendo) ci provo sul mio pc dove ho 16GB di memoria e quindi ho meno problemi.

EDIT: Non funziona, sto cercando il bug

import java.util.*;

import java.sql.Timestamp;

public class EratosteneLong2{

	static HashSet<Long> buffer = new HashSet<Long>();
	
	public static void run(Long limit){
		boolean flag = false;
		Long y;
		buffer.add(2L);
		buffer.add(3L);
		for(Long x = 3L; x<limit; x+=2){
			
			if(x%3L!=0 || (x+1)%6==0 || (x-1)%6==0){	
				buffer.add(x);
				}
		}//ext for
		System.out.println("candidate list ready");
		Long L = Math.round(Math.sqrt(limit));
		
		for(Long x:buffer){
			if(x >= L){
				break;
			}
			//System.out.println(x);
			Iterator<Long> it = buffer.iterator();
			while(it.hasNext()){
				y = it.next();
				if(y<x) continue;
				
				if(x!=y && y%x==0){
					
					it.remove();
				}
			}
		}//ext for

	
	}//run method
	
	public static void main(String[] args){
		Long limit = 80000000000L;
			
		long now = System.currentTimeMillis();
		//run(limit);
		run(limit);
		long final_time = System.currentTimeMillis() - now;
		buffer.clear();
		System.out.println("time spent= "+final_time);
		
	
	}
}

[ingframin]

Dimenticavo che un'altra cosa figa che ho notato, cioé che in C++ (ci ho provato pure in C++...) uso un solo core della CPU mentre Java senza nessun accorgimento da parte mia li usa tutti e 4...
Alla faccia di chi dice che Java non è potente... :-O

[DoctorT]

Originariamente inviato da sottovento

Ottimizzazione: prova a sostituire i Vector con gli ArrayList
Secondo ma vai sotto il minuto

grande sottovento !! solo cambiando il nome della classe sono arrivato a 53 secondi ... tempo dimezzato

adesso prova a fare il crivello con un array al posto del bitset anche se occupa più memoria

[sottovento]

Originariamente inviato da DoctorT

grande sottovento !! solo cambiando il nome della classe sono arrivato a 53 secondi ... tempo dimezzato

adesso prova a fare il crivello con un array al posto del bitset anche se occupa più memoria

Stavo controllando la tua implementazione perche' secondo me' e' possibile fare un tentativo di parallelizzazione.

Il tuo thread e' semplicemente un ciclo, che si ripete per ogni intervallo in cui hai segmentato il crivello. L'idea sarebbe di lanciare una serie di thread, tutti uguali, i quali lavorano su parti diverse del crivello cercando di darsi il minimo fastidio

1) Il numero di thread da lanciare, secondo me, dovrebbe essere pari a

Runtime runtime = Runtime.getRuntime();
int nrOfProcessors = runtime.availableProcessors();
int threadCount = nrOfProcessors * 2;

2) La tua implementazione e' furba, e memorizza solo i numeri primi che hanno dei multipli all'interno del crivello. Gli altri li conta solamente
La cosa mi piace moltissimo, pero' il contatore dovra' essere condiviso fra i vari thread. Sincronizzare tutti i thread solo per incrementare un contatore e' davvero una orribile perdita di tempo. La soluzione quindi e' utilizzare un contatore di tipo java.util.concurrent.atomic.AtomicInteger. Questo ti garantisce che le letture/scritture del contatore siano atomiche con una sincronizzazione molto leggera (anzi, su molti sistemi non serve nemmeno la sincronizzazione).

3) Si potrebbe quindi fare una classe statica che metta a disposizione il contatore di cui sopra e che contenga una coda di tipo java.util.concurrent.ArrayBlockingQueue<>.
In questa coda, calcoli tutti gli intervalli da un thread produttore (quindi li' dentro ci puoi mettere solo l'indice di inizio).
I thread tutti uguali che fanno il calcolo saranno tutti bloccati in attesa su quella coda, estrarranno l'indice di inizio (eventualmente anche quello di fine) e vanno a lavorare su quel pezzo di crivello.

La sincronizzazione dovrebbe quindi essere ridotta al minimo, e c'e' da aspettarsi di una ulteriore riduzione del tempo di esecuzione. L'utilizzo della CPU dovrebbe quindi schizzare ben oltre il 30% della prima implementazione.

Cosa ne pensi?

[DoctorT]

molto interessante, ma con le mie scarse conoscenze di programmazione parallela non credo che riuscirei a tirare fuori niente di buono

[sottovento]

Originariamente inviato da DoctorT

molto interessante, ma con le mie scarse conoscenze di programmazione parallela non credo che riuscirei a tirare fuori niente di buono

Ci possiamo provare, ma credo avro' bisogno del tuo aiuto, perche' parte del tuo
codice non e' ancora chiarissimo.
L'idea e' quella di avere piu' thread, ciascuno dei quali lavora su una parte diversa del crivello. E' ovvio che quindi le singole parti devono essere indipendenti.
Purtroppo non e' proprio vero nel tuo codice, poiche' per ottimizzazione hai fatto in modo che alcune parti dipendano dai cicli precedenti, come per esempio:

                for (; s * s <= high; s++)
                {
                    if (is_prime[(int)s])
                    {
                        primes.add((int)s);
                        next.add((int)(s * s - low));
                    }
                }

ecc.
Queste parti devono quindi essere riscritte in modo da togliere la dipendenza dai cicli precedenti. E' la parte piu' critica poiche' si aggiunge del codice peggiorativo. Il fatto di poter usare il multithreading potrebbe non compensare questo peggioramento, o potrebbe non compensarlo sempre. Non sempre parallelizzare un algoritmo produce benefici.

Comunque, la class piu' semplice e' quella della coda + contatore:

class IntervalProvider
{
    public static AtomicInteger counter = new AtomicInteger(1);
    public static ArrayBlockingQueue<Long> queue = new ArrayBlockingQueue<>(50000);
    public static int numberOfConsumerThreads = Runtime.getRuntime().availableProcessors() * 2;     // Facciamo partire questi thread

    public static void doJob(long limit, int segment_size) throws InterruptedException
    {
        for (long low = 2; low <= limit; low += segment_size)
            queue.put(low);     // La put e' bloccante: se non c'e' spazio nella coda, restera' bloccato in attesa
        for (int i = 0; i < numberOfConsumerThreads; i++)
            queue.put(-1L);  // Questo e' il segnale di terminazione
    }
}

ArrayBlockingQueue<> puo' essere usata come una FIFO: chi scrive resta pero' bloccato in attesa nel caso che la FIFO sia piena. Analogamente chi legge resta bloccato nel caso la FIFO sia vuota.
doJob() verra' lanciato da un thread che scrive tutti gli intervalli in coda.
La coda verra' letta da una serie di thread, tutti uguali, che quindi lavoreranno sui diversi intervalli.

Passiamo alle note dolenti. Il tuo algoritmo deve essere peggiorato facendo in modo che i calcoli di un ciclo non dipendano dal precedente.
Probabilmente se ci dai un'occhiata potresti magari trovare un modo piu' efficente. Ad ogni modo, ho riscritto alcune parti della tua segmented_sieve():

    static void segmented_sieve(long limit, int segment_size) throws InterruptedException
    {
        int sqrt = (int) Math.sqrt((double)limit);
        //long s = 2;
        //long n = 3;

        BitSet sieve = new BitSet(segment_size);
        // togli il commento se vuoi i numeri stampati a video
        // System.out.println(2);
        //System.out.println("sqrtqrt");

        // crea un array di numeri primi < di sqrt(limit)
        boolean[] is_prime = new boolean[sqrt+1];
        for (int i=0;i<=sqrt;i++)
            is_prime[i] = true;
        for (int i = 2; i * i <= sqrt; i++)
        {
            if (is_prime[i])
                for (int j = i * i; j <= sqrt; j += i)
                    is_prime[j] = false;
        }
        // Downwind - al termine di questa operazione, i bit settati in questo array sono 
        // relativi ai numeri primi nell'intervallo 2..sqrt(limit)
        // Questo array non verra' piu' modificato nel corso del programma

        // Ora metto i primi cosi' trovati in lista
        List<Integer> primes = new ArrayList<> ();
        for (int i = 2; i <= sqrt; i++)
        {
            if (is_prime[i])
                primes.add(i);
        }
        // Controllo se e' piu' veloce con gli array, evitando quindi l'unboxing
        int[] vectPrimes = new int[primes.size()];
        for (int i =0; i < vectPrimes.length; i++)
            vectPrimes[i] = primes.get(i);

        // Ok, ora primes contiene i numeri primi che possono avere dei multipli nella sequenza.
        // Gli altri numeri primi non avranno nessun multiplo nel crivello

        // Leggo l'intervallo su cui lavorare
        int count = 0;

        long from = 0L;
        while ((from = IntervalProvider.queue.take()) != -1)     // La take() e' bloccante, quindi il thread stara' qui fintanto che c'e' un intervallo da lavorare
        {
            long to = from + segment_size;
            sieve.set(0, segment_size);		// Downwind - Setta tutti i bit dell'intervallo (estremo destro escluso, e questo dovrebbe essere un errore). Tolgo il -1
            //for (Integer prime : primes)
            for (int i = 0; i < vectPrimes.length; i++)
            {
                int pr = vectPrimes[i];
                if (pr > to)
                    break;
                // Cerco il primo multiplo di pr che sia maggiore o uguale di from
                long m = (from / (long)pr) * pr;
                if (m < from)
                    m += pr;
                else if (m >= to)
                    break;
                while (m < to)
                {
                    sieve.clear((int)(m-from));
                    m += pr;
                }
            }

            for (long i = from + 1; i < to; i += 2) // ricerca i bit con indice dispari
            {
                if (sieve.get((int)(i - from))) // is prime
                {
                    count++;
//                    int c = IntervalProvider.counter.incrementAndGet();
//                    if (c%1000000 == 0)
//                        System.out.println ("-->" + c);
                }
            }
        }
        IntervalProvider.counter.addAndGet(count);
    }

Come vedi, il ciclo for() e' stato sostituito da
from = IntervalProvider.queue.take()
vale a dire, leggo (in maniera bloccante) l'indice di partenza dell'intervallo dalla coda, invece che calcolarmelo. Se l'indice letto e' -1, allora il thread terminera' (infatti nella prima classe puoi vedere che vengono accodati tanti -1 quanti sono i thread consumatori).
Posso ora lanciare parecchie copie di questo metodo.
Il main() semplicemente crea i thread in questione. Praticamente e' quello che hai scritto tu con qualche piccola modifica:

    public static void main(String [] args) throws IOException
    {	
        final int sieve_size = 65536;
        final long n = 8000000000L;

        System.out.println("Attendere...");
        System.out.println("generazioneri primi ..."); 

        long inizio = System.nanoTime();

        ArrayList<Thread> vectThread = new ArrayList<>();
        // create the threads
        System.out.println ("Faccio partire " + IntervalProvider.numberOfConsumerThreads + " consumer thread");
        for (int i = 0; i < IntervalProvider.numberOfConsumerThreads; i++)
            vectThread.add(new Thread( () -> {
                try 
                {
                    Eratostene.segmented_sieve(n, sieve_size);
                } 
                catch (InterruptedException ex) 
                {
                    System.out.println("Interrupted!");
                    System.exit(0);
                }
            }));
        // Start the threads
        vectThread.stream().forEach(t -> t.start());

        // Spedisco tutti gli intervalli
        try
        {
            IntervalProvider.doJob(n, sieve_size);
        }
        catch (InterruptedException ex) 
        {
            System.out.println("Interrupted!");
            System.exit(0);
        }

        // Aspetto che finiscano tutti
        vectThread.stream().forEach(t -> {
            try
            {
                t.join();
            }
            catch (InterruptedException e)
            {
                System.out.println ("Interrupted!");
                System.exit(0);
            }
        });

        long fine = System.nanoTime();
        float time = (fine-inizio)/1000000;
        System.out.println(time + " ms");
        System.out.println ("Trovati " + IntervalProvider.counter.get() + " numeri primi");
    } // end main

Una volta lanciato, la CPU andra' al 100% (non e' piu' limitata al 25/30% come nell'algoritmo single-thread).
Sul mio laptop (quad core) ottengo un debole vantaggio rispetto al tuo algoritmo. E' quindi lecito aspettarsi che il vantaggio diventi piu' consistente su macchine con piu' processori/core mentre si riduca (o diventi peggiorativo) nel caso di macchine con meno core/uniprocessore.

[DoctorT]

non ho aggiornato eclipse a java 8 ci ho messo un bel po' a capire che significava ->

comunque per avere un piccolo miglioramento puoi inserire

if (pr*pr > to)
break;

al posto do

if (pr > to)
break;

[sottovento]

Originariamente inviato da DoctorT

non ho aggiornato eclipse a java 8 ci ho messo un bel po' a capire che significava ->

comunque per avere un piccolo miglioramento puoi inserire

if (pr*pr > to)
break;

al posto do

if (pr > to)
break;

Grande! Immagino che se lo guardi puoi ottimizzare ancora di piu'... e' codice tuo

[Scanca]

Mi scuso se rispindo solo ora ma ho pochissimo tempo.
Ringrazio "ingframin" perchè tramite la "L" sono riuscito ad andare avanti nel progetto e consegnarlo qualche settimana fa.

Il problema sta proprio nell'occupazione in memoria e nel tempo di esecuzione.
Il trucco quel'è?
Intanto usando la struttura dati BitSet si ottimizza molto l'occuazione in memoria, praticamente si sfrutta la posizione per dire se un numero è primo o meno mettendolo a true.
Un secondo trucchetto è quello di memorizzare i candidati primi non multipli di 2, 3 e 5. Così il mio BitSet sarà lungo circa 4.000.000.000 (quindi mi basta un intero ) e posso applicare il crivello di Eratostene su esso trovando in poco tempo quanti sono i numeri primi tra 2 e 8.000.000.000 in circa 30 minuti e occupando 250 Mbyte.

Spero di essere stato chiaro

[Scanca]

Volevo ringraziare tutti di cuore (me ne sono dimenticato prima) dell'aiuto e del tempo che mi avete dedicato. Grazie davvero!!!