Era un uggioso pomeriggio milanese dello scorso anno quando, procrastinando lo studio, venni colpito da un lampo di genio. “Cosa succederebbe se facessi un generatore di numeri casuali alimentato a banane?“. Preso dall’idea andai subito a raccontarla al coinquilino anch’esso elettronico.

Lui mi guardò in faccia e scoppiò a ridere. In quel momento capii che avevo nelle mie mani un grande progetto.

Prima di essere preso per pazzo: ha davvero un senso, e questo post è qui per spiegarvelo. Partiamo dalla radice del problema: i computer sono sistemi deterministici. Detto in modo più semplice, se gli diamo sempre gli stessi dati in ingresso, ci restituiscono sempre gli stessi valori in uscita (Non chiedetemi perché, se aprite quattro volte un programma nel corso della giornata, due volte si apre e due no; quello gli scienziati lo stanno ancora studiando). Che è esattamente quello che ci aspettiamo da un calcolatore. È chiaro da subito, però, che il determinismo e la casualità non sono grandi amici. Anzi, di per sé un computer, non è in grado di fare nulla di casuale.

Se chi legge è tra i fortunati ad aver scritto qualche programma in C, sarà sicuramente incappato in qualche chiamata della funzione rand() . Io stesso ne ho fatto largo uso nei post ([1][2]) riguardo al calcolo di pi greco col metodo Monte Carlo.

Per chi non avesse mai avuto la fortuna, rand() è la funzione che (con molta fantasia nel nome) permette di ottenere numeri casuali in un programma in C. Ma abbiamo appena finito di dire che non c’è nulla di casuale nei computer. Quindi?

Quindi ci hanno ingannato. Purtroppo siamo stati tutti ingannati. Mi prendo anche io la mia parte di responsabilità per averli chiamati, con leggerezza, numeri casuali nei post sul pi greco. Avremmo dovuto chiamarli numeri pseudocasuali. Che è il nome che gli informatici danno a quei set di numeri che hanno la distribuzione in linea con i numeri casuali, ma che in realtà casuali non sono.

Per capirci meglio, per essere definito casuale un set di numeri deve avere almeno (condizione necessaria ma non sufficiente) queste due caratteristiche:

I numeri non devono seguire un qualche pattern o disposizione, e su grandi set di dati La sequenza dei numeri non deve poter essere prevista in anticipo.

È evidente che la difficoltà delle macchine deterministiche è quella di rispondere al punto 2. Il punto 1 viene risolto con facilità, e da luogo a quelli che abbiamo appena chiamato numeri pseudocasuali. Numeri che rispettano la distribuzione uniforme, ma che in realtà non sono veramente casuali.

E insomma, le banane cosa c’entrano?!

Adesso ci arriviamo!

Quando è necessario fornire ad un computer dei veri numeri casuali, si usano dei sistemi hardware, detti appunto true random number generator (TRNG). Esisono tanti tipi di TRNG, che sfruttano diverse grandezze fisiche casuali e le convertono in informazioni digitali che vengono passate al computer.

I più comuni sfruttano fenomeni fisici come il rumore termico dei resistori, l’effetto valanga dei diodi ed altri effetti caotici. Altri sfruttano fenomeni quantistici più complessi come il rumore shot, il decadimento radioattivo o effetti fotonici.

La via percorribile utilizzando le banane è quella del decadimento radioattivo. Le banane infatti sono note per contenere molto potassio, e una piccola ma significativa percentuale del potassio presente in natura è radioattiva. Nello specifico parliamo dell’isotopo 40K, che costituisce lo 0,01% del potassio in natura.

Ora che l’affermazione del “generatore di numeri casuali alimentato a banane” ha un po’ più di senso o almeno un po’ di contesto, resta una domanda: cosa ce ne facciamo dei veri numeri casuali in un computer? La mia risposta vorrebbe essere “non mi interessa, voglio che il focus del mio progetto sia limitato al generarli”, ma sono sicuro che lascerebbe molti insoddisfatti.

Per cui proverò a rispondere anche a questa domanda: la crittografia. Questa è la principale ragione per cui vengono studiati i numeri casuali e il loro rapporto coi calcolatori. I numeri casuali vengono usati per generare le chiavi crittografiche, che sono l’unico fattore a determinare l’efficacia di un sistema di crittografia. Come afferma il Principio di Kerckhoffs: “la sicurezza di un crittosistema non deve dipendere dal tenere celato il crittoalgoritmo ma solo dal tenere celata la chiave”.

È chiaro che, se chi attacca può prevedere in qualche modo la chiave, questa non sarà più nascosta e saremo in presenza di un sistema vulnerabile. Un “buon random” quindi sta alla base di un buon sistema crittografico. [3]

Ok ma… se il random è casuale, come facciamo a distinguere un buon random da un cattivo random?

Per analizzare la qualità dei generatori di numeri casuali, esistono numerosi test statistici. Nel caso dei generatori di numeri casuali, sono disponibili numerose suite di test. Due delle più diffuse sono ent [4] e dieharder [5]. Il primo è stato progettato come test light per i generatori di numeri casuali a decadimento radioattivo, è molto semplice e veloce, necessita di pochi dati ma il risultato è puramente indicativo.

Dieharder di contro è una suite di test che viene considerata lo standard di riferimento per i generatori di numeri casuali, esegue test molto approfonditi ma necessita di diversi gigabyte di dati su cui eseguire i test.

Per ora ci concentreremo solo sui risultati ottenuti con ent.

Prepariamo i dati per eseguire un primo test con ent.

I dati vengono scritti sulla seriale dal generatore sulla seriale, possiamo salvarli su un file da console linux con cat /dev/ttyACM0 >> campionetesto.txt . Sfruttiamo il comando di redirect degli stream di bash in modalità “append“, in questo modo potremo interrompere l’acquisizione e riprenderla più tardi senza sovrascrivere il file.

Il campione raccolto nel corso di due giorni consta di 90628 numeri, uno per riga, compresi tra 0 e 65535. I numeri sono però salvati come file di testo ascii, mentre ent analizza file binari. Si può scrivere un brevissimo programma in C per convertirli in binario:

Fatto questo, possiamo eseguire per la prima volta il test di ent:

valerio@valerio ~/tests $ ./ent campione.txt Entropy = 7.997995 bits per byte. Optimum compression would reduce the size of this 181256 byte file by 0 percent. Chi square distribution for 181256 samples is 498.15, and randomly would exceed this value less than 0.01 percent of the times. Arithmetic mean value of data bytes is 127.4942 (127.5 = random).

Monte Carlo value for Pi is 3.138799695 (error 0.09 percent).

Serial correlation coefficient is 0.005408 (totally uncorrelated = 0.0).

Ent ci resitutisce diversi risultati:

Entropia: l’entropia è la quantità di “casualità” contenuta in una porzione di informazione. La teoria dell’informazione ci dice che la minima dimensione teoricamente ottenibile tramite compressione senza perdite di informazione , è rappresentata dal valore di entropia.

l’entropia è la quantità di “casualità” contenuta in una porzione di informazione. La teoria dell’informazione ci dice che la , è rappresentata dal valore di entropia. Distribuzione chi quadro: questo test viene usato per capire quanto bene la distribuzione dei nostri valori aderisce a una distribuzione teorica (nel nostro caso la distribuzione uniforme). Dal manuale di ent, questo valore deve essere più vicino possibile a 256, con valori percentuali compresi tra 10 e 90%

questo test viene usato per capire quanto bene la distribuzione dei nostri valori aderisce a una distribuzione teorica (nel nostro caso la distribuzione uniforme). Dal manuale di ent, questo valore deve essere più vicino possibile a 256, con valori percentuali compresi tra 10 e 90% Media aritmetica: La semplice media aritmetica dei nostri bit. Essendo i nostri valori compresi tra 0 e 255, dovrebbe essere circa uguale a 127.

La semplice media aritmetica dei nostri bit. Essendo i nostri valori compresi tra 0 e 255, dovrebbe essere circa uguale a 127. Valore di pi greco col metodo Monte Carlo: Il metodo dovrebbe essere ormai fin troppo conosciuto da chi segue i miei post, in questo contesto rappresenta però più un dato simpatico che un dato utile.

Il metodo dovrebbe essere ormai fin troppo conosciuto da chi segue i miei post, in questo contesto rappresenta però più un dato simpatico che un dato utile. Autocorrelazione: Serial correlation in inglese, rappresenta la dipendenza tra i valori della serie. Nel caso ottimo deve essere uguale a zero.

Da questo primo test tutti i valori risultano passabili, tranne il chi quadro. Prima di addentrarci nelle ragioni del cattivo risultato, è necessario comprendere che cosa è questo valore e come la casualità viene “estratta” dall’evento di una radiazione.

Partiamo dal capire cos’è e come funziona il chi-quadro (indicato anche con χ²). Si tratta di un valore utilizzato in statistica per verificare l’aderenza di un set di valori a una distribuzione teoricamente prevista. Affrontiamo prima l’utilizzo tradizionale in statistica, poi volgeremo lo sguardo all’uso del χ² per questa applicazione.

Dato un set di dati, chiamiamo frequenza il numero di volte che si presenta un determinato dato, e gradi di libertà il numero di possibili valori meno uno. Perché meno 1? Quest’operazione, che a prima vista può sembrare abbastanza strana, è in realtà molto logica. Supponiamo di volere analizzare il lancio di una moneta. Abbiamo due possibili outcome: testa e croce. Se ci pensiamo, però, la percentuale di volte che esce testa, è direttamente determinato dalla percentuale di volte che esce croce. Se consideriamo un evento con tre possibili outcome, la percentuale del terzo outcome è direttamente determinata dagli altri due. In questo senso, parliamo di gradi di libertà.

La formula del chi-quadro è:

dove oᵢ è la frequenza osservata e aᵢ è la frequenza attesa teoricamente.

Prendiamo come esempio il classico lancio di uno dado a sei facce. Il lancio del dado può risultare in sei possibili risultati, questo ci da 5 gradi di libertà. Immaginiamo di lanciare il dado 1000 volte, vogliamo verificare quella che in statistica prende il nome di ipotesi zero, ovvero verificare che, entro una certa probabilità, i nostri risultati siano veramente casuali.

Questi sono i dati che ottieniamo dal nostro esperimento:

Sommando i singoli contributi, otteniamo il valore di chi-quadro del nostro esperimento: 3.068. E ora? Prendiamoci un secondo per riflettere sui dati e sulla nostra formula, tenendo a mente che stiamo cercando di misurare, in questo caso, quanto i nostri dati aderiscono ad una distribuzione uniforme. Se tutto fosse perfetto (in un ipotetica dimensione in cui questo può succedere), avremmo le frequenze osservate oᵢ uguali alle frequenze teoriche aᵢ, per cui il nostro denominatore (oᵢ-aᵢ)²→0. Di conseguenza χ²→0. Sappiamo però che la perfezione non appartiene al nostro mondo, altrimenti non avremmo nemmeno bisogno della statistica. In altre parole, è estremamente improbabile che i nostri dati rispecchino esattamente la distribuzione teorica, un valore di chi quadro troppo vicino a zero ci deve far insospettire.

Dall’altro lato, più ci allontaniamo dalla distribuzione teorica, più il numeratore cresce, a parità di denominatore. Il che porta il valore di chi-quadro a crescere. Questo per l’uso normale del chi-quadro è un’ottima notizia, significa poter rifiutare l’ipotesi nulla, e quindi saper di stare lavorando con dati che non sono solo frutto del caso, ma che hanno un qualche significato. Per la nostra applicazione questa è invece una pessima notizia: significa che i nostri dati non sono distribuiti in maniera uniforme.

Questo vuol dire che, per la nostra applicazione, dovremo cercare una via di mezzo. Ok, ma numericamente? Il valore accettabile di chi-quadro varia in base ai gradi di libertà del sistema. Qui inizia la parte difficile, fatta di integrali poco amichevoli e di equazioni non risolvibili algebricamente [6]. Fortunatamente qualcuno ha già fatto il lavoro per noi, disponendo i risultati in tabelle comodamente reperibili su internet.

Le righe della tabella rappresentano i gradi di libertà del sistema, nel nostro caso del dado sono 5. Le colonne rappresentano i livelli di probabilità che il valore calcolato sia maggiore del valore tabulato. Esistono anche tabelle che indicano la probabilità che il valore calcolato sia minore, queste vengono chiamate left tail tables, mentre quella mostrata qua sopra è una right tail table. Questo perché in un caso viene considerata la parte destra del grafico, e nell’altro la parte sinistra. Nel nostro caso abbiamo χ²=3.068, che si colloca tra il 90% e il 25% dei casi. Questo è sufficiente per dire che non ci sono variazioni eccessive da un comportamento che possiamo classificare come random. [7]

Right tail del chi-quadro per d = 5 e chi-quadro = 3,068, da statdistributions.com

Tornando alle banane, prendiamo come percentuali di riferimento il 90% e il 10%. Purtroppo non sono riuscito a trovare tabelle per 255 gradi di libertà, fortunatamente però l’autore di ent mette a disposizione uno strumento per calcolare i valori per ogni probabilità [6].

Per 255 gradi di libertà (2⁸ -1) abbiamo i valori: P(χ²,90%)=226.5 e P(χ²,10%)=284.3. Dai test di ent sui valori registrati dal generatore, avevamo ottenuto un valore di 498.15, decisamente fuori dal range dell’accettabilità. La percentuale di probabilità restituita da ent era < 0.01%.

Ent, tramite l’opzione -c ci restituisce la frequenza di ogni byte, che è riportata qui di seguito, omessi i byte poco rilevanti:

valerio@valerio ~/tests $ ent logbinario.txt -c

Value Char Occurrences Fraction

0 739 0.004077

1 402 0.002218

2 1033 0.005699

3 674 0.003718

... ... ...

254 � 722 0.003983

255 � 691 0.003812 Total: 181256 1.000000 Entropy = 7.997995 bits per byte. Optimum compression would reduce the size

of this 181256 byte file by 0 percent. Chi square distribution for 181256 samples is 498.15, and randomly

would exceed this value less than 0.01 percent of the times. Arithmetic mean value of data bytes is 127.4942 (127.5 = random).

Monte Carlo value for Pi is 3.138799695 (error 0.09 percent).

Serial correlation coefficient is 0.005408 (totally uncorrelated = 0.0).

Si nota subito che il byte 1 ha significativamente meno count degli altri, e che il byte 2 ne ha molti di più. A uno sguardo attento, i count che sembrano “mancare” a 1, sono assegnati a 2. Questo è sicuramente il nostro problema.

Dopo un po’ di prove infruttuose, ho deciso di dividere i byte di posizione pari dai byte di posizione dispari. Questo perché per ogni numero a 16 bit (2 byte) generato, vengono creati due byte, uno di posizione pari e uno di posizione dispari. Ho praticamente diviso il file in Most Significant Byte e Least Significant Byte. Questi sono i risultati:

valerio@valerio ~/tests $ ent logbinariomsb.txt -c

Value Char Occurrences Fraction

0 358 0.003950

1 393 0.004336

2 347 0.003829

3 333 0.003674

... ... ...

254 � 374 0.004127

255 � 330 0.003641 Total: 90628 1.000000

...

Chi square distribution for 90628 samples is 240.88, and randomly would exceed this value 72.82 percent of the times.

... valerio@valerio ~/tests $ ent logbinariolsb.txt -c

Value Char Occurrences Fraction

0 381 0.004204

1 9 0.000099

2 686 0.007569

3 341 0.003763

... ... ...

254 � 348 0.003840

255 � 361 0.003983 Total: 90628 1.000000

...

Chi square distribution for 90628 samples is 891.23, and randomly would exceed this value less than 0.01 percent of the times.

...

Il gruppo dei MSB non riporta nessun problema significativo, mentre il gruppo dei LSB è dove risiede il problema.

Per capire da dove si origina questo problema, dobbiamo prima di tutto capire come i numeri sono generati internamente. Il tubo geiger, tramite un circuito di interfaccia, manda un segnale sul pin 2 (PB2/INT0) di arduino quando viene colpito da una radiazione. Il pin 2 è configurato in modo da generare un interrupt quando riceve un rising edge (quando il segnale passa da 0 a 1): attachInterrupt(digitalPinToInterrupt(2), randomCore, RISING); . L’interrupt chiamerà la funzione randomCore() , che è così definita:

void randomCore() {

Serial.println((unsigned int)(( micros() >> 2 ) & 0xFFFF));

return;

}

La funzione, quando chiamata, chiama a sua volta la funzione micros() di arduino. Questa funzione restituisce un numero a 32 bit che rappresenta il numero di microsecondi trascorsi dall’accensione del sistema. Essendo un numero unsigned a 32 bit, andrà in overflow dopo 4294.96 secondi, ovvero ogni circa 70 minuti. Essendo la velocità del microcontrollore insufficiente per ottenere un aggiornamento più preciso, la micros() si aggiorna a salti di 4 microsecondi, mantenendo sempre i due least significant bit a zero.

Per questo, eseguiamo uno shift del valore restituito dalla micros() di due bit a destra. In questo modo otteniamo un valore di 30 bit. Se usassimo anche i bit più signficativi, otterremmo dei numeri progressivi fino al successivo overflow del timer. Ogni numero sarebbe sicuramente maggiore del precedente e sicuramente minore del successivo, nell’arco dei 70 minuti necessari per l’overflow. Questo non è sicuramente random.

Conserviamo quindi solo i primi 16 byte del valore della micros() . Questo valore avrà overflow ogni 262144 microsecondi, rendendo estremamente improbabile il verificarsi della situazione di cui sopra. (Il contatore geiger registra di fondo circa 30 pulsazioni al minuto).

Il valore problematico è quindi il seguente, dove il MSB è irrilevante:

16 8 0

XXXX XXXX 0000 0001

Ha catturato la mia attenzione il fatto che questo valore si presenta ogni 4*²⁸ = 1024 microsecondi, cioè circa 1 millisecondo, ed è il valore successivo a quello che genera un overflow interrupt. La mia attenzione si è quindi spostata sul codice della funzione millis() del core di arduino [8] [9].

La funzione millis() del core funziona impostando il prescaler del TIMER0 a 64, che per un timer a 8 bit con un clock di 16MHz comporta un overflow del timer ogni 1.024 millisecondi. L’overflow del timer genera un interrupt, con vettore TIMER0_OVF . Se l’impulso del tubo geiger arriva contemporaneamente all’overflow del TIMER0 , avremo due interrupt concorrenti: TIMER0_OVF e INT0 . Questa situazione è gestita dal microcontrollore con un sistema di priorità degli interrupt, il cui ordine di priorità è indicato nel datasheet:

Priorità degli interrupt, dal datasheet dell’ATMega328/P [10]

L’overflow del TIMER0 ha priorità molto inferiore a quella dell’interrupt esterno, la mia ipotesi è quindi la seguente:

Caso 1: l’interrupt INT0 arriva contemporaneamente all’interrupt TIMER0_OVF . Avendo priorità maggiore viene eseguito prima l’interrupt esterno a scapito della millis() , intaccando la precisione della funzione ma senza produrre effetti visibili sui numeri generati; Caso 2: l’interrupt INT0 arriva al ciclo di clock successivo rispetto all’interrupt TIMER0_OVF . Essendo già trascorso un ciclo di clock, l’interrupt TIMER0_OVF sarà già in esecuzione. Quando l’esecuzione sarà terminata, micros() sarà già al valore 2, per cui il numero generato sarà registrato con valore 2.

È possibile che anche l’uso della seriale abbia un’influenza su questo ritardo, ma non ho approfondito. La soluzione a questo problema sarà esposta nel prossimo post, ed implica l’uso del TIMER1 per la generazione dei numeri, invece della funzione micros() di Arduino.

In chiusura, voglio ringraziare Mauro Mombelli e la community di StackExchange cryptography per il prezioso aiuto [11] [12].

Note, fonti e link:

[1] Follia computazionale: pi greco e il metodo Monte Carlo

[2] Follia computazionale: pi greco, i thread e la follia collettiva

[3] inside secure — The importance of true randomness in cryptography — white paper

[4] Ent: A Pseudorandom Number Sequence Test Program

[5] Dieharder: A Random Number Test Suite — Robert G. Brown

[6] Fourmilab — chi-squared calculator

[7] How to evaluate chi-squared results

[8] Millis in wiring.c su github

[9] µC eXperiment blog — Examination of the Arduino millis() Function

[10] Datasheet dell’ATMega328/P

[11] StackExchange Cryptography — Random number generator issue

[12] StackExchange Cryptography — How to evaluate chi squared result?

[13] AVRBeginners — Timers

Originally published at www.valerionappi.it on March 25, 2018.