Dinamica di rilassamento per un polimero sottoposto a forze

Universit`
a degli Studi di Padova
` DI SCIENZE MATEMATICHE, FISICHE E NATURALI
FACOLTA
Corso di Laurea in Fisica
Tesi di laurea triennale
Dinamica di rilassamento
per un polimero sottoposto
a forze estensive
Candidato:
Relatore:
Giacomo Bighin
Prof. Enzo Orlandini
Matricola 564280-SF
Anno Accademico 2008–2009
INDICE
1 Introduzione
1
2 Modellizzazione di un polimero
2.1 La catena ideale
2
2.2 La catena “reale”
4
2
3 La Dinamica Molecolare
5
3.1 La scelta del tipo di simulazione
3.2 La Dinamica Molecolare
5
3.3 L’algoritmo Velocity-Verlet
6
3.4 La scelta del potenziale
7
3.5 La forza di Langevin
9
4 Teoria del rilassamento di un polimero
4.1 Il modello di Rouse
10
4.2 Correzione allo scaling
11
5 La simulazione
11
5.1 Il polimero semplice
12
5.2 Il polimero circolare
13
5.3 Il polimero lineare annodato
5.4 Confronto finale
14
6 Conclusioni e prospettive
7 Ringraziamenti
8 Figure
5
10
14
15
16
17
iii
iv
1
INTRODUZIONE
In un seminario dal titolo “There’s plenty of room at the bottom”, già nel
1959 Richard Feynman intuì la possibilità di manipolare direttamente
atomi, molecole e polimeri:
“[· · · ] sarebbe interessante, credo, in linea principio, per un
fisico poter sintetizzare qualunque sostanza descrittagli da un
chimico. Dai gli ordini e il fisico la sintetizza. Come? Mette
giù gli atomi esattamente dove il chimico dice, e così crea la
sostanza. La comprensione dei problemi della chimica e della
biologia potrebbe essere molto maggiore se sviluppassimo le nostra
capacità di vedere quel che stiamo facendo, e di fare le cose su scala
atomica. É uno sviluppo che non credo possa essere evitato.”
Anche a partire da quel discorso, e comunque solo in tempi piuttosto
recenti, si è sviluppata la nanotecnologia, cioè il controllo della materia
a livello atomico e molecolare. Le nanotecnologie hanno aperto un
campo di indagine totalmente nuovo, dando la possibilità di analizzare
sperimentalmente proprietà fisiche in precedenza non indagabili. In
particolare le tecniche più usate per manipolare una singola molecola
sono l’AFM (Atomic Force Microscope) e l’Optical Tweezer. Entrambe
queste tecniche sono state usate massicciamente per lo studio della fisica
della singola macromolecola, in particolare per esplorare la risposta
elastica di quest’ultima ad una forza di elongazione o compressione. In
particolare si è studiata la dinamica di rilassamento di vari biopolimeri
(DNA, RNA, proteine) a partire da configurazioni fortemente allungate
e quindi di non equilibrio (si veda ad esempio [1] e [2]). É questo il
soggetto di cui si occupa questa tesi.
Di pari passo con gli sviluppi della Fisica Sperimentale si sentì
l’esigenza di nuovi strumenti matematici per descrivere le proprietà
della materia su nanoscala. Risultati fondamentali per il lavoro che
segue sono quelli ricavati dal chimico P. Flory (cui venne assegnato
il Nobel per “le sue fondamentali conquiste, sia teoriche che sperimentali,
nel campo della chimica fisica della macromolecole”), e quelli di P.R. Rouse,
sulla dipendenza del tempo di rilassamento di un polimero dal numero
di molecole che lo compongono.
In questa situazione risulta fondamentale un confronto costante tra la
teoria e i dati sperimentali; il modo più diretto per fare questo confronto
è la simulazione mediante computer del modello teorico proposto. Nel
corso del mio lavoro di tesi ho simulato al computer il rilassamento
e l’approccio all’equilibrio di un polimero in soluzione, sotto varie
configurazioni e topologie iniziali, per poi confrontare i risultati con
alcune previsioni teoriche, di cui ho già accennato sopra e che verranno
illustrate nel corso della tesi, e con alcuni dati sperimentali ottenuti dal
rilassamento di una molecola di DNA (si veda ad esempio [6] e [7],
in particolare quest’ultimo osserva dal punto di vista sperimentale lo
stesso fenomeno che andrò a simulare).
Nel corso dei capitoli seguenti verranno discussi i seguenti punti:
• La modellizzazione di un polimero
• La scelta degli algoritmi più adatti per simulare al computer il
rilassamento di un polimero.
• La dinamica del processo di rilassamento.
1
• Infine verrano discussi, alla luce dei risultati teorici ricavati ottenuti, i dati ottenuti della simulazione.
2
MODELLIZZAZIONE DI UN POLIMERO
Con il termine polimero si indica una molecola formata attraverso
la ripetizione di molte unità elementari, chiamate monomeri, tenuti
insieme da legami covalenti a formare una struttura periodica.
Per simulare al computer un polimero dobbiamo essere in grado
di costruire un modello soddisfacente, che tenga conto delle seguenti
caratteristiche fondamentali:
• La connettività, cioè il forte legame covalente che agisce tra due
monomeri vicini; se si cerca di modificarne la distanza interviene
un’intensissima forza di richiamo (in prima approssimazione si
tratta di una forza elastica). L’intensità della forza è tale che la
distanza tra due monomeri vicini si può considerare costante.
• La repulsione tra due monomeri non vicini, e la conseguente naturale tendenza del polimero a non stare eccessivamente ripiegato
su se stesso.
Il modello di Kuhn
Sviluppiamo dunque un modello coarse, detto modello di Kuhn, nel
quale il polimero è costituito da una serie di “perle” tenute insieme
come a formare una collana: una “perla” rappresenta la lunghezza
massima per la quale il polimero è ancora rigido, mentre la giunzione
libera tra due “perle” rappresenta il fatto che sopra tale distanza caratteristica il polimero è flessibile. Ad esempio una molecola di DNA a
doppia elica può essere modellizzata con “perle” del diametro di circa
50 nm, equivalenti a circa 150 basi nucleotidiche.
Questo approccio verrà usato
regolarmente nel corso di questo lavoro: d’ora in poi con il
termine monomero indicheremo sempre non un monomero
dal punto di vista strettamente chimico, ma un monomero
di Kuhn, che identifica una porzione rigida del polimero, in generale formata da più molecole.
Se introduciamo nell’ordine la
connettività e la repulsione possiamo ricavare un buon modello
Figura 1: Il DNA: è un polimero rigido?
di polimero.
2.1
LA CATENA IDEALE
Il modo più semplice per modellizzare un polimero consiste nel rendere
conto solamente della sua caratteristica più evidente: i legami covalenti
tra monomeri vicini. Come già descritto precedentemente, i legami
covalenti sono di solito molto forti e quindi tale caratteristica viene
rispettatta imponendo un vincolo fisso sulla distanza tra monomeri
2
primi vicini. Non vengono fatte altre assunzioni sull’orientazione dei
legami o sulla torsione. Viene trascurata in questo modello qualunque
interazione che non sia quella che si oppone allo spostamento di due
monomeri vicini dalla loro distanza caratteristica b, che viene anche
detta distanza di Kuhn. Ovviamente questo modello è approssimato,
ma ci sono molti sistemi polimerici che vengono ben descritti anche da
un approccio di questo tipo.
Cerchiamo ora di ricavare qualche quantità rilevante dal modello
appena introdotto. Queste quantità ci saranno utili nell’analisi del
processo di rilassamento. Definiamo l’end-to-end vector1 , che non è
altro che il vettore che va dal primo e all’ultimo monomero. Per una
catena di N + 1 monomeri, indicando con ri il vettore tra l’i-esimo e
l’i + 1-esimo vettore, si ha:
Ree =
End-to-end vector
N
∑ ri
i =1
Denotando con h#i la media di ensemble, cioè la media fatta su tutte
le possibili configurazioni dell’oggetto all’equilibrio termico, si vede
che per l’isotropia dello spazio:
hRee i = 0
Chiediamoci ora quanto vale il quadrato dell’end-to-end vector:
*
hR2ee i = hRee · Ree i =
N
∑ ri
!
N
·
∑ rj
!+
j =1
i =1
=
N
N
∑ ∑ h ri · r j i
i =1 j =1
Ricordando l’interpretazione geometrica del prodotto scalare, come
proiezione di un vettore sull’altro, nel caso di un omopolimero2 i cui
legami abbiamo lunghezza b abbiamo che:
N
hR2ee i = b2 ∑
N
∑ hcos
i =1 j =1
θij i
Ancora per l’isotropia dello spazio la media degli angolo θ di legame
dà 0, ad eccezione degli N legami per cui i = j, per i quali si ha
cos (θ ) = 1. Possiamo quindi concludere:
hR2ee i = Nb2
(2.1)
Il raggio giratore è un altra quantità che useremo, analogamente alla
end-to-end distance, per studiare l’approccio all’equilibrio del polimero.
É definito per un omopolimero nel seguente modo:
R2g =
1
N
N
∑ (Ri − Rcm )2
i =1
dove indichiamo con Ri la posizione dell’ i-esimo polimero e con
Rcm = N1 ∑iN=1 Ri la posizione del centro di massa.
Ricaviamo infine un risultato che ci sarà utile in seguito. Abbiamo
già dimostrato che per una catena ideale valgono le seguenti relazioni:
hRi = 0
hR · Ri = Nb2
1 Il modulo dell’end-to-end vector viene chiamato end-to-end distance
2 Un polimero in cui tutti i monomeri, e di conseguenza le distanza dei legami covalenti,
sono identiche
3
Raggio giratore
Il teorema del limite centrale ci assicura che per N grande la distribuzione di R al variare delle
√ configurazioni sarà gaussiana, con
p
dispersione σ = hR · Ri = Nb2 :
1
R·R
P(R) = √
exp −
2Nb2
2πNb2
Energia libera per una
catena ideale
Ricordando che la probabilità di una configurazione è proporzionale al volume occupato nello spazio delle fasi (Ω(R)), si può trovare
l’entropia del polimero ideale:
S(R) = k ln Ω(R) = k ln (αP(R))
∆S(R) = S(R) − S(0)
per poi ricavare l’espressione dell’energia libera di Helmholtz:
∆F (R) = − T∆S(R) = kT
2.2
Calcolo dell’esponente
di Flory
R
Nb2
(2.2)
LA CATENA “REALE”
L’ultimo aspetto da introdurre per modellizzare realisticamente un
polimero è l’interazione di volume escluso, cioè il fatto che alcune zone
dello spazio non sono accessibili al polimero, in quanto già occupate
da altre parti del polimero stesso. É intuitivo credere che, dal momento
in cui si costringe il polimero a non ritornare su se stesso, il volume
occupato a parità di lunghezza deve aumentare.
Per mostrarlo consideriamo un polimero formato da N monomeri,
immerso in una soluzione, con raggio giratore R. I monomeri che
costituiscono il polimero occupano una sfera di raggio R e di volume proporzionale a R3 . Se indichiamo con v il volume di ciascun
monomero, e con ρ = RN3 la densità, otteniamo che la probabilità Pe
che un secondo monomero si trovi nella zona esclusa dal primo è
proporzionale a:
vρ = v
N
R3
Attribuendo ad ogni interazione di volume escluso un’energia kT,
otteniamo che il costo energetico totale delle interazioni di volume
2
escluso sarà Ee = kTNPe = kTv N
. Ora possiamo ricavare l’energia
R3
totale di una catena reale, inserendo l’espressione per l’energia per
R2
libera di un polimero Fl = kT Nb
2 ricavata nella 2.2, dove si è indicata
con b la lunghezza di Kuhn, cioè la lunghezza di un legame covalente
tra due monomeri. Otteniamo:
N2
R2
Ftot = Ee + Fl = kT v 3 +
R
Nb2
Cerchiamo il valore di R all’equilibrio, quello cioè che minimizza
3
tot
l’energia, imponendo ∂F
∂R = 0, e otteniamo :
1
2
3
RF ≈ v 5 b 5 N 5
(2.3)
Successivamente a questa approssimazione iniziale, dovuta a Flory,
sono stati messi a punto modelli più accurati che mostrano che R F ∝ N ν
con ν ≈ 0.588. L’esponente ν viene appunto chiamato esponente di
3 A meno di costanti moltiplicative: il fatto notevole è come R scala al variare N.
4
Flory4 . Si noti che nel caso di una catena ideale si ricava dalla 2.1 che
ν = 12 .
3
LA DINAMICA MOLECOLARE
Dopo aver ricavato un un modello realistico per descrivere un polimero, analizziamo quali sono gli algoritmi più adatti per simularne
l’approccio all’equilibrio.
3.1
LA SCELTA DEL TIPO DI SIMULAZIONE
I tipi di simulazione usati per studiare una molecola sono in buona
sostanza due: la Dinamica Molecolare e gli algoritmi di tipo Montecarlo. L’idea che sta alla base della Dinamica Molecolare è quella di
integrare passo per passo le equazioni della dinamica del sistema, per
visualizzarne direttamente l’evoluzione: la traiettoria è determinata
dalla forze agenti sul sistema. Al contrario in un approccio di tipo
Montecarlo vengono effettuate delle mosse casuali, create ad hoc (si
veda ad esempio [3], pag. 392). A causa di questa differenza:
• Le simulazioni di tipo Montecarlo permettono di simulare agevolmente polimeri di dimensioni maggiori, ma allo stesso tempo
introducono delle approssimazioni ulteriori che rendono più difficile ricondursi, dai risultati della simulazioni, alla reale fisica del
problema.
• Le simulazioni di Dinamica Molecolare ci forniscono un modo per
calcolare anche le quantità dinamiche del sistema, e ci permettono
di studiare fenomeni quali il trasporto di calore, la diffusione, e,
appunto, il rilassamento.
Questi motivi al momento di scegliere l’algoritmo di simulazione
ci hanno fatto propendere per una Dinamica Molecolare. Inoltre dal
momento che i polimeri simulati non hanno dimensioni troppo grandi
(al più 120 monomeri) la maggior lentezza della Dinamica Molecolare
rispetto al Montecarlo non è stata un fattore determinante.
3.2
LA DINAMICA MOLECOLARE
Con il termine Dinamica Molecolare si indica tutta una serie di tecniche
atte a simulare al computer l’evoluzione temporale di atomi e molecole,
integrando numericamente e ripetutamente le equazioni del moto. In
altre parole andando ad integrare le equazioni del moto, supponendo
di conoscere le forze che agiscono sul sistema e nel sistema per ogni
configurazione, possiamo tabulare le posizioni e le velocità in funzione
del tempo. Per far ciò abbiamo bisogno di:
• Un algoritmo per integrare numericamente le equazioni del moto.
4 É interessante notare che ν è anche l’inverso della dimensione frattale dell’oggetto. In
effetti le catene ideali e reali hanno anche grande interesse dal punto di vista puramente
matematico, dove vengono studiate sotto i nomi di random walk e self avoiding random
walk
5
Dinamica Molecolare vs.
Montecarlo
• Un potenziale da cui ricavare la forza, che fornisca risultati in
accordo con le richieste esposte nel capitolo 2. Esso può essere
esatto in alcuni casi, ma nella stragrande maggioranza delle simulazioni molecolari vengono usati potenziali empirici, come ad
esempio quello di Lennard-Jones.
• Un’ulteriore forza, di tipo stocastico, che tenga conto del moto
Browniano della molecola all’interno del solvente, dovuta alla
presenza di una bagno termico alla temperatura T. Aggiungendo
tale forza ottengo per ogni run della simulazione una traiettoria
differente. Per ottenere una buona statistica vado a fare la media
su queste realizzazioni della forza.1
3.3
L’ALGORITMO VELOCITY-VERLET
Le equazioni del moto sono state integrate numericamente, usando l’algoritmo Velocity-Verlet. Il problema dell’integrazione delle equazioni
del moto, consiste nel ricavare, date al tempo t le posizioni e le velocità
delle particelle e le forze agenti sul sistema, la configurazione al tempo
t + τ. Ricaviamo intanto l’algoritmo di Verlet sviluppando x(t ± τ ) in
serie di Taylor fino al terz’ordine:
(
2
3 ...
x(t + τ ) = x (t) + x˙ (t) + τ2 x¨ (t) + τ6 x (t) + O(τ 4 )
2
3 ...
x(t − τ ) = x (t) − x˙ (t) + τ2 x¨ (t) − τ6 x (t) + O(τ 4 )
Sottraendo membro a membro e riorganizzando i termini, si ottiene:
x(t + τ ) = 2x(t) − x(t − τ ) + τ 2 x¨ (t) + O(τ 4 )
Come si vede l’errore sulla posizione è nell’ordine di τ 4 . La stessa
cosa non si può dire per l’errore sulla velocità, che usando l’approssimazione numerica della derivata:
v(t) ≈
x(t + τ ) − x(t − τ )
+ O(τ 2 )
2τ
risulterebbe essere nell’ordine di τ 2 . Per risolvere questo problema
bisogna modificare leggermente l’algoritmo, ottenendo il cosiddetto
Velocity-Verlet. Per far avanzare il sistema di un tempo τ si calcola
innanzitutto la velocità al tempo t + τ2 e la posizione al tempo t + τ,
svilippando ancora in serie di Taylor:
v(t +
τ
τ F(x(t))
) = v(t) +
+ O(τ 2 )
2
2 m
τ
τ 2 F(x(t))
x(t + τ ) = x(t) + τv t +
= x(t) + v(t)τ +
+ O(τ 3 )
2
2
m
per poi ricalcolare le forze e far avanzare la velocità di un altro mezzo
step temporale:
τ τ F(x(t + τ ))
v(t + τ ) = v t +
+
2
2
m
Con l’algoritmo Velocity-Verlet abbiamo un errore sulle velocità che
è nell’ordine di τ 4 , da confrontare con l’errore nell’ordine di τ 2 per
l’algoritmo di Verlet. L’effetto più evidente di questa maggior precisione
1 Il numero minimo di run per ottenere una buona statistica nel mio caso andava da 150 a
200
6
sulle velocità è da una miglior conservazione dell’energia, in assenza di
forze esterne, quando si usa l’algoritmo Velocity-Verlet. Altri vantaggi
dell’algoritmo sono:
• La sua efficienza computazionale. In particolare la funzione
F (x(t))
x¨ = m che valuta l’accelerazione data dalla forza agente su
ogni particella è in genere molto pesante computazionalmente;
nell’algoritmo di Verlet viene invocata una sola volta per iterazione. Per ottenere una precisione simile con un algoritmo di tipo
Runge-Kutta dovremmo arrivare al quarto ordine, nel quale la
funzione F (x(t)) viene invocata quattro volte. A seconda del tipo
di simulazione il vantaggio in termini di tempo di simulazione
può variare; nel caso della mia simulazione la funzione F (x) è
O(n2 ) e impiega gran parte del tempo totale. Quindi il vantaggio
è considerevole.
• La simmetria rispetto all’inversione temporale. A meno di errori di arrotondamento la simulazione può essere invertita per
ritornare alle condizioni iniziali.
3.4
LA SCELTA DEL POTENZIALE
Data l’estrema complessità delle interazioni tra molecole, la scelta di un
potenziale per simulazioni molecolari spesso ricade su un potenziale
empirico, che possa tenere conto, aggiustando vari parametri, della
complessità delle interazioni tra molecole. Modellizzando un polimero,
è necessario scegliere due tipi differenti di potenziale che implentino
le caratteristiche di connettività e di repulsione descritte a pagina 2.
Abbiamo quindi bisogno di:
• Un potenziale “a molla”, che descriva i legamenti covalenti tra
monomeri consecutivi lungo la catena (connettività).
• Un potenziale repulsivo, che tenga conto del fatto che due monomeri differenti non possono sovrapporsi, e, anzi, preferiscono
stare ad un certa distanza (interazione di volume escluso).
A seconda poi del caso si possono introdurre ulteriori potenziali
per tenere conto della resistenza al ripiegamento o alla torsione del
polimero. In questa tesi si è scelto di simulare un polimero non rigido,
che dunque non presenta resistenza né al ripiegamento né alla torsione.
La scelta del potenziale per la parte repulsiva è ricaduta su un
potenziale di tipo Lennard-Jones (LJ) troncato; indicando con r la
distanza tra due monomeri, l’espressione del potenziale è la seguente2 :
(
6 12
4e σr
− σr
+ e se r < re
VLJ (r ) =
0
se r ≥ re
nel potenziale di Lennard-Jones non troncato e è la rappresenta la
profondità della buca di potenziale, mentre σ e l’ascissa del primo
zero. Nel nostro caso σ e e non hanno più questo significato fisico, dal
momento che il potenziale è stato troncato e traslato lungo l’asse delle
ordinate perché tenda asintoticamente a zero. Con re indichiamo il
1
minimo del potenziale LJ standard, che vale σ · 2 6 .
Per la parte di legame covalente abbiamo scelto un potenziale di tipo
FENE (Finitely Extensibile Nonlinear Elastic potential, si veda anche
2 La e aggiunta rispetto al potenziale di LJ standard serve far in modo che il minimo del
potenziale valga 0.
7
Il potenziale
Lennard-Jones troncato
Il potenziale FENE
Potenziale tra primi vicini
300
Potenziale totale
FENE
Lennard-Jones
250
Potenziale
200
150
100
50
0
0.7
0.8
0.9
1
1.1
1.2
Distanza intermolecolare
1.3
1.4
1.5
Figura 2: Il potenziale tra due primi vicini, scomposto anche nella parte
repulsiva e nella parte attrattiva.
[11]), spesso usato, invece di un più semplice potenziale armonico, per
simulazioni molecolari. Indicando con r la distanza tra due monomeri
consecutivi:
1 2
VFENE (r ) = − krm
log 1 −
2
r
rm
2 !
Qui k è una costante elastica, mentre rm indica la massima distanza
raggiungibile prima della rottura del legame tra due primi vicini. Sommando questi due potenziali, possiamo modellare la forza da applicare
tra due monomeri primi vicini, che consisterà sia di una parte repulsiva
che di una parte attrattiva:
FPV (r ) = −∇ VFENE (r ) + VLJ (r )
L’energia totale di un polimero aperto formato da N monomeri sarà
data quindi dalla somma dell’energia cinetica e di entrambi i potenziali
che abbiamo descritto:
H=
N −1
p2i
1 N N
+
V
(
r
)
+
∑ 2m ∑ FENE i,i+1 2 ∑ ∑ VLJ (ri,j )
i =1
i =1
i =1 j =1
N
j 6 =i
In figura 2 mostriamo i potenziali FENE e Lennard-Jones, e la loro
somma. I parametri usati per disegnare il grafico sono gli stessi della
simulazione3 e cioè:
Lennard-Jones
FENE
e
σ
κ
rm
1.0
1.0
30.0
1.5
3 Solo nel caso della simulazione con polimeri lineari annodati il parametri e è stato
aumentato a 3.35
8
3.5
LA FORZA DI LANGEVIN
Il polimero simulato si muove in un solvente. Vogliamo ricavare per un
monomero considerato singolarmente l’effetto del moto Browniano, dovuto agli urti del monomero con le particelle di solvente a temperatura
T costante in cui è immerso (bagno terminico). Proviamo a ricavare
in maniera euristica un’espressione della forza, denominata forza di
Langevin, alla quale sono soggetti i monomeri per effetto del solvente.
Quando un monomero si sta muovendo in una direzione, incontrerà
in media più molecole di solvente lungo quella direzione, le quali ostacoleranno il suo moto. Questa fatto può essere espresso nel seguente
modo:
m
dv(t)
= −γv(t)
dt
dove γ è un coefficiente legato all’attrito4 . Questa forza però da
sola non basta a render conto del moto browniano del polimero che si
osserva nella realtà, dato che predice che il moto si dovrebbe fermare
in un tempo più o meno breve. Dobbiamo aggiungere un secondo
termine che denoteremo con fL (t) che tenga conto delle collisioni con
le molecole di solvente che avvengono su una scala temporale più breve
rispetto all’effetto generale dell’attrito. Otteniamo:
m
dv(t)
= −γv(t) + fL (t)
dt
(3.1)
Chiediamoci quali condizioni imporre su fL (t):
• Dato che la forza netta risultante dovuta agli urti è già inclusa nel
termine dipendente dalla velocità possiamo richiedere che fL (t)
abbia media nulla:
hfL (t)i = 0
(3.2)
• Data la natura stocastica della forza, i valori a tempi diversi
devono essere scorrelati:
hfL (t) · fL (t + τ )i = 0
τ 6= 0
il che è equivalente ad assumere:
hfL (t1 ) · fL (t2 )i = σ2 δ(t2 − t1 )
(3.3)
La scelta di σ2 è legata a soddisfare le richieste che il sistema, per
t 1 raggiunga l’equilibrio alla temperatura T del bagno termico in
cui è immerso.
L’equazione 3.1 può essere risolta analiticamente, fino ad ottenere la
seguente soluzione:
Z
γt
1 t
γ
v(t) = v(0) exp −
+
exp − (t − τ )
fL (τ )dτ
m
m 0
m
Sfruttando la 3.2, cioè il fatto che la funzione fL ha media nulla,
otteniamo:
γt
hv(t)i = v(0) exp −
m
9
La relazione di
fluttuazionedissipazione
Allo stesso modo dall’integrazione delle equazioni 3.1 e 3.3 possiamo
ottenere:
hv (t)i =
2
v20 exp
γt
−2
m
+
σ2 1 − e−2γt
2γm
Facendo tendere t all’infinito si arriva all’equilibrio, e si ottiene infine:
hv2 (∞)i =
σ2
2γm
All’equilibrio varrà il teorema di equipartizione dell’energia:
3
1
mhv2 (∞)i = k B T
2
2
da cui si ottiene la relazione di fluttuazione-dissipazione che correla
appunto σ2 e γ:
σ2 = 6γk B T
4
T E O R I A D E L R I L A S S A M E N TO D I
UN POLIMERO
Lo scopo di questo lavoro è studiare l’approccio all’equilibrio di un
polimero, applicando una forza esterna e poi rilasciandola, al variare
delle condizioni iniziali: siamo quindi interessati a due aspetti:
• Il tempo che il polimero impiega a raggiungere l’equilibrio, in
funzione del numero di monomeri, della topologia iniziale (lineare, ad anello, annodata) e dell’intensità della forza applicata
inizialmente. Si parlerà di tempo di rilassamento τR .
• Il “modo” in cui il polimero raggiunge l’equilibrio, e quindi la
dipendenza di alcune quantità significative dal tempo.
Abbiamo già introdotto la end-to-end distance e il raggio giratore.
Queste quantità sono la scelta naturale per descrivere il rilassamento di
un polimero, dato che rappresentano una misura diretta dell’elongamento del polimero per effetto della forza esterna. Varieranno quindi
da un massimo, all’ultimo istante in cui viene applicata la forza, ad un
minimo stabile, entro una certa banda di rumore, quando il sistema ha
raggiunto l’equilibrio. É questa proprio una della definizioni classiche
di equlibrio: la situazione per cui gli osservabili macroscopici associati
al sistema non variano. Ricaviamo dunque alcuni risultati teorici sul
tempo di rilassamento, che andranno poi confrontati con i dati ricavati
dalla simulazione.
4.1
IL MODELLO DI ROUSE
Il modello di Rouse descrivere esattamente la dinamica di un polimero
ideale (in assenza quindi di interazioni di volume escluso), in mancanza
4 É da notare la somiglianza di quest’espressione con la legge di Stokes per il moto in un
fluido, nella quale si ha γ = 6πηr
10
di modi idrodinamici. Tra i vari risultati tale modello prevede la seguente legge di scala tra il tempo medio di rilassamento τR , l’estensione
media e il grado di polimerizzazione N:
τR ∝ NR2
(4.1)
Se s’introduce l’interazione di volume escluso le equazioni della
dinamica non sono più risolvibili esattamente. Tuttavia possiamo
inserire nella 4.1 la relazione R ∝ N ν ricavata nella 2.3, ottenendo:
τR ∝ N 2ν+1
(4.2)
Ribadiamo che uest’ultimo passaggio non è del tutto rigoroso, dal
momento che l’equazione 2.3 vale per un polimero all’equilibrio, mentre, al contrario, l’equazione 4.2 vale per un polimero non all’equilibrio.
Ciononostante, come vedremo nella sezione successiva, questo argomento di scala riesce comunque a prevedere con buona accuratezza
l’approccio all’equilibrio di un polimero reale.
Notiamo che l’argomento appena esposto vale tanto per la end-toend distance quanto per il raggio giratore, dal momento che tra le due
quantità c’è una relazione lineare (si veda ad esempio [3], pag. 63).
4.2
CORREZIONE ALLO SCALING
L’ultima osservazione da fare prima di completare il nostro modello di
polimero, è osservare che i calcoli presentati nei capitoli sono validi per
un polimero con un numero N di monomeri arbitrariamente grande.
Dal momento che il numero di monomeri è finito dobbiamo introdurre
le cosiddette correzioni allo scaling. Se ad esempio nel caso di un
polimero arbitrariamente lungo vale la relazione:
h Ri = AN ν
nel caso di un polimero finito dobbiamo applicare la seguente correzione
B
C
D
(4.3)
h Ri = AN ν 1 + 1 + + 3 + · · ·
N
N2
N2
1
In realtà già fermandosi al termine di ordine N − 2 si ottiene in
parecchi casi un ottimo accordo con i dati sperimentali.
5
LA SIMULAZIONE
Siamo ora in possesso di tutti gli strumenti teorici per analizzare i
risultati della simulazione al computer dei polimeri. In particolare
analizziamo l’approccio all’equilibrio per tre topologie differenti:
• Polimero lineare semplice.
• Polimero lineare con un nodo.
• Polimero circolare.
Nelle simulazioni è stata applicata una forza esterna ai due capi del
polimero, nel caso di un polimero lineare; nel caso di un polimero
circolare la forza è stata applicata a due punti diametralmente opposti.
11
I parametri della
simulazione
L’intensità della forza è stata accuratamente scelta in modo tale che le
catene polimeriche fossero quasi al limite di rottura dei legami covalenti.
La simulazione inizia con la forza esterna applicata, ad un tempo
ti < 0, scelto in modo che il polimero raggiunga entro t = 0 una
situazione di equilibrio con la forza esterna. A t = 0 la forza esterna
viene rimossa e il polimero è libero di raggiungere un nuovo equilibrio.
Riportiamo qui i parametri usati per la simulazioni:
Parametro
T
(temperatura del bagno termico)
γ
(viscosità di Langevin)
dt
(timestep di integrazione)
m
(massa di un monomero)
F
(modulo della forza esterna)
Valore
1
8
0.01
1
100
Tutti i valori riportati sono in unità arbitrarie. D’ora in poi quando ci
riferiremo ai tempi sarà sempre in unità di integrazione.
5.1
IL POLIMERO SEMPLICE
Per quanto riguarda il polimero lineare ci aspettiamo che, quando viene
raggiunto l’equilibrio con la forza esterna, la sua end-to-end distance sia
all’incirca Nb, dove con b indichiamo la distanza di Kuhn. Se prestiamo
fede all’equazione 2.3, quando il polimero ha raggiunto l’equilibrio la
end-to-end distance dovrebbe, invece, scalare come N ν dove ν ≈ 0.588.
Inoltre secondo le stime di Rouse i tempi di rilassamento dovrebbero
scalare come N 2ν+1 .
Per verificare tutte queste ipotesi possiamo mettere in grafico le endto-end distance (o equivalentemente i raggi giratori) in funzione del
tempo, per polimeri di varia lunghezza, ed effettuare dei riscalamenti:
• Se riscaliamo le distanze dividendo per un fattore N i grafici
relativi ai diversi polimeri, dovrebbero sovrapporsi quando i
polimeri sono in equilibrio con la forza estena, e dunque per t ≤ 0
(figura 5).
• Se riscaliamo invece le distanze dividendo per un fattore N ν i
grafici dovrebbero tendere all’equilibrio in tempi differenti, ma la
end-to-end distance di equilibrio, ora riscalata, dovrebbe essere la
stessa (figura 6, tutte le figure relative al presente capitolo sono
riportate nel capitolo 8).
• Se in aggiunta al punto precedente riscaliamo i tempi dividendo
per un fattore N 2ν+1 dovrebbero coincidere anche i tempi di rilassamento. I polimeri possono comunque approcciare l’equlibrio in
maniera differente, ma il tempo entro il quale l’equilibrio viene
raggiunto, perché le stime di Rouse siano verificate, dovrà essere
lo stesso (figura 7).
Come si può osservare nelle figure 4, 5, 6 e 7 queste tre previsioni sono
verificate con buona accuratezza e il polimero effettivamente rispetta
l’argomento di Flory e le stime di Rouse. Come previsto, però abbiamo
dovuto applicare le correzioni allo scaling. I polimeri raggiungono lo
12
stesso equilibrio se, al momento di effettuare i riscalamenti, definiamo
un valore efficace della costante di Flory, νefficace = 0.675. É più corretto
però esprimere il risultato ottenuto nei termini mostrati nella sezione
4.2; mettiamo a sistema i dati trovati con l’equazione 4.3:

νefficace = A0 · 20ν 1 + √B + · · ·

A
·
20


20

..
.



 A · 80νefficace = A0 · 80ν 1 + √B + · · ·
80
possiamo così ricavare il primo coefficiente di correzione allo scaling:
B = −0.782
5.2
IL POLIMERO CIRCOLARE
Consideriamo ora un polimero circolare. Se viene applicata una forza
estensiva a due monomeri diametralmente opposti, i rimanenti monomeri si disporranno come a formare due catene, ciascuna lunga N2
monomeri. A questo punto, se consideriamo le due catene come sostanzialmente distinte, o comunque non interagenti, ci possiamo aspettare
che:
ν
• Le end-to-end distance di equilibrio scalino come N2 .
• I tempi di rilassamento scalino come
2ν+1
N
2
.
In altre parole, trattando un polimero circolare come due catene
distinte formate da N2 monomeri, ci aspettiamo che possano valere
gli stessi argomenti di scaling usati per l’analisi del rilassamento dei
polimeri lineari1 .
Mettiamo quindi in grafico in figura 8 le stesse quantità della figura
6, riscalando come appena illustrato, e osserviamo che tale previsione è
ben rispettata. L’unico caso in cui i dati della simulazione si discostano
sensibilmente da questa previsione è nel caso del polimero circolare
più piccolo, quello formato da 10 monomeri. Evidentemente in questo
caso la forza repulsiva tra monomeri non primi vicini agisce in maniera
considerevole anche tra monomeri appartenenti alle due catene, che
negli altri casi avevamo legittimamente considerato come distinte.
É anche intuibile perché le nostre previsioni sovrastimino la end-toend distance all’equilibrio finale: le forze repulsive impediscono che
il polimero, anche quando in equilibrio con la forza esterna, assuma
una configurazione troppo schiacciata. Di conseguenza i due punti di
applicazione della forza rimagono più vicini e la end-to-end distance è
minore. Rimandiamo alla sezione 5.4 per un confronto con il tempo di
rilassamento per un polimero lineare.
1 Nel caso del polimero circolare, non esistendo una prima e un’ultima molecola, la endto-end distance è definita come la distanza tra i due monomeri diametralmente opposti
scelti per l’applicazione della forza.
13
5.3
IL POLIMERO LINEARE ANNODATO
Analizziamo ora il caso del polimero annodato. Si tratta di un nodo trifoglio (3, 1) tagliato in punto. Ai due capi risultanti viene
applicata la forza estensiva per
t < 0, esattamente come negli
altri casi. In questo caso, però,
possiamo supporre che la end-toend distance quando la forza viene applicata sia proporzionale a
N − 8, dato che circa 8 monomeri
vengono “persi” nel nodo. Mettendo in grafico la end-to-end distance in funzione del tempo con Figura 3: Un corda annodata, topoi consueti riscalamenti, in figura
logicamente equivalente al
9, vediamo che in prima approspolimero lineare annodato
simulato.
simazione l’equilibrio viene raggiunto contemporaneamente; tuttavia guardando più in dettaglio l’accordo non è così buono come nei
due casi precedenti. Questo comportamento è da attribuire ad una serie
di fenomeni:
• Soprattutto per valori di N piccoli, a volte durante l’approccio
all’equilibrio il nodo si può sciogliere. Aumenta dunque la endto-end distance finale, fino ad arrivare ad un valore che tende ad
avvicinarsi alla end-to-end distance di un polimero identico ma
non annodato, si veda la figura 10. Questa tendenza ovviamente
è più o meno marcata a seconda della frazione di polimeri che si
scioglie2 .
• Il nodo costringe la molecola a ritornare su se stessa, dando modo
anche alle forze repulsive che agiscono tra monomeri non primi
vicini di agire. Questo aspetto non è completamente considerato
nel nostro modello, che dunque fornisce risultati sensibilmente
divergenti rispetto all’osservazione.
5.4
Polimero lineare
semplice vs. polimero
circolare
CONFRONTO FINALE
Come già accennato nelle sezioni precedenti, facciamo ora un confronto
tra le varie configurazione studiate, proponendoci di osservare le eventuali differenze nell’approccio all’equilibrio tra un polimero circolare,
uno annodato e uno lineare.
Mettiamo a confronto l’approccio all’equilibrio di una molecolare
circolare composta da N = 38 monomeri e di una molecola lineare
composta da N 0 = N2 = 19 monomeri. Come già mostrato nella sezione
5.2 il tempo di rilassamento scala in maniera identica, ma ci chiediamo
se c’è qualche differenza nel “modo” in cui l’equilibrio viene raggiunto.
Calcoliamo dai dati della simulazione l’approssimazione numerica
della derivata prima, nel seguente modo:
f 0 ( xi ) =
f ( x i +1 ) − f ( x x i −1 )
x i +1 − x i −1
2 Ovviamente nel caso reale il nodo può o sciogliersi oppure non sciogliersi. Osserviamo
una tendenza a raggiungere la end-to-end distance del nodo non sciolto perché i dati in
grafico sono le medie su 200 simulazioni
14
Le derivate della end-to-end distance per il polimero circolare e per
il polimero lineare praticamente coincidono, e quindi se ci limitiamo
a considerare questo parametro un polimero lineare di lunghezza N2
è effettivamente equivalente ad un polimero circolare di lunghezza
N (figura 11). Considerando invece la derivata del raggio giratore si
osserva che questo parametro ci permette di distinguere la topologia
del polimero. Il raggio giratore nel caso del polimero lineare varia in
maniera molto più veloce, rispetto al polimero circolare (figura 12).
Ripetiamo ora la stessa analisi, applicandola però ad un polimero
lineare semplice e ad un polimero lineare annodato. Si è considerato un
polimero lineare annodato, di lunghezza N = 38; per poter studiare i
primi istanti dell’approccio all’equilibrio dobbiamo confrontarlo con un
polimero lineare di eguale lunghezza all’equilibrio con la forza esterna.
Per questo motivo abbiamo scelto un polimero lineare con grado di
polimerizzazione3 N 0 = 38 − 9 = 29.
Vogliamo anche in questo caso capire se dal “modo” in cui si arriva
all’equilibrio è possibile distinguere la topologia. I risultati, in figura
13 e 14, sono sostanzialmente analoghi a quanto visto per il polimero
circolare. Se consideriamo solo la end-to-end distance questa varia
in maniera sostanzialmente identica, indipendentemente dal fatto che
il polimero sia annodato o meno. Se consideriamo, invece, il raggio
giratore, questo ci fornisce qualche informazione sulla topologia: il
polimero lineare semplice ha un approccio all’equilibrio inizialmente
più veloce rispetto al polimero lineare annodato. É ragionevole pensare
quindi che il raggio giratore contenga più informazioni sulla struttura
fisica della molecola rispetto alla end-to-end distance.
Per indagare in maniera più approfondita le informazioni contenute
nel raggio giratore, confrontiamo un polimero circolare, un polimero
lineare annodato e un polimero lineare semplice, tutti di eguale lunghezza N = 38. Ovviamente quando viene raggiunto l’equilibrio con
la forza esterna la topologia del polimero influisce marcatamente sul
raggio giratore; il polimero circolare avrà il raggio giratore minore,
seguito dal polimero lineare annodato e dal polimero lineare semplice.
Al contrario, mano a mano che i tre polimeri approcciano l’equilibrio
finale le informazioni sulla topologia contenute nel raggio giratore
vengono progressivamente perse. All’equilibrio finale il raggio giratore
è funzione solo del grado di polimerizzazione N e non dipende dalla
topologia, come si può vedere in figura 15. É importante notare che la
curva relativa al polimero lineare annodato non raggiunge in maniera
monotona l’equilibrio, dato che il nodo per un tempo abbastanza lungo
si scioglie.
6
CONCLUSIONI E PROSPETTIVE
Il lavoro di questa tesi è consistito nel simulare al computer polimeri
con differenti gradi di polimerizzazione e differenti topologie, nel ricavare un semplice modello teorico per descrivere il comportamento dei
polimeri e nella successiva analisi dei risultati ricavati. Tale analisi ha
confermato che il modello della catena “reale” descrive accuratamente
il comportamento dei polimeri nella maggior parte dei casi che abbiamo
preso in considerazione.
3 In altre parole nel nodo vengono “persi” 9 monomeri
15
Polimero lineare
semplice vs. polimero
lineare annodato
Raggio giratore per le
tre topologie
Alcune prospettive per uno studio successivo, che non sono state
indagate per limiti di tempo, potrebbero essere:
• Continuare il lavoro di analisi prendendo in considerazione anche topologie più complicate, e osservare se compaiono nuove
dipendenze dalla topologia.
• Lo studio del processo inverso, cioè l’approccio all’equilibrio con
la forza esterna, partendo da un polimero all’equilibrio in un
solvente. Ho avuto modo di verificare che il processo non è
del tutto simmetrico rispetto a quello da me studiato; sarebbe
interessante verificare se valgono comunque le stesse leggi di
scala che si sono dimostrate efficaci nello studio del tempo di
rilassamento.
7
RINGRAZIAMENTI
Desidero ringraziare il prof. Enzo Orlandini per la disponibilità e la
pazienza con cui mi ha seguito durante il lavoro di tesi.
16
8
FIGURE
Figura 4: Approccio all’equilibrio per quattro polimeri lineari. Il grado di
polimerizzazione è diverso, e quindi lo sono anche le distanze di
equilibrio e i tempi di rilassamento
Figura 5: All’equilibrio con la forza esterna le end-to-end distance sono
proporzionali a N.
17
Figura 6: Approccio all’equilibrio per polimeri lineari. Riscalando le end-to-end
distance di un fattore N ν i polimeri raggiungono, in tempi diversi, la
stessa distanza di equilibrio.
Figura 7: Approccio all’equilibrio per polimeri lineari, riscalamento sui tempi
secondo Rouse (N 2ν+1 ), sulle distanze secondo Flory (N ν ): i varie
polimeri lineari raggiungono lo stesso equilibrio, nello stesso tempo,
in modi differenti
18
Figura 8: Approccio all’equilibrio per polimeri circolari. Trattando il polimero
circolare come due catene indipendenti con N2 polimeri si ottengono
buoni risultati, tranne che nel caso N = 10.
Figura 9: Approccio all’equilibrio per polimeri annodati.
19
Figura 10: Approccio all’equilibrio per polimeri annodati. Per N = 33 alcuni
polimeri si sciolgono.
Figura 11: Considerando la end-to-end distance (e la sua derivata) il rilassamento di un polimero circolare formato da N monomeri è indistinguibile
da quello di un polimero lineare semplice formato da N2 monomeri.
20
Figura 12: Se consideriamo il raggio giratore, possiamo invece distinguere tra
il polimero circolare e il polimero lineare semplice
Approccio all’equilibrio per un polimero lineare e per uno chiuso (derivata prima)
0.005
Lineare annodato, N=38
Lineare semplice, N=29
Derivata prima della end-to-end distance
0
-0.005
-0.01
-0.015
-0.02
-0.025
-0.03
-0.035
-0.04
0
200
400
600
800
1000
Tempo
Figura 13: Derivata prima della end-to-end distance, per il rilassamento di un
polimero lineare semplice e di un polimero lineare annodato.
21
Approccio all’equilibrio per un polimero lineare e per uno chiuso (derivata prima)
0.001
Lineare annodato, N=38
Lineare semplice, N=29
Derivata prima del raggio giratore
0
-0.001
-0.002
-0.003
-0.004
-0.005
-0.006
-0.007
-0.008
-0.009
0
200
400
600
800
1000
Tempo
Figura 14: Derivata prima del raggio giratore, per il rilassamento di un
polimero lineare semplice e di un polimero lineare annodato.
Figura 15: All’equilibrio con la forza esterna il raggio giratore dipende fortemente dalla topologia; man mano che si raggiunge l’equilibrio finale
tale dipendenza gradualmente scompare
22
BIBLIOGRAFIA
[1] Ueda M. et al. (1999), Direct measurement of DNA by means of AFM,
Nucleic Acids Symposium Series No. 42 245—246, Oxford University
Press
[2] Wang M. D. et al. (2007), Stretching DNA with Optical Tweezers,
Biophysical Journal 72:1335-1346
[3] Rubinstein, M. e Colby R. H. (2003), Polymer Physics, Oxford
University Press, New York
[4] Peliti, L. (2003), Appunti di Meccanica Statistica, Bollati Boringhieri,
Torino
[5] Thijssen, J.M. (1999), Computational Physics, Cambridge University
Press, Cambridge
[6] Ferre, S. e Blanch H. W. (2004), The Hydrodynamics of DNA Electrophoretic Stretch and Relaxation in a Polymer Solution, Biophysical Journal,
87:468-465
[7] Perkins, T. T. et al (1994), Science 264, 819:822
[8] Obermayer B. et al. (2007), Stretching dynamics of semiflexible polymers,
Eur. Phys. J. E 23, 375-388
[9] Humphrey W., Dalke A. and Schulten K. (1996), VMD - Visual
Molecular Dynamics, J. Molec. Graphics, vol. 14, pp. 33-38.
[10] Huang K. (1997), Meccanica Statistica, Zanichelli, Bologna
[11] Kremer K. e Grest, G.S. (1990), J. Chem. Phys. 92, 5057
23

Download Report