Ugrás a tartalomhoz

A kvantumkémia alapjai és alkalmazása

Veszprémi Tamás, Fehér Miklós

Educatio Társadalmi Szolgáltató Nonprofit Kft.

7.3. Korrelációs számítások

7.3. Korrelációs számítások

Ahhoz, hogy az előző pontban bemutatott eljárások hatékonyságát és használhatóságát érzékelni tudjuk, hasonlítsuk először össze egy FCI számítás eredményeit a HF, CI, CC, MP, CAS, MR és CASPT2 módszerek eredményeivel. Célpontunk a vízmolekula, az alkalmazott bázis minden esetben a cc-pVDZ szint. A számításokat R. J. Bartlett és munkatársai végezték el állandó OH távolság (Re=0,9755 Å), és HOH kötésszög (110,60) mellett.

J. Olsen, P. Jorgensen, H. Koch, A. Balkova, R.J. Bartlett, J. Chem. Phys., 104 , 8007 (1996).

A disszociáció hatásának tanulmányozására a számításokat megismételték 2Re és 3Re kötéstávolságokkal is. Az eredményeket a 7.1. táblázatban foglaltuk össze. A táblázat tartalmazza a teljes energiákat, és az FCI eredményekhez viszonyított energiakülönbségeket (Δ). Minthogy a molekula geometriáját nem optimálták, az eredmények nem mérhetők össze a 6.10. táblázat adataival.

Jól látható, hogy az eltérés (Δ) a gerjesztési szint növekedtével – a várakozásnak megfelelően – csökken. A csökkenés a CC módszer esetén a leginkább rohamos, így pl. Re mellett a CCSDT módszer hibája 8-szor kisebb, mint a CCSD-é, a CCSDTQ módszeré pedig mintegy 26-szor kisebb a CCSDT hibájánál. A konvergencia sokkal gyorsabb, mint a CI-konvergencia: a CCSD energia mélyebb, mint a CISDT-vel számított, majdnem eléri a CISDTQ értékét. Az is jól látható azonban, hogy noha disszociáló molekula esetén is működik a módszer, a konvergencia már jóval lassabb.

A M  øller–Plesset perturbációs számítások eredményei komoly mértékben függenek a kiindulási geometriától. Az egyensúlyi Re távolságnál a Hartree–Fock-determináns dominál, ezért a perturbációs sor gyorsan konvergál. A molekula disszociációja során azonban a helyzet változik, az alapállapot determinánsának súlya csökken és már nem nevezhető az eltérés „kis perturbációnak”. Ennélfogva a perturbációs sor konvergenciájában zavarok lépnek fel: 2Re-nél oszcillációt tapasztalunk, 3Re-nél pedig egy értelmetlen divergáló sort kapunk. Meg kell jegyeznünk, hogy az oszcilláció jelensége számos esetben fellép, még stabilis molekulák esetében is. Mivel a módszer nem variációs, nem meglepő, hogy a számított energia esetenként az FCI szint alá is mehet – természetesen ez sem tartozik a módszer előnyei közé.

 

7.1. táblázat. Energiák és energiakülönbségek az FCI és más módszerek között (a.e.-ben) különböző OH távolságok mellett.

 

   

A táblázat a teljes energiák 1 -szeresét tartalmazza.

Két különböző aktív tér hatását láthatjuk a CAS, MR és CASPT2 számítások eredményeinek az összehasonlítása során. A kisebb (a-val jelölt) aktív tér három betöltött és három virtuális pályából áll, a nagyobb (b) 3 betöltött és 4 virtuális pályát tartalmaz. A viszonylag gyenge CAS eredmények azt sugalják, hogy a dinamikus korreláció megfelelő kezeléséhez jóval nagyobb aktív tér lenne szükséges. A korrekció egyik lehetséges módja az MR számítás. Az MRSDCIa és MRSDCIb számításokban az eredeti a és b „referencia” teret kiegészítették egyszeresen és kétszeresen gerjesztett determinánsok sorozatával. Ezáltal igen jóminőségű hullámfüggvényeket kaptak, melyekkel nemcsak az alapállapotot, hanem a disszociáció is jól leírható. A b referencia a hibát az a-hoz képest mintegy a harmadára csökkenti.

A CAS számításokban hiányolt dinamikus korreláció figyelembe vételének másik lehetséges módja a perturbációs számítás. A CASPT2 számítások eredményei jók, noha nem érik el az MR számítások minőségét. Érdekességük, hogy nem javulnak az eredmények lényegesen a nagyobb (b) aktív tér hatására. Az MR és a CASPT2 módszerek esetén is megfigyelhető, hogy a számítások hibája nagyjából azonos a teljes potenciálfelületen.

 

7.2. táblázat. A vízmolekula néhány számított és mért tulajdonsága. #+

 

   

# Hajgató Balázs számítási eredményei.

+ A számításokban a cc-pVTZ bázist használtuk.

+ + a = adiabatikus, v = vertikális ionizációs energia (l. 6 . fejezet 10 . megjegyzése).

A CAS-számításokban az aktív tér mérete (8,8).

A CASPT2-, CASPT3- és MRSDCI-számításokban a CAs (8,8) hullámfüggvényt használtuk.

A 7.1. táblázat alapján már lehet némi elképzelésünk arról, hogy a korrelációs energia mekkora hányadát tudjuk figyelembe venni az egyes módszerekkel, illetve az adott módszer bizonyos szintjén. A 7.2. és 7.3. táblázatokban molekuláris sajátságok korrelációs módszerekkel számított értékeit hasonlítjuk össze. Mindkét táblázat folytatás: a 7.2 táblázatban a vízmolekulára végzett számítások eredményei találhatók és összevethetők a „sima” HF számítások 6.11-beli eredményeivel. Az adatokból látható, hogy a kísérleti energiától még eléggé távol vagyunk, de a számított geometriai adatok igen pontosak. (Ismételten hangsúlyozzuk, hogy az energiák nem összemérhetők a 7.1. táblázat adataival.) A kötésszög 1-2 fokos eltérése a kísérleti értéktől elsősorban a kísérleti és számított szerkezeti adatok különböző jelentésével értelmezhető. Az ionizációs energiákat nem a Koopmans-tétellel számítottuk, hanem az alap- és ionállapot energiájának különbségéből. (Mivel post-HF számításokban a molekulapálya nem értelmezhető, a Koopmans-tételt sem használhatjuk.) Külön számítottuk a vertikális ionizációs energiákat (amikor az ionállapot geometriája nem változik az alapállapotéhoz képest), és az adiabatikus ionizációs energiákat (amikor az ionállapotban a molekula relaxálódik). Nem vettük viszont figyelembe az alap- és ionállapotokban a zéruspont energiát, a kísérleti értékektől való eltérésnek ez a fő forrása.

A 7.3. táblázat értelmezéséhez idézzük emlékezetünkbe a 6.5. táblázat következtetéseit: amíg a metil-izocianát molekula szerkezetét HF szinten kis bázis segítségével is megnyugtatóan tisztázni tudtuk, a CH3NCO leírásához nagy bázisra, és mindenek előtt polarizációs függvényekre volt szükség. A nagyon flexibilis CH3NCS térszerkezetének modellezésében a HF módszer kudarcot vallott. Még a legnagyobb bázis és több sorozat polarizációs függvény használata mellett is C3v szimmetria, azaz 180 fokos C–N–C szög adódott, noha a kísérletek a hajlott szerkezetet egyértelművé tették. A korreláció figyelembe vétele ugrásszerűen javítja az eredményeket és már viszonylag kis bázissal is elfogadható egyezést kapunk.

Mindez ideig azt vizsgáltuk, hogy a számítás szintjének emelése hogyan hat a számítási eredményekre. Leszögezhetjük a nyilvánvaló tényt: általában javítja azt.

Aki az ellenkezőjére tippelt, itt abbahagyhatja...

Persze, mint minden válasz, ez sem egészen pontos. Az energia (variációs módszer esetén) valóban mindig javul, de más molekuláris sajátságokra ez már nem mindig áll fenn. Fordítsunk most egyet a kérdésen. Az elmélet egy adott szintjén dolgozva, vajon milyen messze vagyunk az egzakt megoldástól, vajon mekkora az elkövetett hiba? A válasz bizonytalan, hiszen az elektronkorreláció függ a vizsgált rendszer méretétől, szerkezetétől, valamint a kérdéses tulajdonságtól és nem tudjuk előre megmondani, hogy pl. egy adott molekula esetében újabb gerjesztések figyelembe vétele nem változtatja-e meg lényegesen az addigi következtetéseinket. Úgy tűnik (de az ellenvéleményt nem zárhatjuk ki), hogy nem túl „faramuci” molekulák esetében az S és D típusú gerjesztések figyelembe vétele a molekula legtöbb tulajdonságának leírásához elegendő. Általában minél flexibilisebb egy molekula, annál nehezebb a kezelése. Ugyancsak magas szintű módszerekre van szükség gyenge kötések, illetve gyengén kölcsönható rendszerek (pl. van der Waals-kötések, van der Waals-molekulák) tárgyalásához. E könyv második részében látni fogjuk, hogy nehéz a gerjesztett állapotok tárgyalása, és a molekulaszerkezet megváltozásának a kvantumkémiai kezelése, pl. kémiai reakciók leírása. Annyi viszont biztos, hogy emelve a számítás szintjét az eredmények egyre kevésbé változnak, és szükségképpen konvergálnak az FCI szintjéhez. Az elmondottak illusztrálására lássunk egy példát!

 

7.3. táblázat. A CH3NCO és CH3NCS számított és mért szerkezeti adatai.

 

   

 

7.4. táblázat. A BrCNO néhány számított szerkezeti adata. ++

 

   

+ Geometriai adatok  Å -ben és fokban, energia atoegységben, Δ E c m 1 -ben, az alkalmazott bázis 6-31G*.

+ + T. Pasinszki, N. Westwood, J. Phys. Chem. , 99 , 6401–6409 (1995) nyomán.

+ + + A linearitáshoz vezető energiagát.

A BrCNO rövid élettartamú reaktív molekula, mely ezért csak speciális körülmények között és speciális módszerekkel tanulmányozható. A molekula és rokon molekulák térszerkezetének tanulmányozása egyértelműen azt jósolja, hogy a molekula egyensúlyi szerkezete lineáris, vagy ha hajlott is, a linearitáshoz vezető gát alacsony, várhatóan kisebb mint 50 cm1. A molekula térszerkezetének kvantumkémiai vizsgálata nagyon tanulságos példát szolgáltat a különböző számítási módszerek hatékonyságára és hiányosságaira. A 7.4. táblázatból látható, hogy a HF, MP3, QCISD, CCSD számítások a molekula szerkezetét lineárisnak, a többi módszer hajlottnak jósolja. Ha a perturbációs számításokat nézzük (HF = MP1), lineáris-hajlott-lineáris-hajlott váltakozás figyelhető meg. A páros rendű tagoknál (MP2, MP4...) mindig újfajta gerjesztéseket veszünk figyelembe (az MP2-nél belépnek a kétszeres, az MP4-nél az egyszeres, háromszoros és négyszeres gerjesztések), az utánuk következő páratlan tagoknál pedig ezek az új tagok (pl. MP3-nál a kétszeres gerjesztések) csatolódnak egymással. Mindennek az a következménye, hogy a páros tagok az újonnan figyelembe vett gerjesztések hatását kissé eltúlozzák, míg a páratlan tagok ezek hatását alulbecsülik (pl. az MP2 módszer eltúlozhatja, az MP3 módszer alábecsülheti a kétszeres gerjesztések hatását). Tulajdonképpen ez látható a táblázatban. A háromszoros gerjesztések hatását a különböző módszerek eltérően veszik figyelembe és leírásuk egyre pontosabb az MP4–QCISD(T)–CCSD(T) sorban. Ez nagyon jó összhangban van a számított gátmagasság csökkenésével. Jól látszik azonban, hogy a gátmagasság pontos leírásához az alkalmazott és meglehetősen magas szint még mindig nem elegendő, további gerjesztések figyelembe vételére, illetve az alkalmazottnál jóval nagyobb bázis használatára lenne szükség. Ezzel várhatóan még tovább csökkenne a gát magassága. Mindenesetre a feltételes mód nem vígasztalhat bennünket: az adott szinten a feltett kérdésre egyértelmű választ nem tudunk adni.

Az elmondottakból következik, hogy ismerjük az elvi algoritmust a Schrödinger-egyenlet egzakt megoldására. Semmi mást nem kell tennünk, csak a bázist és a korrelációs szintet növelni a végtelenségig (7.5. táblázat). Miután pedig ezen kimosolyogtuk magunkat, gondoljuk végig, milyen tényleges lehetőségek állnak rendelkezésünkre!

 

7.5. táblázat. Bázis és a korreláció szintje: a számítási szint definiálásának két fő faktora.

 

   

Az egyre magasabb szintek egyre költségesebbek. A számítógép szükségletet becsülendő nézzük az egyes módszerek műveletigényét. Ha n a betöltött és N a virtuális spinpályák száma, akkor – figyelembe véve az elvileg elvégzendő integrálások számát, a mátrix diagonalizálás, és a fellépő egyéb műveletek költségét – a műveletigényt az egyes módszerek alkalmazása során meg tudjuk becsülni (7.6. táblázat). Amint látható, a számítási szint emelése a műveletigény és így a gépidő szükséglet nagyságrendekkel való növekedését eredményezi. (A memória- és háttértár-szükségletről még nem is beszéltünk.) Ráadásul, amíg egy HF számításhoz általában DZP minőségű bázis már megfelelő méretű, korrelációs számításokhoz ennél nagyobb bázisra van szükség. (Ez abból a tényből következik, hogy e módszerekben a virtuális pályák szerepe azonos a betöltött pályákéival.) A bázis kiválasztásának két kritériuma, a minőség és a méret nyilvánvalóan sugalja, hogy a bázist is a korrelációs módszerekhez kell illeszteni. A 6.5. pontban említett „korreláció-konzisztens” bázis már ilyen.

 

7.6. táblázat. A pontosság és használhatóság egyensúlya.

 

   

+ Csak abban az esetben érvényes, ha N sokkal nagyobb, mint n.

+ + Nem-hidrogés atomokat figyelembe véve.

Természetesen lehetőségeinket az alkalmazott számítógép kapacitása határozza meg. E határokon belül kell megtalálni a feladat megoldásához az optimális bázist és számítási szintet, és ezek fényében gondolkodhatunk el eredményeink megbízhatóságán is. Mindemellett bízhatunk a számítástechnika töretlen fejlődésében is, amint azt a 7.7. táblázatbeli összehasonlítás mutatja. A triamino-trinitro-benzol RHF/6-31G** számításához szükséges gépidő igényt látjuk különböző számítógépeken, különböző programok használatával. Összesen 300 bázisfüggvényünk van és az összehasonlítás egyetlen HF számításra vonatkozik.

Micsoda döbbenetes fejlődés! Amihez nem is olyan régen még sok-sok évnyi munka kellett (volna), ma néhány másodperc alatt lenyelhető. És még mi minden várható! Az ún. Moore-törvény szerint a számítógépek úgy másfél-kétévenként megduplázódik. A táblázat adatai azt is bizonyítják, hogy nem csupán a hardware fejlődik rendületlenül, hanem a számítási módszerek és algoritmusok hatékonysága is rohamosan nő. Remek kilátások, rózsaszín a világ! Jelenleg rutin HF-számításokat végezhetünk kb. 100 atomos rendszereken, még néhány év, és akár fehérjéket is számíthatunk CCSDT szinten! Vagy talán mégsem? Tételezzük fel, hogy a mainál 1000-szer gyorsabb számítógép áll a rendelkezésünkre. Mivel a CCSDT módszer a pályák számának 8. hatványától függ, ezzel a géppel az adott szinten a mainál – 100082 – kb. kétszer nagyobb molekulát számíthatunk.

Eredetileg Gordon Moore a mikroelektronikai eszközök tárolókapacitásának növekedését becsülte a fenti módon.

 

7.7. táblázat. A triamino-trinitro-benzol RHF/6-31G* számításához szükséges processzoridők. +

 

   

+ Gaussian News, Summer 1995, p.9.

+ + Szieberth Dénes számításai.

Mintha szertefoszlana az előbbi eufória! Ha egy fehérje mondjuk 10000 atomot tartalmaz... . bizony, ilyen számítás belátható időn belül nem végezhető el.

Szerencsére a fejlődés minden irányban töretlen. Nézzük például, hogyan lehet a 7.6. táblázat elvileg korrekt adatait számunkra előnyösen „megviccelni”. A független részecske modellből következik, hogy a számítások szűk keresztmetszete az integrálok előállítása. Az egyelektronos integrálok száma elvileg ugyan a molekulaméret (ill. a bázisfüggvények számának) négyzetével arányos, ám a Gauss-függvények szorzata is Gauss-függvény. A szorzat-Gauss-függvény „mérete” pedig exponenciálisan csökken a kiindulási függvények távolságával (6.5. ábra), és úgy 8-10  Å távolságban az eredő függvény integrálja már elhanyagolható. Ebből következően minél nagyobb a molekuláris rendszer, annál több integrál esik ki és így – nem varázslat – az egyelektronos integrálok száma máris csupán lineárisan függ a bázisfüggvényektől, és a többi integrálás költségén csupán egy-egy szorzást kell elvégezni. Hasonló gondolatmenetet követve kiderül, hogy a kiértékelendő kételektron integrálok száma a bázisfüggvények számának negyedik hatványáról aszimptotikusan csökken a második hatványra. Az „aszimptotikus” jelző azt jelenti, hogy minél nagyobb a rendszer, annál inkább közeledik a második hatvány felé, mivel annál több nagyon távoli kölcsönhatást, illetve az ezeket képviselő integrált lehet elhanyagolni. Ez számunkra különösen kedvező, hiszen nem is a kis rendszerekkel van bajunk. Különböző trükkökkel még további komoly sebességnövekedés érhető el. Az ún. FMM módszer (Fast Multipole Method) a közeli és távoli kölcsönhatásokat megkülönbözteti, a közelieket analitikai integrálással számítja, a távoliakat viszont bizonyos algoritmus szerint közelíti.4 Gyorsítani lehet a mátrix diagonalizálás alapvetően n3-ös méretfüggését is, pl. annak a kihasználásával, hogy a Hamilton mátrix (főként post-HF módszerek esetén) nagyon szellős, azaz nagyon sok zérus elemet tartalmaz. Mindezek nyomán ma már léteznek olyan algoritmusok HF számításokra, melyek számítógép-igénye a bázisfüggvények számával (aszimptotikusan) csupán lineárisan változik. Máris változik a kép, hiszen a 7.6. táblázatban megadott maximális 100-200 atomos rendszer egyszeriben ezres nagyságrendűvé alakul! Ugyancsak csökkenthető a korrelációs módszerek gépidő-igénye a 7.6. táblázat adataihoz képest. Pulay bebizonyította, hogy lokalizált pályákat használva az MP2 módszer gépidő-igénye a pályák számával eredeti 5. hatvány helyett aszimptotikusan a lineáris felé változik. Bíztató kísérletek léteznek a CCSD módszer 1-2 nagyságrenddel való gyorsítására is.4 Jó reményeink vannak tehát akár biomolekulák kielégítő kvantumkémiai vizsgálatára is – ha nem is FCI szinten.