Potenziali interatomici da machine learning: accuratezza DFT a costo classico

La dinamica molecolare (MD) è uno strumento potentissimo per studiare i materiali a livello atomico: diffusione, transizioni di fase, meccaniche dei difetti. Il suo tallone d’Achille è il campo di forze: i potenziali classici (Lennard-Jones, Tersoff, MEAM) sono veloci ma poco accurati per legami complessi; la DFT è precisa ma permette solo sistemi di poche centinaia di atomi per pochi picosecondi. I potenziali interatomici da machine learning (MLIP o NNP) rompono questo compromesso.

Questo articolo spiega come si costruisce un MLIP, perché può avvicinarsi all’accuratezza della DFT a costo di un potenziale classico, e quali applicazioni sono già state dimostrate.

Il problema dei campi di forze classici

Un potenziale classico descrive l’energia in funzione delle posizioni atomiche con una forma funzionale fissa (somme di termini a due corpi, tre corpi, ecc.) con parametri ottenuti fittando dati di equilibrio. Il problema è l’estrapolazione: il potenziale è calibrato vicino all’equilibrio e fallisce in stati lontani da esso — stati di transizione, superfici, difetti, alta pressione. I MLIP invece apprendono l’energia in tutto lo spazio di configurazione rappresentato nel training set, senza presumere una forma funzionale.

L’approccio: energia come somma di contributi atomici

L’idea fondamentale, introdotta da Behler e Parrinello nel 2007, è di scrivere l’energia totale come somma di contributi atomici, ciascuno funzione solo dell’ambiente locale (entro un raggio di cutoff R_c di solito 6–8 Å):

E_tot = ∑_i ε(d_i)    con d_i = fingerprint dell’atomo i

Il fingerprint d_i descrive l’ambiente dell’atomo i (distanze, angoli, tipi degli atomi vicini). La rete neurale — o il kernel — ε impara la funzione che mappa d_i in un’energia atomica. Sommando i contributi di tutti gli atomi si ottiene l’energia totale; le forze si ottengono come gradiente negativo: F_i = −∂E_tot/∂R_i.

Quattro passi per costruire un MLIP

Il flusso di lavoro standard comprende quattro fasi. 1. Dataset di riferimento: si generano configurazioni diverse (bulk cristallino, superfici, difetti, strutture amorfizzate, stati di alta pressione) e si calcola l’energia e le forze atomiche di ciascuna con DFT. Per un potenziale del carbonio si includono tipicamente bulk privo di difetti, superfici, vacanze, cluster isolati, grani di confine e dislocazioni. 2. Fingerprinting: ogni configurazione viene convertita in fingerprint invarianti. 3. Selezione del training set: si evitano ridondanze raggruppando i fingerprint in cluster (k-means) o per ampiezza della forza (force-binning) e si campiona da ciascun gruppo. 4. Regressione: si addestra il modello (rete neurale o kernel) a predire energia e forze.

Tre architetture principali: kernel, lineare, rete neurale

La kernel ridge regression (KRR) predice la forza su un atomo i come combinazione ponderata di kernel di similarità tra il suo fingerprint e quelli del training set:

F_i = ∑_j=1^N_t a_j exp[−(d_ij/σ)²/2]   (kernel ridge regression)

dove d_ij è la distanza nel fingerprint space tra la configurazione i e il riferimento j, e σ è il parametro di lunghezza di scala. KRR raggiunge la precisione chimica (errore < 0.1 eV/Å sulle forze) per metalli come Al, Cu, Ti, W, Si e C. Il modello lineare è più veloce: F = V⋅w, dove V è la matrice dei fingerprint e w il vettore dei pesi, determinato con minimi quadrati in pochi secondi anche per training set grandi. Le reti neurali profonde (NNP) scalano meglio su sistemi complessi: per ogni elemento si usa una rete indipendente che mappa il fingerprint G_i (simmetry functions) in un’energia atomica; i pesi della rete si ottimizzano minimizzando l’errore su energia e forze simultaneamente.

GAP per il carbonio: un caso paradigmatico

Il Gaussian Approximation Potential (GAP) per il carbonio di Deringer, Csányi e colleghi è considerato un riferimento. Addestrato su calcoli DFT di fasi amorfe e liquide del carbonio, il GAP predice correttamente le energie di formazione di grafite, diamante, fullereni e nanotubi con errori di pochi meV/atomo. La dispersione dei fononi del diamante, calcolata con GAP, si sovrappone quasi perfettamente ai dati DFT. Rispetto alla DFT, il GAP è più veloce di diversi ordini di grandezza (tipicamente 10³–10⁵×), pur mantenendo errori entro 2 kJ/mol sui percorsi di annealing simulati, come confrontato con spettroscopia EELS sperimentale.

Accuratezza, velocità e trasferibilità

Il principale limite dei MLIP è la trasferibilità: il potenziale funziona bene solo per configurazioni simili a quelle del training set. Se si usano in una MD a temperatura molto più alta o con composizioni diverse, possono produrre forze fisicamente assurde. Per questo il training set deve essere “diversificato con intenzione”: includere stati lontani dall’equilibrio, superfici, difetti, stati di alta pressione — tutti i contesti fisici rilevanti per l’applicazione prevista.

Domande frequenti

Che differenza c’è tra un potenziale classico e un MLIP?

Il potenziale classico ha una forma funzionale fissa (es. Lennard-Jones, Tersoff) con pochi parametri calibrati vicino all’equilibrio. Un MLIP non ha forma funzionale: usa una rete neurale o un kernel per apprendere direttamente la relazione struttura-energia dal training set DFT. I MLIP sono più accurati, specialmente lontano dall’equilibrio, ma richiedono un training set diversificato e sono più lenti dei potenziali classici (anche se molto più veloci della DFT).

Come si ottiene il training set per un potenziale ML?

Si generano configurazioni diverse del sistema (bulk, superfici, difetti, stati di alta pressione, strutture perturbate casualmente) e si calcolano energia e forze con DFT per ciascuna. La chiave è la diversità: campionare sia stati di equilibrio che stati lontani dall’equilibrio. Si usano tecniche di clustering o force-binning per evitare ridondanze e ridurre il costo computazionale mantenendo la copertura dello spazio.

Il GAP può simulare la transizione grafite-diamante?

Sì: potenziali GAP per il carbonio addestrati su fasi amorfe e liquide sono stati usati per simulare la nucleazione della transizione grafite-diamante, riproducendo il meccanismo microscopico e concordando con i dati sperimentali. Il GAP predice anche correttamente le energie di formazione di nanotubi e fullereni, che non erano nel training set originale.

Quante configurazioni DFT servono per un buon MLIP?

Tipicamente alcune migliaia fino a qualche decina di migliaia di strutture, a seconda della complessità del sistema. La qualità conta più della quantità: 2 000 strutture ben scelte (diverse, rappresentative dei diversi regimi fisici) valgono più di 10 000 strutture simili tra loro. La validazione incrociata è essenziale per verificare che il potenziale non abbia memorizzato il training set.

Perché le forze sono più utili delle energie nel training?

Ogni calcolo DFT produce una sola energia ma 3N forze (tre componenti per ogni atomo). Includere le forze nel training moltiplica di fatto il numero di punti di training disponibili, accelerando la convergenza del potenziale. Alcuni studi mostrano che aggiungere le forze riduce l’errore sull’energia di un fattore 2–5 a parità di strutture DFT calcolate.