Eratostene in Haskell

Oggi vediamo come si possa calcolare la sequenza dei primi m numeri primi utilizzando il crivello di Eratostene con Haskell. In realtà, faremo un post molto breve che non è altro che una piccola esercitazione di programmazione in Haskell. Se davvero voleste scrivere un algoritmo super efficiente per calcolare una sequenza di numeri primi, dovete guardarvi questo post e perdere un bel po’ di tempo a studiarvi tutte le varianti proposte.

Ma poiché il tempo a nostra disposizione è quello che è; e oggi è una bella giornata e c’è il sole; i figli richiedono la mia presenza; per tutti questi motivi sarò breve e commenterò solo la prima soluzione semplice per realizzare un crivello di Eratostene.

L’algoritmo è semplicissimo: prendete la lista dei primi m numeri naturali (a partire da 2), e toglieteci tutti i multipli di 2; poi prendete il secondo elemento della lista, cioè il 3, e togliete tutti i multipli di 3; poi prendete il terzo elemento della lista, cioè il 5, e togliete tutti multipli di 5; e così via. Ne ha parlato anche il nostro amico zar in un suo post.

E’ molto semplice scrivere un algoritmo così in Haskell. Supponiamo di saper già fare la sottrazione tra due liste con una funzione che si chiama `minus`. Allora, i primi m primi si calcolano così:

primesTo m = 2 : eratostene [3,5..m] where
  eratostene []     = []
  eratostene (p:xs) = p : eratostene (xs `minus` [p, p+2*p..m])

Vediamo un po’ di capirci qualcosa. Definiamo al solito una funzione per calcolare la lista dei primi fino a m. Tale funzione è pari alla lista contenente il numero 2 in testa, seguita dalla lista prodotta da un’altra funzione locale chiamata eratostene. Alla quale passiamo la lista dei numeri dispari fino a m.

La funzione eratostene è ricorsiva: se gli diamo una lista vuota, ritorna una lista vuota. Se gli diamo una lista con in testa il numero p (che assumiamo essere già primo) e una lista xs, ci da come risultato una lista che contiene in testa il numero p, e in coda la lista ottenuta cancellando da xs la lista contenente tutti i multipli dispari di p fino a m (scritta come [p, p+2*p..m]).

Sintentico ed elegante come al solito! Vi torna? Se avete qualche dubbio scrivetelo nei commenti

Rimane da capire come si realizza la funzione `minus` in maniera elegante e soprattutto efficiente. Fare la sottrazione non è difficile in generale, il difficile è farlo in maniera efficiente. In effetti, la prima semplice soluzione sarebbe di prendere ogni elemento da xs e vedere se si trova nella lista [p, p+2*p..m]. Ma questo significa che per ogni elemento di xs dobbiamo scorrere tutta la lista [p, p+2*p..m]. Quindi, la complessità sarebbe pari al prodotto delle dimensioni delle due liste. Non c’è una maniera migliore di farlo?

La risposta è sì se assumiamo che le due liste siano ordinate in ordine crescente. Ed ecco l’algoritmo:

minus (x:xs) (y:ys) = case (compare x y) of
  LT -> x : minus xs (y:ys)
  EQ -> minus xs ys
  GT -> minus (x:xs) ys

minus xs _ = xs

Immanzitutto, stiamo utilizzando un nuovo costrutto che si chiama case of, ed è piuttosto simile allo switch case del C. Qui l’argomento del case è il risultato della funzione compare, che è una delle tre costanti LT (less than), EQ (equal) o GT (greater than). In particolare, se x è minore di y, allora il risultato di minus è x seguito da xs `minus` ys. Se invece x e y sono uguali, allora bisogna rimuovere x, e quindi il risultato è semplicemente xs `minus` ys. Infine, se x è maggiore di y, allora devo semplicemente scorrere la seconda lista e quindi il risultato è (x:xs) `minus` ys.

L’ultimo pezzo serve nel caso in cui la seconda lista sia vuota. In quel caso, il risultato della “sottrazione” tra liste è semplicemente la prima lista. La complessità di questo algoritmo è pari, nel caso peggiore, al numero totale di elementi nelle due liste.

Ed ecco il codice complessivo:

module Main where
import System.Environment

main :: IO()
main = do args <- getArgs
          putStrLn ( show $ primesTo (read $ args!!0) )

primesTo m = 2 : eratostene [3,5..m] where
  eratostene []     = []
  eratostene (p:xs) = p : eratostene (xs `minus` [p, p+2*p..m])

minus (x:xs) (y:ys) = case (compare x y) of
  LT -> x : minus xs (y:ys)
  EQ ->     minus xs ys
  GT ->     minus (x:xs) ys

minus xs _ = xs

Ed eccoci arrivati in fondo. Con questo post si conclude una prima parte di questo minicorso sui linguaggi funzionali e sul linguaggio Haskell. Per completare ci servono: strutture dati e programmazione sequenziale con Haskell (sì, si può!) e con le monadi. Alla prossima!

Annunci
Post a comment or leave a trackback: Trackback URL.

Commenti

  • zar  On 4 marzo 2012 at 22:57

    Bello! E’ anche veloce?

    • glipari  On 4 marzo 2012 at 23:03

      A sentir loro:

      http://www.haskell.org/haskellwiki/Prime_numbers

      non particolarmente, dato che ne propongono un numero piuttosto elevato di varianti, sempre più complicate. Però era quello più semplice da spiegare (e da studiare, soprattutto!). Bisognerebbe aver il tempo di fare qualche analisi comparativa. Prima o poi mi metto a “tradurre” anche quello di Atkin, se mi riesce di trovare un po’ di tempo per farlo… 😦

      • zar  On 4 marzo 2012 at 23:06

        Uh, quello di Atkin è difficile…

  • compscient  On 5 marzo 2012 at 08:35

    Guardando velocemente i tempi di esecuzione sembrerebbe quasi peggio delle divisioni successive. (E, inoltre, siccome ad ogni passo ricorsivo crea una nuova lista, credo che la memoria occupata sia davvero troppa, anche se non ho controllato)

    D’altra parte, le liste immutabili non sono il tipo di dato più adatto ad un algoritmo di questo tipo.

    • glipari  On 5 marzo 2012 at 08:44

      Sulla velocità hai probabilmente ragione. Sulla memoria bisognerebbe controllare, se non sbaglio Haskell ha un meccanismo per fare “caching” delle liste immutabili. Comunque, sul link riportato qui su ci sono anche delle versioni con liste mutabili, sembra siano un po’ più veloci (a sentire loro), ma non ho fatto alcuna prova.

      • compscient  On 5 marzo 2012 at 18:39

        Credo che un moderato miglioramento nei tempi di esecuzione si possa ottenere molto facilmente aggiungendo una condizione del tipo p*p > m in erastotene (e fra l’altro, sarebbe anche più rispondente alle reali specifiche del crivello).

      • compscient  On 5 marzo 2012 at 18:40

        Ah, e per quanto riguarda la memoria, in effetti non ne usa moltissima. Pensavo molto peggio, ma probabilmente per un linguaggio di questo tipo il garbage collector risulta molto più semplice (e migliore).

Rispondi

Inserisci i tuoi dati qui sotto o clicca su un'icona per effettuare l'accesso:

Logo WordPress.com

Stai commentando usando il tuo account WordPress.com. Chiudi sessione / Modifica )

Foto Twitter

Stai commentando usando il tuo account Twitter. Chiudi sessione / Modifica )

Foto di Facebook

Stai commentando usando il tuo account Facebook. Chiudi sessione / Modifica )

Google+ photo

Stai commentando usando il tuo account Google+. Chiudi sessione / Modifica )

Connessione a %s...

%d blogger hanno fatto clic su Mi Piace per questo: