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
© Copyright 2025 ExpyDoc