Chimica fisica superiore - Compendio (Parte 3 di 3)

8. Simulazione di Sistemi in Fase Condensata

Le tecniche computazionali classiche impiegate per lo studio di sistemi in fase condensata sono: ↪ atomi come particelle classiche

1. il metodo Monte Carlo → studio delle proprietà strutturali e spettroscopiche
2. le simulazioni di Dinamica Molecolare Classica → studio di proprietà strutturali, spettroscopiche e dinamiche

Metodo Monte Carlo

Nel metodo Monte Carlo (MC) si parte da una geometria iniziale della molecola e si genera, ad ogni step, una perturbazione random che sposta un solo atomo.

Nella procedura di Monte Carlo “Metropolis”, che impiega la versione “Markov Chain”, la nuova configurazione generata ad ogni step viene accettata solo se

• ha un’energia inferiore rispetto alla configurazione precedente
• se obbedisce alla distribuzione di Boltzmann

Il metodo Monte Carlo è non deterministico => due simulazioni MC che partono dalla stessa configurazione generano configurazioni differenti

Il metodo MC possiede capacità di “tunneling” = superare barriere energetiche.

La simulazione può essere eseguita impiegando sia coordinate cartesiane che coordinate interne.

Lo svantaggio del metodo MC è che non è adatto allo studio di fenomeni o proprietà dipendenti dal tempo.

DINAMICA MOLECOLARE CLASSICA

Nelle simulazioni di Dinamica Molecolare Classica, si parte da una configurazione iniziale della molecola e si genera, ad ogni time-step (intervallo temporale) una nuova configurazione in cui tutti gli atomi sono stati spostati.

Le nuove configurazioni vengono ottenute risolvendo numericamente un set di equazioni del moto di Newton, cioè un set di equazioni differenziali del 2° ordine accoppiate:

F_i = m_i d²r_i/dt²

Le forze agenti sui singoli atomi sono anche pari a

F_i = - ∂V(r₁, r₂, ...r_N) / ∂r_i i = 1,...,N

dove V(r₁,...,r_N) è detto "Force Field" = campo di forza, e rappresenta il potenziale d'interazione intra- e intermolecolare.

Inizialmente si inseriscono nel programma i dati di input: coordinate e velocità atomiche + force field. Con quest'ultimo si calcolano le forze agenti sugli atomi, e dalle F_i si ottengono a loro volta le nuove coordinate e velocità atomiche. Con tali dati si può reiterare il calcolo, di modo da ottenere la "traiettoria", cioè l'insieme delle configurazioni ottenute in sequenza.

Il time step dev'essere scelto in modo che

• sia sufficientemente grande di modo da effettuare calcoli rapidi
• sia sufficientemente piccolo per permettere di studiare l'evento più veloce

=> time step dell'ordine dei femtosecondi (st ~10^-15s)

Le simulazioni di Dinamica Molecolare Classica sono deterministiche

Formalismo di Lagrange e di Hamilton

Le equazioni del moto possono essere espresse con il formalismo di Lagrange o con quello di Hamilton → portano agli stessi risultati ottenuti risolvendo numericamente le equazioni di Newton.

Vantaggi:

1. si può utilizzare qualsiasi insieme di variabili non ridondanti
2. si possono introdurre variabili addizionali → utili per imporre dei vincoli e per simulazioni in insiemi NVT o NPT

Differenza fra le formulazioni:

• Lagrange → 1 sola eq. differenziale del 2° ordine
• Hamilton → 2 eq. differenziali accoppiate del 1° ordine

Scelta del time step

Il time step è scelto in modo tale da essere di un ordine di grandezza inferiore rispetto al processo più veloce.

Motivi molecolari (vibrazioni e rotazioni) → ~10¹³-10¹⁴ s^-1 → time step ~ 10^-15 s (femtosecondi)

Per studiare processi chimici lunghi (~ millisecondi e secondi) si devono utilizzare degli algoritmi.

Vincoli

Per alleggerire i calcoli, è possibile imporre vincoli su

1. lunghezze di legame atomo-H → simulazioni più lunghe
2. angoli → poco usati perché possono avere influenza sulle proprietà del sistema
3. geometrie molecolari → es. modelli di acqua rigidi

In questo approccio, occorre trattare le osservabili sperimentali come risultanti da proprieta' del sistema, soggette a fluttuazioni casuali che sono dovute alle interazioni del sistema stesso.

Consideriamo una variabile "A":

Valore medio di A: <A> = ∫ A P(A) dA

dove P(A) = "funzione di distribuzione" → probabilita' di misurare un valore compreso fra A e A+dA

Valore quadratico medio di A: <A²> = ∫ A² P(A) dA

Varianza (= indeterminazione sulla misura di A) : σ_A = <A²> - <A>²

Funzione di correlazione di due osservabili "A" e "B" :

C_AB = <AB> - <A><B>

La funzione di correlazione definisce la correlazione fra A e B, cioe’ se variano indipendentemente oppure se e quale dipendenza sussiste.

Coefficiente di correlazione :

ρ = C_AB / σ_Aσ_B     (-1 ≤ ρ ≤ +1)    disuguaglianza di Schwartz

• -se ρ = +1 ⇒ A e B sono completamente correlate (variano insieme)
• -se ρ = 0 ⇒ A e B sono completamente scorrelate (sono indipend.)
• -se ρ = -1 ⇒ A e B sono anticorrelate (sono uguali e opposte)

Funzione di “autocorrelazione temporale” = ρ indica la correlazione di una variabile con se stessa in tempi diversi :

C_AA (t,t’) = <A(t) A(t’)>

Essa e' utile per analizzare la presenza di periodicita’, tendenze o cicli.

le interazioni intermolecolari

3) non sono descritti correttamente gli effetti di polarizzazione

4) non si può studiare la reattività chimica

Per superare tali problemi, si può ricorrere alla Dinamica Molecolare Ab Initio; tale metodo computazionale non impiega campi di forza o parametri empirici per descrivere le interazioni fra atomi, ma calcola le forze agenti sugli atomi al volo ("on the fly") a partire dalla struttura elettronica, che viene ottenuta sfruttando metodi quantomeccanici (generalmente DFT).

La AIMD si basa sull'Approssimazione di Born-Oppenheimer, per cui, in virtù della notevole differenza fra le masse dell'e⁻ e dei nuclei, il moto dei nuclei è indipendente dal moto degli elettroni, i moti degli e⁻ dipendono solo in modo parametrico dalle coordinate nucleari:

V(R) = <Φ₀ |Ĥ_elet ({r_i; R})|Φ₀> + E_nuc (R)

In tal modo

1. il moto dei nuclei è trattato classicamente => si risolvono le equazioni del moto di Newton
2. il moto degli e⁻ è trattato quantisticamente => si risolve l'equazione di Schrödinger elettronica, non relativistica, indipendente dal tempo e a molti corpi

Dunque in una simulazione di AIMD si calcola la struttura elettronica, si ottiene l'energia potenziale, si calcolano le forze agenti sugli atomi come F=-∇V(R) e infine si risolvono le eq. di Newton per ottenere la traiettoria.

3) per sviluppare metodi “linear-scaling”

Il difetto è che non si ottiene una definizione unica degli orbitali localizzati.

Generalmente, si utilizza un algoritmo di localizzazione che trasforma le funzioni di Bloch in "Funzioni di Wannier massimamente localizzate" (MLWF), che hanno la proprietà di decadere esponenzialmente con la distanza dal proprio centro.

- Calcolo di proprietà vibrazionali

Per calcolare le proprietà vibrazionali di cristalli (quali modi normali e frequenze di vibrazione) generalmente si effettuano calcoli di dinamica reticolare nell’ambito dell’approssimazione armonica, svolte su celle unitarie ottimizzate alla temperatura di 0 K.

La dinamica molecolare ab-initio permette di superare questi limiti, in quanto :

1. si possono studiare sistemi in fase condensata (non solo cristalli)
2. si supera l'approssimazione armonica
3. si introducono gli effetti della temperatura

Inoltre, a partire dalla traiettoria calcolata in una simulazione di AIMD, si ottengono

1. il vettore di polarizzazione ⇒ correlato allo spettro IR
2. il tensore di polarizzabilità ⇒ correlato allo spettro Raman