Lattice QCD per principianti

Digressione numerica - 1/2

Sì, va bene, tutto fico, tutto bello, ma alla fine noi vogliamo sapere come si fanno le simulazioni numeriche. Che ci faccio con tutto questo apparato teorico se poi non so farci i conti? Ok, come volete. Però dobbiamo ripartire dal modello di Ising. Abbiamo visto (lo abbiamo visto nel caso della magnetizzazione) che possiamo scrivere il valore di aspettazione di un certo osservabile $A$ in questo modo: $$\langle A \rangle = \frac{1}{Z}\sum_{\cal C} A[{\cal C}] \text{e}^{-\beta H[{\cal C}]}$$ Ricapitoliamo:

${\cal C}$ è una generica configurazione di spin;
$\sum_{\cal C}$ è la somma su tutte le possibili configurazioni;
$A[{\cal C}]$ è il valore che assume l'osservabile $A$ nella configurazione ${\cal C}$;
$H[{\cal C}]$ è il valore dell'energia nella stessa configurazione;
$Z = \sum_{\cal C} \text{e}^{-\beta H[{\cal C}]}$ è la funzione di partizione canonica del sistema;
$\beta = 1/T$ (ci siamo messi in unità $k_{\text{B}} = 1$) in cui $T$ è la temperatura del sistema.

Ora siamo pronti per

L'approccio ingenuo ingenuo

Il mio piccolo lettore romanista dirà subito: "Be', ma che ce vò, ciabbiàmo supercomputer fichissimi, alla fine quelli sò spin, valgono $\pm 1$, generiamo tutte le configurazioni, e che ce vorrà mai, calcoliamo $A$ e $H$ per ogni configurazione, fàmo le somme, il rapporto, ecco fatto!" Sì, ok. 'Spetta. Un sistema di Ising con $N$ spin ha un numero totale di configurazioni pari a $2^N$. Tenete conto che per calcolare numericamente, con una buona precisione statistica e sistematica, gli esponenti critici del modello di Ising in $d=3$ è stato necessario considerare sistemi con un numero di spin pari a circa $10^9$ ($1000$ punti per lato per $3$ dimensioni). Un sistema del genere ha un numero di configurazioni pari a $\exp(10^9\ln 2) \simeq 10^{300000000}$ (se non ho sbagliato i conti, il che è sempre possibile; controllate). Immaginiamo di avere un supercomputer veramente fantastico, che genera una configurazione, aggiorna le somme eccetera nel tempo di $10^{-15}$ secondi. Sono veramente ottimista. Ciò significa che in un secondo otteniamo $10^{15}$ configurazioni ("mammagara!" esclama il piccolo simulatore numerico che è in me), e a sua volta implica che per ottenere tutte le configurazioni dobbiamo aspettare solo per un tempo pari a circa $10^{299999985}$ secondi. "Be', no, scusa, io a una cert'ora ciò judo".

Ok, non funziona: possiamo pensare però a

L'approccio ingenuo

Appena tornato da judo, il mio piccolo lettore romanista mi dice: "Vabe', ho capito che come dicevo prima non se pò fà. Però ciò la soluzione in tasca. Fàmo come i sondaggi, che se intervistano 1000 italiani i sondaggisti ciànno i risultati delle elezioni col margine del 3%... E vabe', ciavrémo un errore statistico, lo calcoleremo: estraiamo un po' di configurazioni a caso... 'Na mijonata ne potremo estrarre, no? Ma fàmo pure de più, 'na mijardata... Poi le pesàmo cor peso de Boltzmann, fàmo le somme, i rapporti... Nooo? Mo che c'è che non va..."

Per spiegargli cosa c'è che non va mi tocca riscrivere la funzione di partizione canonica in un altro modo e fare almeno un disegnetto. Notiamo che nella nostra definizione originaria, la funzione di partizione è la somma su tutte le possibili configurazioni di spin. Possiamo cambiare punto di vista e scriverla come una somma su tutti i possibili livelli di energia del sistema: $$ Z = \sum_{E} g(E) \text{e}^{-\beta E} $$ $g(E)$ è la cosiddetta densità degli stati: conta quante configurazioni diverse di spin hanno la stessa energia $E$, e stati con la stessa energia sono equiprobabili. Questo è il famoso postulato delle uguali probabilità a priori. Se avete seguito Meccanica Statistica lo sapete già, e sennò lo sapete adesso.

Notate che, a meno di un fattore $k_{\text{B}}$, che però abbiamo posto uguale a $1$, il logaritmo di $g(E)$ è per definizione l'entropia microcanonica di un sistema con energia fissata $E$: $S(E) = \log g(E)$. Quindi potremmo scrivere $$ Z = \sum_{E}\text{e}^{S(E)-\beta E}$$ Ma guarda tu che combinazione, all'esponente compare una cosa che assomiglia tanto all'energia libera... (non è un caso, se proprio ci tenete a saperlo).

Ora chiediamoci quale sia la probabilità $P(E)$ che una data configurazione abbia proprio energia $E$: la risposta è facile: $$ P(E) = g(E) \text{e}^{-\beta E} $$ Se considerassimo un sistema con energia continua dovremmo scrivere $ P(E)\text{d}E = g(E) \text{e}^{-\beta E}\text{d}E $ per ovvi motivi. $g(E)$ è una funzione che cresce circa a potenza con $N$ (numero di spin nel sistema), mentre l'esponenziale negativo ammazza tutto. Il risultato globale, almeno in senso qualitativo, è quello mostrato nella figura che segue:

Ok, la figura fa un po' cacare, mi sono anche dimenticato di colorare una freccia, le ombre sono sbagliate e le funzioni fanno schifo; però spero di aver reso l'idea. Nel limite termodinamico, cioè per numero di spin $N$ molto molto grande, tendente a infinito, l'esponenziale scende molto molto mooolto velocemente ($E$ è una quantità estensiva, quindi cresce linearmente con $N$), la $g(E)$ sale molto molto velocemente e la $P(E)$, almeno intorno al massimo, va verso una gaussiana talmente piccata da diventare... una delta di Dirac.

Da questo piccolo esercizio impariamo che estrarre configurazioni a caso è una perdita di tempo. Noi vorremmo estrarre configurazioni significative, ossia vicino al picco, perché sono quelle che dànno un contributo rilevante. Se estraiamo configurazioni completamente a caso, con distribuzione di probabilità piatta, la probabilità di beccarne una vicino al picco è circa zero, e diventa esattamente zero nel limite termodinamico: la stragrande maggioranza delle configurazioni che estrarremmo sarebbe assolutamente inutile, e il nostro valore d'aspettazione non avrebbe senso alcuno. E noi siamo proprio interessati a fare simulazioni col maggior numero di spin che sia possibile, perché vogliamo evitare effetti di volume finito e poi in fondo perché la Meccanica Statistica funziona appunto nel limite termodinamico. Quindi mi dispiace, piccolo e simpatico amico romanista, ma neanche questo approccio funziona.

"E vabe', mai 'na gioia..."

L'approccio furbo

No, non temere, piccolo e simpatico amico romanista. Abbiamo una soluzione. Questa soluzione si chiama importance sampling. Significa che dobbiamo continuare a perseguire l'idea di estrarre configurazioni a caso, ma non completamente a caso, cioè, non con distribuzione piatta. Le vogliamo estrarre proprio seguendo la distribuzione di Boltzmann. Cioè vogliamo che la probabilità di estrarre una configurazione con energia $E$ sia proporzionale al peso di Boltzmann. Una formula vale più di mille parole: $$ p(E) \propto \text{e}^{-\beta E} $$ Immaginate di riuscire a fare una cosa del genere: non sappiamo ancora come, anche se lo scopriremo presto, ma per il momento immaginate di riuscire a farlo. Estraete $M$ configurazioni ${\cal C}_i,\;\,i = 1,\cdots M$, per ogni configurazione calcolate il valore dell'osservabile $A$, quindi avremo $A_i,\;\,i = 1,\cdots M$ e il valore d'aspettazione del vostro osservabile sarà semplicemente $$ \langle A \rangle = \frac{1}{M} \sum_{i=1}^{M}A_i $$ E se proprio vi prude di conoscere l'errore statistico, ammettendo per un momento che le vostre configurazioni siano statisticamente non-correlate (non so se tornerò su questo argomento), potete sempre calcolare $$ \sigma_A = \sqrt{\frac{\langle A^2\rangle - \langle A \rangle^2}{M-1}} $$. "Ammazza che ficata... Maaa... ciò un po' de domande". Chiedi pure, piccolo e simpatico amico romanista.

"Scusa, ma ner valore d'aspettazione $\langle A \rangle$ che hai scritto prima, che fine ha fatto er peso de Boltzmann".
È implicito nel fatto che le configurazioni sono estratte proprio seguendo la distribuzione di Boltzmann, non c'è bisogno di rimettercelo dentro.
"Ah... Eeee... la $g(E)$?".
Configurazioni con energia uguale sono equiprobabili. E poi estraggo configurazioni, quindi vedi la prima definizione del valore d'aspettazione, la $g((E)$ non c'è.
"Hai capito... Ok, ma $Z$? La funzione de partizione che fine ha fatto."
Si cancella tra numeratore e denominatore. In qualche senso è inclusa nel prezzo.
"Anvedi... Sì ma ancora non m'hai spiegato come faccio a estrarre le configurazioni in questo modo..."
Questo lo vediamo nella prossima puntata, perché sennò poi la pagina diventa troppo lunga e i lettori s'annoiano. E io ho bisogno di andare a fare una passeggiata. Comunque serve una cosa chiamata "catena di Markov".

Vai alla puntata precedente

Cerca nel blog

Azione Parallela