C7790 Počítačová chemie a molekulové modelování -1Petr Kulhánek kulhanek@chemi.muni.cz Národní centrum pro výzkum biomolekul, Přírodovědecká fakulta Masarykova univerzita, Kotlářská 2, CZ-61137 Brno 10. Plochy potenciální energie I C7790 Počítačová chemie a molekulové modelování I C7800 Počítačová chemie a molekulové modelování I - cvičení C7790 Počítačová chemie a molekulové modelování -2- Opakování t t itH    ),( ),(ˆ x x    časově závislá Schrödingerova rovnice C7790 Počítačová chemie a molekulové modelování -3- Opakování t t itH    ),( ),(ˆ x x    časově závislá Schrödingerova rovnice )()(ˆ xx kkk EH   časově nezávislá Schrödingerova rovnice )()(),( tft xx   systém se může nacházet v několika kvantových stavech stav je popsán vlnovou funkcí Y a má energii E C7790 Počítačová chemie a molekulové modelování -4- Opakování t t itH    ),( ),(ˆ x x    časově závislá Schrödingerova rovnice )()(ˆ xx kkk EH   časově nezávislá Schrödingerova rovnice )()(),( tft xx   ),()(),(ˆ RrRRr mmme EH YY )()(ˆ , RR llVRTlR EH   )(),()( RRrx  Y pohyb elektronů ve statickém poli jader elektronické vlastnosti systému pohyb jader v efektivním poli elektronů vibrace, rotace, translace Bornova- Oppenheimerova approximace C7790 Počítačová chemie a molekulové modelování -5- Opakování ),()(),(ˆ RrRRr mmme EH YY )()(ˆ , RR llVRTlR EH   pohyb elektronů ve statickém poli jader elektronické vlastnosti systému pohyb jader v efektivním poli elektronů vibrace, rotace, translace lVRTmoptmk EREE ,, )(  výsledná energie stavu elektronická složka energie vibračně, rotačně, translační složka energie optimální geometrie jader, při které je Em minimální C7790 Počítačová chemie a molekulové modelování -6- Opakování ),()(),(ˆ RrRRr mmme EH YY )()(ˆ , RR llVRTlR EH   pohyb elektronů ve statickém poli jader elektronické vlastnosti systému pohyb jader v efektivním poli elektronů vibrace, rotace, translace je možné obdobným způsobem dále rozdělit na samostatné příspěvky vibrační, rotační a translační kTjRiVlVRT EEEE ,,,,  C7790 Počítačová chemie a molekulové modelování -7- Nápověda ▪ vibrační energie je kvantována ▪rotační a translační stavy nebudeme uvažovat Molekula vodíku hvEV        2 1 Určete základní energie stavů 1s-1s a 1s+2s kvantové vibrační číslo 0,1,2,... Energie základního stavu 1s+1s: )0()(  vErEE Vo C7790 Počítačová chemie a molekulové modelování -8Struktura vs stav systému ✓ ✓ )0()(  vErEE Vo pouze část kvantově mechanického popisu stavu systému C7790 Počítačová chemie a molekulové modelování -9Domácí úkol I 1. Navrhněte vhodné energetické referenční stavy pro atom vodíku. 2. Navrhněte energetický referenční stav tak, aby byl stejný pro libovolný atom. 3. Odvoďte vztah pro excitační energii ze stavu n do stavu n+1 pro atom vodíku. 4. Navrhněte vhodný energetický referenční stav pro harmonický oscilátor. 5. Odvoďte vztah pro excitační energii ze stavu v do stavu v+1 pro harmonický oscilátor. 6. Navrhněte vhodný energetický referenční stav pro anharmonický oscilátor. 7. Může mít anharmonický oscilátor energii větší než De? 8. Proč může mít tuhý rotátor nulovou energii a harmonický oscilátor a částice v potenciálové jámě ne? 9. Srovnejte energie pro základní stav translačního, rotačního a vibračního pohybu molekuly vodíku. V případě translačního pohybu uvažujte objem boxu, který pojme 1 mol ideálního plynu za standardních podmínek. 10. Z jakého důvodu je Bornova-Oppenheimerova aproximace použitelná? 11. Kolikrát se prodlouží výpočet energie pokud se porovná výpočet pro molekulu benzenu s výpočtem molekuly bifenylu metodou CCSD(T). Každý vnitřní molekulový orbital (dva elektrony) je popsán jednou bázovou funkcí. Každý valenční molekulový orbital (dva elektrony nebo neobsazený) je popsán dvěma bázovými funkcemi. C7790 Počítačová chemie a molekulové modelování -10Domácí úkol II 1. Zamyslete se nad výhodou a nevýhodou reprezentace geometrie systému pomocí kartézských a interních souřadnic. C7790 Počítačová chemie a molekulové modelování -11- PES Plochy potenciální energie (Potential Energy Surface) Vlastnosti Vizualizace Významné body na PES C7790 Počítačová chemie a molekulové modelování -12Konfigurační prostor )(RE R = bod v 3N rozměrném prostoru (N je počet atomů) ),,,....,,,,,,( 222111 NNN zyxzyxzyxR kartézské souřadnice prvního jádra (atomu) Jednotlivé body tvoří konfigurační prostor. Každý bod v konfiguračním prostoru pak představuje unikátní strukturu daného systému. jednotlivé souřadnice jsou stupni volnosti daného systému C7790 Počítačová chemie a molekulové modelování -13Výpočet potenciální energie Výpočet potenciální energie E(R) je možný: ▪ aproximativním řešením Schrödingerovy rovnice (kvantová mechanika, QM) ▪ pomocí empirických silových polí (molekulová mechanika, MM) ▪ hybridním QM/MM přístupem ▪ pomocí zhrubených modelů C7790 Počítačová chemie a molekulové modelování -14Výpočet potenciální energie Výpočet potenciální energie E(R) je možný: ▪ aproximativním řešením Schrödingerovy rovnice (kvantová mechanika, QM) ▪ pomocí empirických silových polí (molekulová mechanika, MM) ▪ hybridním QM/MM přístupem ▪ pomocí zhrubených modelů pouze kategorie metod C7790 Počítačová chemie a molekulové modelování -15Výpočet potenciální energie Výpočet potenciální energie E(R) je možný: ▪ aproximativním řešením Schrödingerovy rovnice (kvantová mechanika, QM) ▪ HF metoda ▪ post HF metody (MPn, CI, CC) ▪ DFT metody (různé funkcionály) ▪ pomocí empirických silových polí (molekulová mechanika, MM) ▪ formy a parametry silových polí ▪ hybridním QM/MM přístupem ▪ rozhraní, typ QM-MM interakce, link atomy, ... ▪ pomocí zhrubených modelů stovky metod lišící se použitými aproximacemi neovlivňují obecné zákonitosti/vlastnosti E(R) )(RE C7790 Počítačová chemie a molekulové modelování -16Grafické zobrazení funkcí f(x,y) x y f(x) x funkce jedné proměnné funkce dvou proměnných C7790 Počítačová chemie a molekulové modelování -17Grafické zobrazení funkcí f(x,y) x y f(x) x funkce jedné proměnné funkce dvou proměnných z x y funkce tří proměnných barva je reprezentací funkční hodnoty f(x,y,z) volumetrické zobrazení C7790 Počítačová chemie a molekulové modelování -18Příklady volumetrického zobrazení Renderings of the atomic orbitals for a molecule (LiH, H is the white atom) with boundary enhanced volume contours. (a) shows the atomic orbitals rendered with only 1 s, 2 s, and 2 px of Li. (b) shows the atomic orbitals rendered with 2 py of Li, and 1 s of H on top of (a). Yun Jang; Varetto, U. Interactive Volume Rendering of Functional Representations in Quantum Chemistry. IEEE Transactions on Visualization and Computer Graphics 2009, 15, 1579–5186. C7790 Počítačová chemie a molekulové modelování -19Zobrazení E(R) 1 atom 2 atomy N atomů ),,( 111 zyxE ),,,,,( 222111 zyxzyxE ),,,....,,,,,,( 222111 NNN zyxzyxzyxE nezobrazitelné pouze volumetricky Příklad: enzym BsoBI má ~ 10000 atomů => 30000 stupňů volnosti, pro vizualizaci by bylo nutné použít 30000+1 dimenzionální prostor C7790 Počítačová chemie a molekulové modelování -20Vlastnost E(R) Potenciální energie je invariantní vůči: • posunutí (translaci) celého systému • natočení (rotaci) celého systému } bez působení vnějších silových polí (např. elektrostatické, magnetické, atd.) C7790 Počítačová chemie a molekulové modelování -21Invariance vůči posunutí HCHO )( 1RE )( 2RE )()( 21 RERE  21 RTR  ,....},,,,,{ TTTTTT zyxzyxT  vektor posunutí )( 1RE )( 2RE )()( 21 RERE  !!! Neplatí pro posun v silovém poli !!! C7790 Počítačová chemie a molekulové modelování -22Invariance vůči natočení HCHO )( 1RE )()( 21 RERE  21 RR Θ rotační matice )( 1RE )()( 21 RERE  !!! Neplatí pro rotaci v silovém poli !!! )( 2RE )( 2RE C7790 Počítačová chemie a molekulové modelování -23Dvouatomová molekula ),,,,,( 222111 zyxzyxE molekula vodíku ▪ tři translační stupně volnosti ▪ dva rotační stupně volnosti (molekula je lineární) C7790 Počítačová chemie a molekulové modelování -24Dvouatomová molekula ),,,,,( 222111 zyxzyxE )(rE meziatomová vzdálenost molekula vodíku 6-5=1 molekula vodíku ▪ tři translační stupně volnosti ▪ dva rotační stupně volnosti (molekula je lineární) C7790 Počítačová chemie a molekulové modelování -25Tříatomová molekula ),,,,,,,,( 333222111 zyxzyxzyxE molekula vody ▪ tři translační stupně volnosti ▪ tři rotační stupně volnosti C7790 Počítačová chemie a molekulové modelování -26Tříatomová molekula ),,,,,,,,( 333222111 zyxzyxzyxE molekula vody ▪ tři translační stupně volnosti ▪ tři rotační stupně volnosti ),,( 21 rrE 9-6=3 r1 r2  Interní souřadnice r1,r2, zobrazitelné pouze volumetricky C7790 Počítačová chemie a molekulové modelování -27Grafické zobrazení E(R) Funkce E(R) je běžně nezobrazitelnou funkci. Zobrazuje se z ní tedy pouze relevantní část v dvou či trojdimenzionálním prostoru, který co nejlépe vystihne studovaný problém. E(x,y) x y E(x) x )(1 Rfx  )(1 Rfx  )(2 Rfy  transformační (projekční) funkce C7790 Počítačová chemie a molekulové modelování -28Grafické zobrazení E(R) Funkce E(R) je běžně nezobrazitelnou funkci. Zobrazuje se z ní tedy pouze relevantní část v dvou či trojdimenzionálním prostoru, který co nejlépe vystihne studovaný problém. E(x,y) x y E(x) x )(1 Rfx  )(1 Rfx  )(2 Rfy  transformační (projekční) funkce Co ostatní stupně volnosti? C7790 Počítačová chemie a molekulové modelování -29Grafické zobrazení E(R) Funkce E(R) je běžně nezobrazitelnou funkci. Zobrazuje se z ní tedy pouze relevantní část v dvou či trojdimenzionálním prostoru, který co nejlépe vystihne studovaný problém. E(x,y) x y E(x) x )(1 Rfx  )(1 Rfx  )(2 Rfy  transformační (projekční) funkce Co ostatní stupně volnosti? )(2 Rrc f )(3 Rrc f0 )(    cr RE hodnota E(R) je vůči rc minimální C7790 Počítačová chemie a molekulové modelování -30Ilustrativní příklad C7790 Počítačová chemie a molekulové modelování -31Stacionární body x1 x2 x3 x1, x2 a x3 jsou stacionární body (lokální extrémy funkce). Odvozují se od nich vlastnosti kvantových stavů systému. lVRTk ExEE ,1)(  C7790 Počítačová chemie a molekulové modelování -32Stacionární body x3 tečna funkce E(x) ve stacionárním bodě má nulovou směrnici (a=0) Směrnice funkce je dána gradientem funkce (tj. první derivací funkce) x xE    )( )tan(a a Podmínka nutná pro stacionární bod 0 )(    x xE C7790 Počítačová chemie a molekulové modelování -33Stacionární body x3 tečna funkce E(x) ve stacionárním bodě má nulovou směrnici (a=0) Směrnice funkce je dána gradientem funkce (tj. první derivací funkce) x xE    )( )tan(a a Podmínka nutná pro stacionární bod 0 )(    x xE 0 )( 3    xx xE C7790 Počítačová chemie a molekulové modelování -34Typy stacionárních bodů lokální minima lokální maximum C7790 Počítačová chemie a molekulové modelování -35Určení typu stacionárního bodu   ... )( 2 1)( )()( 2 2 2        x x xE x x xE xExxE Taylorova řada: C7790 Počítačová chemie a molekulové modelování -36Určení typu stacionárního bodu   ... )( 2 1)( )()( 2 2 2        x x xE x x xE xExxE Taylorova řada: Ve stacionárním bodě:   ... )( 2 1 )()( 2 2 2     x x xE xExxE gradient je nulový kvadrát odchylky je vždy kladný C7790 Počítačová chemie a molekulové modelování -37Určení typu stacionárního bodu   ... )( 2 1)( )()( 2 2 2        x x xE x x xE xExxE Taylorova řada: Ve stacionárním bodě:   ... )( 2 1 )()( 2 2 2     x x xE xExxE gradient je nulový kvadrát odchylky je vždy kladný Hodnota funkce roste při vychýlení ze stacionárního bodu, pokud má druhá derivace v daném bodě kladnou hodnotu. Stacionární bod je pak lokálním minimem. Hodnota funkce klesá při vychýlení ze stacionárního bodu, pokud má druhá derivace v daném bodě zápornou hodnotu. Stacionární bod je pak lokálním maximem. 0 )( 2 2    x xE 0 )( 2 2    x xE C7790 Počítačová chemie a molekulové modelování -38Stacionární body 0 )( 2 2    x xE Lokální minimum: 0 )(    x xE Lokální maximum: 0 )( 2 2    x xE 0 )(    x xE !!! podmínka nutná !!! C7790 Počítačová chemie a molekulové modelování -39Dvourozměrný případ E(x,y) x y lokální minima C7790 Počítačová chemie a molekulové modelování -40Dvourozměrný případ E(x,y) x y lokální minima lokální maxima C7790 Počítačová chemie a molekulové modelování -41- sedlové body Dvourozměrný případ E(x,y) x y lokální minima lokální maxima C7790 Počítačová chemie a molekulové modelování -42Zobecnění pro E(R) Stacionární bod: 0 )(    R RE podmínka nutná, každá složka gradientu musí být nulová gradient má 3N složek Typ stacionárního bodu: )( )()( )()()( )()()( )()()()( 2 2 1 2 2 1 2 11 2 11 2 11 2 2 1 2 11 2 1 2 11 2 11 2 2 1 2 RH                                                   NN N z RE xz RE z RE yz RE xz RE zy RE y RE xy RE zx RE zx RE yx RE x RE      N počet atomů Charakter stacionárního bodu určuje Hessian, což je matice druhých derivací potenciální energie. Nezaměňovat s Hamiltonianem!!! C7790 Počítačová chemie a molekulové modelování -43Vlastnosti Hessianu kkk cHc  vlastní vektor (eigenvector) vlastní číslo (eigenvalue) Nk 3,...,1 N počet atomů • 6 (5) vlastních čísel je nulových – odpovídá translaci a rotaci systému • zbylá vlastní čísla: • všechna kladná – lokální minimum • jedno záporné, ostatní kladná – sedlový bod prvního řádu • dvě záporná, ostatní kladná – sedlový bod druhého řádu • ..... • všechna záporná – lokální maximum Diagonalizace Hessianu je způsob hledání vlastních čísel a vektorů. C7790 Počítačová chemie a molekulové modelování -44Vlastnosti Hessianu kkk cHc  vlastní vektor (eigenvector) vlastní číslo (eigenvalue) Nk 3,...,1 N počet atomů • 6 (5) vlastních čísel je nulových – odpovídá translaci a rotaci systému • zbylá vlastní čísla: • všechna kladná – lokální minimum • jedno záporné, ostatní kladná – sedlový bod prvního řádu • dvě záporná, ostatní kladná – sedlový bod druhého řádu • ..... • všechna záporná – lokální maximum Diagonalizace Hessianu je způsob hledání vlastních čísel a vektorů. chemicky významné stacionární body C7790 Počítačová chemie a molekulové modelování -45Diagonalizace Hessianu )( )()( )()()( )()()( )()()()( 2 2 1 2 2 1 2 11 2 11 2 11 2 2 1 2 11 2 1 2 11 2 11 2 2 1 2 RH                                                   NN N z RE xz RE z RE yz RE xz RE zy RE y RE xy RE zx RE zx RE yx RE x RE      )( )( 0 )( 00 0 )( 0 000 )( 2 2 2 3 2 2 2 2 2 1 2 Rλ                                   Nc RE c RE c RE c RE      Diagonalizace Hessianu je operace, při kterém se hledá takové natočení souřadného systému, při kterém jsou smíšené druhé derivace energie nulové. Nenulové mohou být pouze diagonální prvky matice. vlastní čísla Vlastní čísla Hessianu pak určují zakřivení funkce ve směru os nového souřadného systému. Tyto osy jsou určeny vlastními vektory, které jsou ortonormální. 1kcijji cc . ortogonální normalizační podmínka C7790 Počítačová chemie a molekulové modelování -46- Energie, Gradient, Hessian C7790 Počítačová chemie a molekulové modelování -47Výpočet potenciální energie Výpočet potenciální energie E(R) je možný: ▪ aproximativním řešením Schrödingerovy rovnice (kvantová mechanika, QM) ▪ HF metoda ▪ post HF metody (MPn, CI, CC) ▪ DFT metody (různé funkcionály) ▪ pomocí empirických silových polí (molekulová mechanika, MM) ▪ formy a parametry silových polí ▪ hybridním QM/MM přístupem ▪ rozhraní, typ QM-MM interakce, link atomy, ... ▪ pomocí zhrubených modelů stovky metod lišící se použitými aproximacemi C7790 Počítačová chemie a molekulové modelování -48Výpočet gradientu energie R R   )(E Gradient energie: )(RE                Nzzyx ,...,,, 111 jedná se o vektor, počet složek 3N               Nz E z E y E x E )( ,..., )( , )( , )( 111 RRRR Výpočet gradientu může být uskutečněn: • analyticky • numericky C7790 Počítačová chemie a molekulové modelování -49Analytický výpočet gradientu Analytický výpočet gradientu je preferovaný způsob výpočtu v případech, kdy je vyjádření a následný výpočet derivací energie snadný. Příklad: 2 0 )( 2 1 )( rrKE R energie dvouatomové molekuly v harmonické aproximaci, K a r0 jsou parametry modelu, r je meziatomová vzdálenost Řešení: viz domácí úkol C7790 Počítačová chemie a molekulové modelování -50Numerický výpočet gradientu Numerický výpočet gradientu je využíván tehdy, pokud není analytický gradient dostupný například z důvodu složitosti implementace algoritmu jeho výpočtu. K výpočtu numerického gradientu lze použít buď metodu dopředných diferencí (FD – forward differences) nebo centrálních diferencí (CD – central differences). V ojedinělých případech je možné použít i vícebodové metody. Metoda centrálních diferencí je přesnější než FD a tudíž preferovaným způsobem výpočtu gradientu. C7790 Počítačová chemie a molekulové modelování -51Numerický výpočet gradientu Dopředné diference Centrální diference x0x-1 x1x0 x1 h xEhxE xx xExE x E x )()()()()( 00 01 01 0        R h hxEhxE xx xExE x E x 2 )()()()()( 00 11 11 0         R C7790 Počítačová chemie a molekulové modelování -52Numerický výpočet gradientu Dopředné diference Centrální diference x0x-1 x1x0 x1 h xEhxE xx xExE x E x )()()()()( 00 01 01 0        R h hxEhxE xx xExE x E x 2 )()()()()( 00 11 11 0         R počítá se jednoupočítá se pro každou složku gradientu celkem 3N+1 výpočtů energie počítá se pro každou složku gradientu celkem 6N výpočtů energie C7790 Počítačová chemie a molekulové modelování -53Výpočet Hessianu )( )()( )()()( )()()( )()()()( 2 2 1 2 2 1 2 11 2 11 2 11 2 2 1 2 11 2 1 2 11 2 11 2 2 1 2 RH                                                   NN N z RE xz RE z RE yz RE xz RE zy RE y RE xy RE zx RE zx RE yx RE x RE      Hessian energie: jedná se o matici, počet složek 3Nx3N Výpočet Hessianu může být uskutečněn: • analyticky (paměťově a výpočetně náročné) • numericky (metodou centrálních diferencí) • z energií (3 x N x 3 x N x 2 výpočtů) • z gradientů (3 x N x 2 výpočtů) C7790 Počítačová chemie a molekulové modelování -54Hledání optimálních geometrií Hledání lokálních minim C7790 Počítačová chemie a molekulové modelování -55Optimální geometrie C7790 Počítačová chemie a molekulové modelování -56Optimalizační metody Metody optimalizace geometrie I. nultého řádu (pouze energie) ▪ simplexová metoda II. prvního řádu (pouze energie a gradient) ▪ metoda největšího spádu ▪ metoda konjugovaných gradientů III. druhého řádu (energie, gradient a Hessian) ▪ Newtonova metoda IV.pseudodruhého řádu (energie, gradient a aproximativní Hessian) ▪ Broydenova–Fletcherova–Goldfarbova–Shannova metoda (BFGS) C7790 Počítačová chemie a molekulové modelování -57- 1) Pro níže uvedenou funkci určete charakter bodů s hodnotami: a) x=1 b) x=0 c) x=-1 2) Za jaké situace může být druhá derivace funkce nulová? 3) Jaký je vztah mezi rozsahem reakce x a reakční koordinátou rc (také označovanou jako x)? 4) Vysvětlete proč hledáme tranzitní stavy reakcí jako sedlové body prvního řádu. Domácí úkol 33015)( 2  xxxE C7790 Počítačová chemie a molekulové modelování -58Domácí úkol 1. Vyjádřete gradient funkce E(R) podle kartézských souřadnic obou atomů. 2. Systém obsahuje 300 atomů. Výpočet jeho energie kvantově–chemickou metodou trvá 15 minut. Výpočet energie a analytického gradientu pak 20 minut. 1. Určete dobu výpočtu numerického gradientu a srovnejte ji s délkou výpočtu analytického gradientu. 2. Určete dobu výpočtu numerického Hessianu, který je počítán a) z energií a b) z analytických gradientů. 3. Navrhněte způsob urychlení výpočtu numerického gradientu a Hessianu. 2 0 )( 2 1 )( rrKE R C7790 Počítačová chemie a molekulové modelování -59Domácí úkol 1. Nastudujte uvedené optimalizační metody. Zaměřte se na jejich princip, výhody a nevýhody vůči ostatním optimalizačním metodám. 2. Srovnejte rychlost optimalizace geometrie molekuly (počet kroků) za použití stejné optimalizační metody pro dvě různé vyjádření geometrie pomocí interních a kartézských souřadnic. Literatura: (1) Leach, A. R. Molecular Modelling: Principles and Applications, 2nd ed.; Prentice Hall: Harlow, England; New York, 2001. (2) Jensen, F. Introduction to Computational Chemistry, 2nd ed.; John Wiley & Sons: Chichester, England; Hoboken, NJ, 2007.