PLOCHA POTENCIÁLNÍ ENERGIE
Zero point energy - Energie nulového bodu • Molekula o určitou část své energie nikdy nemůže přijít • Tzv. Zbytková energie (ZPE) – „vnitřní“ energie molekuly,
která je přítomna vždy (i za 0 K) • Pro každou vibraci
Potenciální energie molekuly • Epot - jedna ze složek celkové energie molekuly • Závisí na geometrii dané molekuly • Charakterizuje stabilitu určité konformace (geometrického
uspořádání) molekuly → čím je potenciální energie menší, tím je molekula stabilnější • Umožňuje určit nejstabilnější (nejpravděpodobnější)
konformaci molekuly
Plocha potenciální energie • Potential energy surface (PES) • Energie molekuly jako funkce její geometrie • Energie vertikální osa, geometrické vlastnosti jako
horizontální osy
Významné body - Stacionární body Body PES, jejichž gradient je nulový vektor. • Minima a maxima (lokální a globální) • Sedlové body
∂E =0 ∂X i
Stacionární body na PES • Během optimalizace se prohledává PES a hledají se
lokální/globální minima. • Přes tranzitní stavy (bariéry) nás optimalizační algoritmy nedostanou. Výsledkem je, že nalezená minima jsou většinou lokální.
M = minimum, TS =1st order saddle point, S = 2nd order saddle point
Konkávní funkce
• F’’(x) > 0 • F’’(x) < 0
Konvexní funkce
• F’’(x) = 0
Inflexní bod
Hessova matice (Hessián, H) obsahuje druhé derivace energie
Minima na PES Problém minimalizace → nalezení pouze nejbližšího lokálního minima
Minima na PES – jak je najít Nederivační Metoda kroku pevné délky Simplexová metoda Derivační • První derivace Metoda největšího spádu Metoda konjugovaných gradientů • Druhá derivace Newton-Raphsonova metoda Quasi-Newtonova metoda
Minima na PES - konvergenční kritéria Všechny výše zmíněné metody se přibližují k minimu tak dlouho, dokud nejsou splněna konvergenční kritéria. Nastavení konvergenčních kritérií Počet kroků výpočtu Rozdíl mezi 2 kroky výpočtu: • ve struktuře • v energii
Metoda kroku pevné délky • Metoda systematicky prochází prostor kolem každé
souřadnice V každém kroku jsou pro každou souřadnici xi vygenerovány další dvě souřadnice (xi´ = xi + δxi, xi´´= xi + 2δxi). Pro obě nové souřadnice je spočtena Epot. Ze získaných bodů [xi, Epot(xi)], [xi ´, Epot(xi´)] a [xi´´, Epot(xi´´)] je vygenerována parabola. Je nalezeno minimum paraboly a pro další krok je xi nahrazeno tímto nalezeným minimem. • Jednoduchá na implementaci • Nejhorší minimalizační metoda zahrnující hodně kroků • Může zcela minout minimum • Příliš se nepoužívá
Simplexová metoda Nejedná se o simplexovou metodou v lineárním programování.
Simplex = geometrický útvar, který má v prostoru o M rozměrech M+1 vrcholů.
V okolí vstupního bodu vytvoříme iniciální simplex. Pro vrcholy simplexu vypočteme Epot a na základě porovnání této energie pro všechny vrcholy simplexu, simplex dále geometricky transformujeme (přetáčíme podle roviny, protahujeme, zkracujeme,… ). Cílem transformací je přesunout simplex z oblastí s velkou Epot do oblastí blíže minimu.
Simplexová metoda • Nejvýhodnější nederivační metoda • Vhodná speciálně v oblastech daleko od minima (slouží k
přiblížení k minimu) • V oblastech poblíž minima konverguje pomalu
Metoda největšího spádu Steepest descent method • Modeluje chování míčku na svahu
(míček padá tam, kde je největší spád → kde na něj nejvíce působí gravitační síla) • Postupuje tím směrem, kam směřuje opačný vektor gradientu (ve směru síly, působící v daném bodě).
Metoda největšího spádu Steepest descent method • Jednoduchá implementace i použití • Pro malé molekuly nejrychleji konvergující metoda • Poblíž minima konverguje pomalu • Může se otočit zpět nebo přejít minimum • Na dlouhých rovných plochách osciluje
Metoda konjugovaných gradientů Conjugate gradient method • Prohledávání není jen funkcí gradientu gk posledního (k-
tého) bodu cesty, ale také funkcí gradientu gk-1 předcházejícího bodu a směru sk-1 posledního přesunu. • Spojení informací o současném a předcházejícím stavu zabraňuje nepříjemnému efektu oscilování, který nám dělal problémy při využití metody největšího spádu.
Metoda konjugovaných gradientů Conjugate gradient method • Výpočetně náročnější než metoda největšího spádu, ale
spolehlivější a vhodná i v oblastech poblíž minima.
Newton-Raphsonova metoda • Nejjednodušší metoda, využívající druhé derivace.
Vyjadřuje potenciálovou funkci f pomocí Taylorova polynomu (pro bod xk): f(x) = f(xk) + (x - xk).f´(xk) + (x - xk)2.f´´(xk)/2 + ...
Tranzitní stavy na PES • Reaktanty a produkty chemické reakce jsou v rámci PES
sousedními lokálními minimy. • Reakční koordináta je v PES nejkratší cestou z minima, odpovídajícího reaktantům do minima, odpovídajícího produktům (IRC - intrinsic reaction coordinate). • Bod A je tranzitním stavem PES, pokud je v tomto bodě gradient funkce f nulový
Tranzitní stavy na PES Studium mechanismu chemické reakce: • Struktura aktivního komplexu • Potenciální energie aktivního komplexu • Průběh reakční koordináty Nalezení reakčních cest mezi dvěma nesousedícími minimy
Tranzitní stavy na PES Metody: • Linear synchronous transit (LST) • Quadratic synchronous transit (QST) • Saddle optimization method (SOM) • Locally updated planes (LUP) • Self penalty walk (SPW)
Tranzitní stavy na PES Linear synchronous transit (LST): • Vytvoří úsečku z reaktantu (R) do produktu (P). • Vypočítá Epot pro některé body (souřadnice atomů) této úsečky. • Jako tranzitní stav (TS) označí ten bod úsečky, pro který je hodnota Epot největší. Quadratic synchronous transit (QST): • Nejdříve pracuje stejně jako LST. • Z bodů R, P a TS sestrojí parabolu a poté na této parabole vyhledává opět maximum Epot (nový TS).
Saddle optimization method (SOM): • Vychází ze struktur P a R. • Zkouší na základě R podle vzoru P vygenergovat mezistrukturu (strukturu s co největší Epot), která se velmi podobá R, ale obsahuje již malé strukturní změny směrem k P. Touto strukturou R2 poté nahradí R. • Analogicky nalezne pro P odpovídající P2 a tím P nahradí. • Opakuje tento proces až do splnění konvergenčních podmínek.
Locally updated planes (LUP): • Na spojnici R a P nalezne maximum a toto maximum relaxuje • Při relaxaci využívá kolmici k nadrovině nalezeném maximu. Self penalty walk (SPW): • Reakční cesta je vyhledávána pomocí minimalizace průměrné energie podél dané cesty.
Globální minimum na PES • Lokální minimum s nejmenší hodnotou Epot na celé PES • Nejstabilnější uspořádání atomů a elektronů daného systému (nejvyšší pravděpodobnost výskytu v reálném chemickém prostředí). Metody • Systematické prohledávání (grid search) • Molekulová dynamika • Stochastické a Monte Carlo metody • Genetické algoritmy • Difuzní metody (DFT)
Systematické prohledávání (grid search): • Proloží hyperplochou mřížku a v jejích vrcholech vypočítá Epot. Takto zmapuje polohy lokálních minim a mezi nimi pak najde globální minimum. Molekulová dynamika Stochastické a Monte Carlo metody: • Začínají v nějaké vstupní geometrii (nejčastěji v lokálním minimu. Nové konfigurace generují náhodným posunem jednoho nebo více atomů (random kick).
Genetické algoritmy • Algoritmus je postaven na myšlence, že existují populace objektů, z nichž každý má svou množinu genů. Rodičovské objekty mohou tvořit potomky kombinací svých genů (přičemž může docházet i k mutacím – únik z lokálního minima). Nejplodnější jedinci z populace jsou vybíráni a přenášeni do další generace. Tito jedinci jsou také nejplodnější. Hledá se globální minimum účelové funkce. Difuzní metody • Potenciálová funkce je postupně měněna tak, že ubývají lokální minima, až zaniknou všechna s vyjímkou globálního (příspěvky ve směru kolmém k hyperploše => pro minima vzrůstá energie a pro maxima a sedlové body energie klesá)
MOLEKULOVÁ MECHANIKA
Proč molekulová mechanika? • Tam, kde selhává kvantová mechanika • Modely biomolekul, PES • Velmi efektivní pro studium dynamiky biosystémů →
molekulová dynamika (MD) • Oproti kvantové mechanice jednodušší → omezení na pouhou „Newtonovskou fyziku“ • Pomocí MM nelze studovat chemické reakce
30
Model atomu a vazby v MM • Atomové jádro a elektrony shrnuty do
jedné částice tvaru koule • Částice (koule) nese určitý parciální náboj, který je lokalizován v jejím středu • Částice mají určitý poloměr • Vazby, kterými jsou spojeny jednotlivé atomy jsou reprezentovány „pružinami“ • Atomy jsou v MM rozděleny do tříd (atomové typy): • chemická značka atomu • vaznost atomu a geometrické uspořádání vázaných atomů
Atomové typy aneb není uhlík jako uhlík Příklady atomových typů (silové pole AMBER): • C • • • • • • •
sp2 uhlík v uhlovodících CA sp2 uhlík v aromatickém kruhu CT sp3 uhlík v uhlovodících H vodík, vázaný na dusík HO vodík, vázaný na kyslík N sp2 dusík v aminoskupině NA sp2 dusík v pětičlenném cyklu N3 sp3 dusík v kladně nabitých aminoskupinách
Vazba jako pružina • Chemické vazby mají různou pevnost (disociační energii) • Analogie s pružinami, které mají různou tuhost
→ silová konstanta k [N.m-1] • Vazby mají různou délku → rovnovážná meziatomová vzdálenost r0 [pm] • k a r0 představují vazebné parametry pro biatomickou molekulu
Molekula vodíku
Principy Model prostředí: Silové pole (force field), které je tvořeno dvěmi složkami: • potenciálová funkce - funkce pro výpočet Epot (f: Rn → R)* • parametry - hodnoty konstant, vyskytujících se v potenciálové funkci * n je počet souřadnic, nutných pro popis geometrie molekuly. Pro molekulu s N atomy, N > 2 se jedná o 3N - 6 souřadnic.
Potenciálová funkce • Vytvořena na základě klasické Newtonovy mechaniky • Potenciální energie - součet energií všech vazebných i nevazebných interakcí v rámci molekuly
} }
Vazebné interakce
Nevazebné interakce 35
Vazebné a nevazebné interakce Vazebné interakce • popisují, jak moc se geometrie molekuly liší od „ideálního stavu“ • v rámci silového pole jsou definovány ideální délky vazeb, vazebné úhly a torzní úhly • čím více se určitý reálný parametr liší od ideálního, tím větší má daná interakce energii Nevazebné interakce • Elektrostatické - pro každou dvojici atomů určují, jak na sebe vzájemně působí náboje těchto atomů • Ne-elektrostatické - zbývající
Energie pnutí vazeb k 2 E (bond ) = (r − r0 ) 2
Harmonická aproximace
37
Deformace vazebných úhlů kθ 2 E (ang ) = (θ − θ0 ) 2
θ
Harmonická aproximace
39
Energie rotace vazeb (torze) E (tor ) =
kt (1 + cos(nτ − φ )) 2
τ
Van der Waalsovské interakce • Objeveny na základě studia chemických vlastností
zkapalněných vzácných plynů (helium, argon …). Tyto plyny jsou tvořeny izolovanými neutrálními atomy → je zřejmé, že mezi jejich atomy nepůsobí žádné vazebné nebo elektrostatické interakce. Ale přesto jsou tyto plyny kapalné (→ nějaké interakce je drží pohromadě). • Van der Waalsovské interakce jsou způsobeny třemi typy sil: coulombickými, indukčními a disperznímí
Van der Waalsova interakce • disperzní a repulzní složka • Lennard-Jonesův potenciál
⎛ ⎛ σ ⎞12 ⎛ σ ⎞6 ⎞ E ( LJ ) = 4ε ⎜ ⎜ ⎟ − ⎜ ⎟ ⎟ ⎜ ⎝ r ⎠ ⎝ r ⎠ ⎟ ⎝ ⎠ Repulze atomů
r σ ε
ε
Atrakce atomů
vzdálenost atomů vzdálenost, v níž na sebe atomy začínají silově působit nejmenší energie, kterou může daná dvojice atomů
Elektrostatická interakce • interakce nábojů dvou atomů • Coulombův zákon
E (coul ) =
qi q j 4πε 0rij
r vzdálenost atomů qi, qj vzdálenost, v níž na sebe atomy začínají silově působit ε0 permitivita vakua
• Parciální náboje • RESP (Restrained ElectroStatic Potential fit) • Mulliken
Parametry silového pole Každému atomovému typu přísluší parametry: Charakteristiky atomu - vdW poloměr, náboj atomu, relativní hmotnost atomu, vaznost a geometrie vazeb. Ideální parametry pro interakce atomu - délky vazeb, vazebné úhly a torzní úhly. Konstrainy - Silové konstanty, určující, jak moc je odchylka od ideálních parametrů nevýhodná. Původ parametrů • z experimentu - RTG a neutronová difrakce, NMR, rotační spektroskopie, vibrační spektroskopie (silové konstanty) • výpočtem - fitování energetických hyperploch
Molekula methanu
Software pro MM • AMBER http://amber.scripps.edu/
• CHARMM • GROMACS • OPLS
http://www.charmm.org/
• ENCAD • MOE
http://www.gromacs.org/
Optimalizace cyklohexanu (Hyperchem) a) Proveďte optimalizaci geometrie cyklohexanu v
židličkové konformaci. b) Proveďte optimalizaci geometrie cyklohexanu ve vaničkové konformaci. c) Porovnejte energetické rozdíly mezi oběma optimalizovanými konformacemi s publikovanými experimentálními daty. J.B. Hendrickson, J. Am. Chem. Soc. 83, 4537 (1961). J.B. Hendrickson, J. Am. Chem. Soc. 89, 7036 (1969).
Optimalizace cyklohexanu – židličková konformace • •
• • •
• • • • •
Volba silového pole AMBER, Distance Dependent Dielectric (škálovací faktor na 1), 1-4 škálovací faktory pro „Electrostatic“ a „van der Waals“ budou rovny 0.5. Volba „Cuttofs“ bude nastavena na „None“. Optimalizace budou provedena s „amber2“ parametry silového pole (Select Parameter Set). Stavba modelu cyklohexanu v židličkové konformaci Vytvořte šestiúhelník na jehož vrcholech jsou atomy uhlíku (může být nepravidelný, C spojeny jednoduchou vazbou!) a poté aktivujte „Add H and Model Build“. Změřte vazebné délky, úhly a torze a zaznamenejte. Výpočet Single Point energie cyklohexanu s následnou optimalizací Optimalizační algoritmus bude Conjugate Gradient (Polak-Ribiere), „RMS gradient“ nastavte na 0.1 kcal/mol Å, počet cyklů dejte na 300. Bude se jednat o „In vacuo“ optimalizaci. Zaznamenejte po optimalizaci energii i RMS gradient.
Optimalizace cyklohexanu – vaničková konformace • • • • • • • • • • • •
Volba silového pole Nechejte stejné silové pole jako v předchozím případě. Stavba modelu cyklohexanu ve vaničkové konformaci Zapněte „Multiple Selections“. Vyberte čtyři atomy C, které leží v jedné rovině (nástroj „Select“). Přejděte do „Name Selection“ a zvolte „PLANE“. Výběr rozšiřte o levou nebo pravou část cyklohexanu (zahrňte i vodíky). Aktivujte „Reflect“ v „Edit menu“. Změřte vzdálenost mezi nejbližšími axiálními vodíky. Výpočet Single Point energie cyklohexanu a následná optimalizace Proveďte výpočet energie pro danou vaničkovou konformaci cyklohexanu (Single Point) a zaznamenejte energii i RMS gradient. Optimalizační protokol ponechejte nastavený tak, jako v předchozím případě. Proveďte optimalizaci a zaznamenejte opět energii i RMS gradient.
