Integrator

Nonostante la rassicurante immagine natalizia con gli auguri da parte dell'ex governatore della California qui a fianco, il post di questo mese non è dedicato a un nuovo film in uscita nelle sale, ma al fantastico mondo degli integratori, ovvero i motori di un algoritmo di dinamica molecolare, dei quali mai avrei pensato di poter scrivere. Quando ho deciso di mettere su questo blog, circa un anno fa, pensavo fosse relativamente facile spiegare quello che, da circa 10 anni a questa parte, ho scelto di portare avanti come mia principale attività di ricerca. La verità è che non è per niente facile, soprattutto senza l'ausilio delle mie amate equazioni, tuttavia è molto istruttivo cercare di farlo e, a tratti, anche piuttosto divertente, almeno per me. Per chi legge, non so. Tentiamo dunque l'impresa e fatemi sapere come è andata!

La biofisica computazionale comprende diversi metodi di simulazione: tra questi la dinamica molecolare emerge e si impone quasi prepotentemente per diversi fattori.

Il primo è la possibilità di descrivere il dettaglio atomico degli elementi che compongono le molecole biologiche, il secondo è la relativa facilità con cui è possibile scrivere un codice di simulazione che sostanzialmente non fa altro che risolvere la celeberrima equazione di Newton qui a lato, determinando le forze che agiscono su ciascun atomo e, da queste, le accelerazioni corrispondenti e quindi le traiettorie di tutti gli atomi del sistema. Ho scritto "relativa facilità" perché sembra semplice, ma non lo è. Primo, perché le equazioni da risolvere sono tantissime; secondo, perché non si possono risolvere con carta e penna; terzo, perché anche disponendo del più potente computer, quest'ultimo non è proprio in grado di svolgere un calcolo in modo esatto. I numeri su un computer sono infatti rappresentati in forma binaria come sequenza di 1 e 0, con gli 1 che rappresentano potenze del numero 2, ad esempio il numero 110101 sta per 1x32+1x16+0x8+1x4+0x2+1x1 cioè 53. Qualunque numero non intero sarà dunque sempre affetto da errori di arrotondamento per il semplice fatto di essere stato "costretto" su un calcolatore binario. Risolvere l'equazione di Newton significa dunque trovare una procedura per determinare le posizioni degli atomi nel futuro a partire da quelle nel presente, conoscendo tutte le forze e le velocità, con la migliore approssimazione possibile.

Questo è il compito degli integratori, il cui nome discende proprio dal fatto che "integrano" l'equazione di Newton che contiene l'accelerazione, ovvero una derivata seconda della posizione rispetto al tempo: l'integrale infatti altro non è che l'operazione inversa rispetto alla derivata. Chi di noi si è cimentato con gli integrali ricorderà grandi calcoli che dipendevano fortemente dalla funzione che andava integrata: alcune erano piuttosto semplici, altre decisamente più complicate e per alcune ci si metteva le mani nei capelli, quando c'erano ancora! Tutto questo perché facevamo i calcoli con carta e penna: per un computer invece non importa quanto la funzione sia complicata perché le operazioni saranno sempre le stesse e saranno tutte operazioni numeriche. Il computer è infatti una macchina non troppo intelligente che fa esattamente ciò che le chiediamo di fare e lo fa milioni, miliardi di volte senza protestare, ma solo consumando energia elettrica. L'equazione di Newton in realtà non è proprio scritta così come l'ho riportata, ma in una forma un po' più elegante (sul senso estetico che hanno i fisici per le equazioni e non per i capi di abbigliamento si potrebbe scrivere un intero numero di Vanity Fair con reportage fotografico). Queste equazioni "eleganti" della meccanica classica sono le equazioni di Hamilton-Jacobi o anche solo Hamilton o equazioni della meccanica hamiltoniana. Le equazioni hamiltoniane hanno un'importante proprietà: conservano invariato il valore di una quantità che, con la nota fantasia dei fisici, si chiama Hamiltoniano.

L'Hamiltoniano, di fatto, è l'energia totale del sistema, quella che si conserva sempre quando ci sono in gioco forze "conservative" (sempre per la nota fantasia dei fisici) e non si conserva quando ci sono forze cosiddette dissipative (come ad esempio gli attriti). Infatti le equazioni di Hamilton-Jacobi conservano l'energia totale, che resta a un valore costante. O meglio, dovrebbero conservarla costante perché ci sono sempre in agguato i famosi errori di arrotondamento dovuti al fatto che un computer continua a rappresentare in forma binaria (e quindi con numeri interi) numeri che interi proprio non sono. L'algoritmo che conserva meglio il valore dell'energia, a meno degli errori di arrotondamento, si chiama algoritmo di Verlet e fu ideato da (indovinate chi?) Verlet nel 1957 e successivamente riaggiustato con varianti del tutto equivalenti. In realtà l'algoritmo di Verlet non è il più preciso: ci sarebbero algoritmi ben più complicati che permettono di integrare le equazioni di Hamilton-Jacobi con maggiore precisione. Eppure resta il più utilizzato in tutte le simulazioni di dinamica molecolare. Il motivo? Pur non essendo estremamente preciso, ha una serie di proprietà tra cui quella di restare sempre vicino alla vera traiettoria di un sistema.

In altre parole, sì, sappiamo che fa errori, ma questi errori non sono mai troppo grandi e restiamo sempre sufficientemente vicini alla vera traiettoria, quella che otterremmo se fossimo in grado di fare tutti quegli integrali a mano. Paradossalmente un algoritmo più preciso del Verlet non ha questa proprietà: può essere molto preciso nel breve termine e prendere grossissime cantonate al passare del tempo. Un po' come nella vita: si tratta di decidere se fare un piccolo passo incerto restando sempre abbastanza vicini alla perfezione o cercare a tutti i costi di essere perfetti e finire per perdere la strada. Quanti insegnamenti si traggono dalla fisica computazionale!
L'algoritmo di Verlet non conserva solo l'energia, ma conserva anche altre due quantità. La prima è il volume: resta sempre lo stesso durante tutta la simulazione. La seconda è il numero di particelle: non si prevedono creazioni e distruzioni di particelle durante una simulazione, quindi il numero di atomi non può cambiare. Se però vogliamo confrontare i risultati di una simulazione con una situazione sperimentale, ci rendiamo conto che siamo davvero lontani. In un esperimento infatti conservare l'energia non è per niente banale: di solito quello che facciamo è mantenere costante la temperatura che è legata all'energia cinetica delle particelle (cioè quella associata al fatto che hanno una velocità: è zero se sono ferme), ma non è l'energia totale.

Inoltre neanche il volume è mantenuto costante, se non in una pentola a pressione (che infatti dovrebbe chiamarsi pentola a volume costante): normalmente teniamo ferma la pressione e lasciamo che il volume vari come meglio crede. Infine non è propriamente vero che manteniamo costante il numero di particelle, anche se possiamo assumere con buona approssimazione che questo vari poco. Il nostro mitico Verlet va quindi modificato per ottenere almeno un algoritmo che mantenga costante la temperatura e la pressione: sul numero variabile di particelle si sta lavorando, ma al momento almeno io preferisco lavorare in questi casi con un metodo diverso dalla dinamica molecolare, ovvero il metodo Monte Carlo.
Cosa possiamo inserire per mantenere costante la temperatura? Esattamente lo stesso dispositivo che abbiamo in casa: un termostato, ovvero un pezzo di programma al computer che calcola la temperatura, la confronta con una temperatura di riferimento e apporta le necessarie modifiche alle velocità delle particelle per mantenere la temperatura costante.

Esistono diverse possibilità e quindi diversi termostati, con vantaggi e svantaggi, ma tutti (e sì, proprio tutti) hanno in comune la caratteristica che, a differenza del Verlet, fanno un po' innamorare le nostre traiettorie facendo perdere (almeno in parte) la cognizione del tempo! L'algoritmo di Verlet, infatti, pur generando una traiettoria approssimativa, genera una traiettoria che ha anche un senso temporale, come se fosse un film proiettato in avanti o all'indietro. Se inseriamo un termostato questa proprietà viene, almeno in parte, persa: non siamo più sicuri che la traiettoria rappresenti davvero l'evoluzione del sistema nel tempo. Possiamo soltanto affermare che le varie fotografie che otteniamo a istanti successivi sono effettivamente fotografie di uno stato corrispondente alle condizioni della simulazione e che quindi il sistema che stiamo simulando prima o poi si comporterà come abbiamo previsto dalla simulazione, ma è difficile dire se lo farà prima o poi.
Avendo superato questo dolore, risulta ora più facile digerire il passo successivo, ovvero il fatto che per tenere costante la pressione, dobbiamo far variare il volume: inseriamo un barostato che calcola la pressione ad ogni istante della simulazione (tenendo presente che quell'istante non ha molto a che vedere con il tempo perché ci siamo già innamorati con il termostato) e la confronta con quella che vorremmo mantenere (di solito una atmosfera, esattamente la pressione alla quale normalmente viviamo) e fa le opportune modifiche al volume perché questa resti costante, gonfiando o sgonfiando il sistema che stiamo simulando come se fosse un palloncino.

In definitiva, più gli integratori cercano di simulare le condizioni di un esperimento, più si allontanano da quella rassicurante eleganza della meccanica hamiltoniana: come in amore, perdiamo sì la cognizione del tempo ma impariamo ad accontentarci di quello che riusciamo a generare e forse ci confrontiamo meglio con la realtà. Forse, eh?