C7800 Počítačová chemie a molekulové modelování I - cvičení Petr Kulhánek, Jakub Štěpán kulhanek@chemi.muni.cz Národní centrum pro výzkum biomolekul, Přírodovědecká fakulta Masarykova univerzita, Kotlářská 2, CZ-61137 Brno Obsah > Požadavky na zpracování výsledků > Tématické okruhy > Molekula vody > Dimer molekuly vody > Interakce halidů s bambus[6]urilem > Referenční manuál > Gaussian > Infinity > Nemesis > Avogadro > VMD > Extrapolace na CBS [Požadavky na zpracování výsledků Výsledky jednotlivých cvičení budou zpracovány do protokolu, který bude mít následující náležitosti: • Jméno a příjmení, název cvičení a datum • Pro každý tematický okruh: • Stručné shrnutí tématu včetně reakčního schématu, pokud je to vhodné • Použitý software včetně verzí • Výsledky (tabulky) • Diskuze výsledků dle zadání • Použitá literatura (např. u experimentálních hodnot) Protokol ve formátu pdf je nutné odevzdat do 9.11. 24:00 do odevzdávárny CV1 očítačová chemie a molekulové modelování I - cvičení -3- Molekula vody očítačová chemie a molekulové modelování I - cvičení Úkoly 1) Vytvořte model molekuly vody a jeho geometrii zoptimalizujte pomocí molekulové mechaniky. 2) Zoptimalizujte geometrii molekuly vody pomocí metody HF/cc-pVDZ 3) Změřte významné geometrické parametry optimalizované geometrie a srovnejte je s výchozím modelem. 4) Ověřte, že nalezená geometrie odpovídá lokálnímu minimu na PES, pomoci vibrační analýzy. 5) Pro optimalizovanou geometrii proveďte výpočet energie včetně výpisu průběhu SCF výpočtu, dipólového momentu a Mullikenových a MK (Merz-Singh-Kollman) atomových nábojů pomocí metody HF/cc-pVDZ 6) Na stejné geometrii opakujte výpočet uvedený v bodě 4 pro báze: cc-pVTZ, cc-pVQZ a cc-pV5Z 7) Energii molekuly vody extrapolujte na CBS (Complete Basis Set), tj. určete hodnotu energie pro metodu HF/CBS 8) Extrapolujte dipólový moment na CBS hodnotu. 9) Srovnejte hodnotu vypočteného dipólového momentu s experimentální hodnotou, případný rozdíl diskutujte. Řešení Štábní kultura 00.input 01.opt 02. freq 03. props 1) Počáteční geometrii molekuly vody vytvořte v programu Avogadro nebo Nemesis. Geometrii modelu předoptimalizujeme pomocí molekulové mechaniky (silového pole). Zvolte takové silové pole, které dle vašeho názoru nejlépe vystihne geometrii molekuly vody. Předoptimalizované geometrie uložte ve formátu xyz pod názvem water.xyz do složky OO.input 2) Soubor water.xyz překopírujte do adresáře Ol.opt a přejmenujte jej na opt.com. Přejmenovaný soubor otevřete v textovém editoru a změňte jej na vstupní soubor pro optimalizaci geometrie v programu Gaussian. Spusťte výpočet. očítačová chemie a molekulové modelování I - cvičení 01. cc-pVDZ 02. cc-pVTZ 03. cc-pVQZ 04. cc-pV5Z Řešení 3) V adresáři Ol.opt otevřete soubor opt.log v textovém editoru a analyzujte jeho obsah. Ověřte, že výpočet proběhl v pořádku a nalezněte optimalizovanou geometrii. Soubor zavřete. Ze souboru opt.log vyextrahujte informace o změně energie, dále odpovídající geometrie a optimalizovanou geometrii do souboru s názvem last.xyz pomocí skriptů z modulu qmutil. Soubor last.xyz otevřete v programu vmd, Avogadro, či Nemesis a změřte významné geometrické parametry. 4) Soubor last.xyz překopírujte do adresáře 04.freq a přejmenujte jej na freq.com. Přejmenovaný soubor otevřete v textovém editoru a upravte jeho obsah pro výpočet normálních vibrací na úrovni teorie HF/cc-pVDZ. Spusťte výpočet. Výstupní soubor freq.log otevřete v textovém editoru a určete typ stacionárního bodu. Normální vibrace zobrazte v programu Avogadro nebo Nemesis. 5) Soubor last.xyz překopírujte do adresáře 03.props/01.cc-pVDZ a přejmenujte jej na props.com. Přejmenovaný soubor otevřete v textovém editoru a upravte jeho obsah pro výpočet energie metodou HF/cc-pVDZ. Zvolte úplný výpis (#P) a výpočet ESP atomových nábojů metodou Merz-Singh-Kollman (Pop=MK). Spusťte výpočet. Výstupní soubor props.log otevřete v textovém editoru a vyextrahujte z něj data do dále uvedené tabulky. 6) Postupujte analogicky jako v bodě 5, použijte postupně metody: HF/cc-pVTZ, HF/cc-pVQZ, HF/cc-pV5Z ) Počítačová chemie a molekulové modelování I - cvičení Řešení 6) (aug-)cc-pVXZ jsou konzistentně korelační báze, které vykazují exponenciální pokles energie systému se vzrůstající velikostí báze. Tento fakt lze využít k extrapolaci energie na úplnou bázi (CBS - complete basis set neboli HF-limitu). E(x) = ECBS + Ae kde x je číselné vyjádření pro X: 2 (D), 3 (T), 4 (Q) ECBS, A, B jsou neznámé konstanty, které lze učit z hodnot energie pro tři hodnoty x (báze). 7) Podobná závislost platí i pro jiné vlastnosti systému, např. dipólový moment, frekvence vibrace, geometrický parametr, apod. P(x) = P_ + Ae ' -B x CBS C7800 Počítačová chemie a molekulové modelování I - cvičení Výsledky Geometrie molekuly vody doplňte použité silové pole Metoda MM HF/cc-pVTZ Rozdíl d(HO) [Á] 0(HOH) □ Výsledky Molekula vody - průběh SCF metody HF/cc-pVDZ SCF Krok E [au] AEpre„ [au] 2 3 4 doplňte podle skutečného počtu SCF kroků změna energie vůči předchozímu kroku Poznámka: Hodnoty energie vyextrahujte ze souboru manuálně nebo použijte příkaz grep. -10- očítačová chemie a molekulové modelování I - cvičení Výsledky Molekula vody Mul liken MK ESP Báze E [au] Er [kcal/mol] u-[D] q(H) q(O) q(H) q(O) cc-pVDZ rr.n\/T7 0,0 LL |J V 1 Ĺ. cc-pVQZ cc-pV5Z CBS t ---- 4 ---- ---- výsledek výpočtu absolutní energie (E(RHF)) relativní energie vůči bázi cc-pVDZ atomové náboje dipólový elektrostatický moment C7800 Počítačová chemie a molekulové modelování I - cvičení Dimer molekuly vody Úkoly 1) Vytvořte model dimeru molekuly vody a jeho geometrii zoptimalizujte pomocí molekulové mechaniky. 8 Zoptimalizujte geometrii dimeru molekuly vody pomocí metody HF/cc-pVDZ Změřte významné geometrické parametry optimalizované geometrie a srovnejte je s výchozím modelem. Pozorované rozdíly se pokuste zdůvodnit. Ověřte, že nalezená geometrie odpovídá lokálnímu minimu na PES, pomocí vibrační analýzy. Pro optimalizovanou geometrii proveďte výpočet energie se započítáním BSSE korekce. Na stejné geometrii opakujte výpočet uvedený v bodě 4 pro báze: cc-pVTZ, cc-pVQZ a cc-pV5Z Energii dimeru molekuly vody extrapolujte na CBS (Complete Basis Set), tj. určete hodnotu energie pro metodu HF/CBS Pro každou bázi určete interakční energii mezi molekulami vod. Srovnejte ji s hodnotou vypočtenou s použitím úplné báze. Srovnejte hodnotu interakční energie s hodnotou BSSE chyby. Výsledky diskutujte. očítačová chemie a molekulové modelování I - cvičení -13- Řešení Štábní kultura 00.input 01.opt 02. freq 03. props 01. cc-pVDZ 02. cc-pVTZ 03. cc-pVQZ 04. cc-pV5Z Postupuje se analogicky jako v předchozím případě. Do výpočtu energie je však nutné zahrnout korekci na BSSE (Counterpoise=2 + definice fragmentů). očítačová chemie a molekulové modelování I - cvičení -14- Výsledek Geometrie dimeru molekuly vody doplňte použité silové pole | Metoda MM HF/cc-pVTZ Rozdíl d(HO) [Á] 0(HOH) □ d(H...O) [Ä] vodíková vazba Dle vlastního uvážení uveďte další geometrické parametry, které nejlépe postihnou rozdíl mezi oběma geometriemi. očítačová chemie a molekulové modelování I - cvičení -15- Výsledek Dimer molekuly vody Báze E ^BSSE ^BSSE E, [au] [kcal/mol] [au] [kcal/mol] [kcal/mol] cc-pVDZ 0,0 cc-pVTZ cc-pVQZ cc-pV5Z CBS výsledek výpočtu absolutní energie včetně korekcef na BSSE interakční energie mezi molekulami vody E. = E -2* E i dun er monomer relativní energie vůči bázi cc-pVDZ Vazba aniontů do kavity BU[6] očítačová chemie a molekulové modelování I - cvičení -17- Úkoly Namodelujte geometrii komplexu BU[6] s chloridovým aniontem. Zoptimalizujte geometrii komplexu pomocí semiempirické kvantově chemické metody PM3. Za použití optimalizované geometrie připravte model komplexu BU[6] s bromidovým aniontem a samotný BU[6]. Zoptimalizujte geometrie komplexu Br@BU[6] a BU[6] metodou PM3. Ověřte, že nalezené geometrie odpovídají lokálním minimům na PES, pomocí vibrační analýzy. Vypočítejte vazebnou energii chloridového a bromidového aniontu do BU[6]. Výsledky srovnejte s experimentálními hodnotami vazebných afinit Zobrazte elektrostatický potenciál namapovaný na povrch elektronové hustoty pro molekulu BU[6]. Zobrazení použijte k vysvětlení vazby aniontů do kavity bambusurilu. Výpočet interakční energie zopakujte pomocí metody PM6. Získané výsledky srovnejte s metodou PM3. očítačová chemie a molekulové modelování I - cvičení -18- 7800 Počítačová chemie a molekulové modelování I - cvičení Řešení 7) Z adresáře 02.PM3/03.BU6/01.opt překopírujte soubor last.xyz do adresáře 02.PM3/03.BU6/03.surface a přejmenujte jej na surface.com. Přejmenovaný soubor otevřete v textovém editoru a upravte jeho obsah pro výpočet energie na úrovni teorie HF/cc-pVDZ. V souboru uveďte jméno archivu (checkpointu) jako surface.chk. Spusťte výpočet. Po dokončení výpočtu vypočítejte elektronovou hustotu a elektrostatický potenciál pomocí utility cubegen z modulu gaussian. $ module add gaussian $ formchk surface.chk $ cubegen 0 density=scf surface.fchk density.cube 0 h $ cubegen 0 potential=scf surface.fchk esp.cube 0 h Soubory density.cube a esp.cube zobrazte v programu VMD. Vytvoří formátovaný archiv, který je využíván v dalších analýzách. očítačová chemie a molekulové modelování I - cvičení -20- Výsledek očítačová chemie a molekulové modelování I - cvičení -21- Gaussian http://www.gaussian.com Nápověda: Tech Support -> Gaussian 09 Help: • Table of Contents • Keyword List !!! Všechny výstupy níže uvedené jsou pouze ukázkami!!! Jejich obsah nesouvisíš řešenými úkoly. Vstupní soubor LinkO příkazy - uvozeny %, lze použít pro specifikaci jména archivu (checkpoint), množství paměti, počtu procesorů určených pro výpočet (všechny položky jsou volitelné) Příkazová řádka - uvozena #, specifikuje typ výpisu (viz Route (#)), metodu a typ výpočtu ►metoda báze vlnové funkce náboj multiplicita M=2S+1 / celkový spin elektronů %chkj^i2o. chk # RHF/cc-pVDZ Opt 7^ optimalizace geometrie optimalizace geometrie molekuly vody H H 0.000000 0.000000 0.000000 -\- 0.000000 0.854766 -0.854766 -0.154167 0.538096 0.538096 prázdný řádek prázdný řádek prázdný řádek specifikace geometrie systému (kartézské souřadnice vÁ) zakončena prázdným řádkem C7800 Počítačová chemie a molekulové modelování I - cvičení Výpočet energie a optimalizace %chk=h2o.chk # RHF/cc-pVDZ Pop=MK <- vypočet energie molekuly vody 0 1 0 H H <— 0.000000 0.000000 0.000000 0.000000 0.854766 -0.854766 -0.154167 0.538096 0.538096 %chk=h2o.chk # RHF/cc-pVDZ Opt Pop=ChelpG optimalizace geometrie molekuly vody 0 1 O H H <- 0.000000 0.000000 0.000000 0.000000 0.854766 ■0.854766 -0.154167 0.538096 0.538096 bez klíčového slova Opt pouze výpočet energie, vlnové funkce a vlastností systému (single point výpočet) zakončeno prázdným řádkem zakončeno prázdným řádkem očítačová chemie a molekulové modelování I - cvičení -24- Průběh výpočtu Výpočet energie Vstupní geometrie O Odhad vlnové funkce (Guess) iO Výpočet E a vlnové funkce (SCF) O Populační analýza a výpočet vlastností systému Populační analýza a výpočet vlastností systému Optimalizace geometrie Q pro vstupní geometrii Vstupní geometrie Q Odhad vlnové funkce (Guess) Výpočet E a vlnové funkce (SCF) Si výpočet gradientu E lokální minimum? ano X 6 opti mahzovana eometrie Populační analýza a výpočet vlastností systému odhad vlnové funkce pro následující krok se bere z kroku předchozího změna geometrie ne očítačová chemie a molekulové modelování I - cvičení -25- [vibrační analýza %chk=h2o. chk # RHF/cc- pVDZ Freq vibracni analýza 0 1 0 0. 0000000 0. 0000000 -0. 0780430 H 0. 0000000 0. 7492420 0. 5000340 H 0. 0000000 -0. 7492420 0. 5000340 končeno prázdným řádkem Pokud je vstupní geometrie stacionárním bodem na PES, lze z počtu imaginárních frekvencí normálních vibrací určit typ stacionárního bodu. 0 imaginárních frekvencí = lokální minimum 1 imaginární (záporná) frekvence = tranzitní stav prvního řádu Vibrační analýza vyžaduje výpočet Hessianu, což může být výpočetně velmi náročné. očítačová chemie a molekulové modelování I - cvičení -26- Průběh výpočtu Výpočet energie Vstupní geometrie O Odhad vlnové funkce (Guess) lO Je použit model ideálního plynu! Harmonická aproximace Harmonická aproximace je platná pouze pro stacionární bod na PES Test na stacionární bod Výpočet E a vlnové funkce (SCF) 8. výpočet gradientu E M očítačová chemie a molekulové modelování I - cvičení -27- CP korekce na BSSE %Chk=ene.chk # RHF/cc-pVDZ Counterpoise=2^- vypocet energie s korekci na BSSE 0 1 0(Fragment=l H(Fragment=l) H(Fragment=l) 0(Fragment=2) H(Fragment=2) H(Fragment=2) <- 231135 374693 074728 648273 .166742 6494 0 0 -0 1 2 063539 253491 864502 893276 355543 1.311764 0.203312 1.116647 0.130817 ■1.229780 ■1.894982 ■0.798971 Výpočet energie s korekcí na BSSE zakončeno prázdným řádkem Fragmenty představují interagující molekuly. Counterpoise korekce na BSSE vyžaduje celkem 5 výpočtů energie[AB, A(B), (A)B, A, B], které program Gaussian provede automaticky. C7800 Počítačová chemie a molekulové modelování I - cvičení Průběh výpočtu Vstupní geometrie Odhad vlnové funkce™ (Guess) Odhad vlnové funkce™ (Guess) Odhad vlnové funkce™ (Guess) I Výpočet E a vlnové funkce (SCF) Populační analýza a výpočet vlastností systému AB Výpočet E a vlnové funkce (SCF) 0 ■ DStl Populační analýza a výpočet vlastností systému 3 Výpočet E a vlnové | funkce (SCF) 1 ™^-™7 í analýza a É vlastností 1 ému 1 Populačn výpočet syst Výpočet E a vlnové funkce (SCF) Populační analýza a výpočet vlastností systému 0 Výpočet E a vlnové funkce (SCF) Populační analýza a výpočet vlastností systému A(B) (A)B B ^cp _ ^ab + (EA EA(B)) + (EB E(A)B) Korigovaná energie BSSE C7800 Počítačová chemie a molekulové modelování I - cvičení Spuštění výpočtu $ module add gaussian $ g0 9 soubor.com aktivace modulu gaussian pouze jednou v daném terminálu \ jméno vstupního souboru program gaussian Po doběhnutí úlohy bude výsledek výpočtu uložen v souboru soubor.log, na posledním řádku souboru musí být uvedeno: Normál termination of Gaussian 09 at Sun v opačném případě je výpočet neúspěšný. Důvod předčasného ukončení výpočtu je nutné hledat ve výstupním souboru. Na klastru WOLF je doporučeno výpočty zadávat do dávkového systému pomocí rozhraní Infinity. Tento postup je podrobně popsán v následující kapitole. očítačová chemie a molekulové modelování I - cvičení -30- Výstupní soubor I a II 1. Vstupní geometrie v interních souřadnicích ! Initial Parameters ! ! (Angstroms and Degrees) ! ! Name Definition Value Derivative Info. ! ! Rl R(l,2) 1.0999 estimate D2E/DX2 ! ! R2 R(l,3) 1.0999 estimate D2E/DX2 ! ! Al A(2,l,3) 101.9929 estimate D2E/DX2 ! 2. Odhad vlnové funkce počet bázových funkci = počet koeficientů c, které je nutné nalézt během SCF procedury Two-electj^i^integral symmetry is turned off. 24 bereis functions, 47 primitive gaussians, 25 cartesian basis functions 5 alpha electrons 5 beta electrons nuclear repulsion energy 8.0071357792 Hartrees. NAtoms= 3 NActive= 3 NUniq= 3 SFac= 1.00D+00 NAtFMM= 60 NAOKFM=F Big=F One-electron integrals computed using PRISM. NBasis= 24 RedAO= T NBF= 24 NBsUse= 24 1.00D-06 NBFU= 24 Harris functional with IExCor= 205 diagonalized for initial guess. počáteční odhad vlnové funkce očítačová chemie a molekulové modelování I - cvičení -31- Výstupní soubor 3. Výpočet energie a vlnové funkce viriálový teorém SCF Done : E (RHF) = ^ Convg = Výsledná energie -75.9901319773 0.5258D-08 A.U. after 10 cycles -V/T = 2.0057 / / Počet SCF kroků potřebných k nalezení E a vlnové funkce Změna energie v posledním kroku SCF procedury Energie je v atomových jednotkách. Úplný průběh SCF procedury je vypsán, pokud je ve vstupním souboru uvedeno #P očítačová chemie a molekulové modelování I - cvičení -32- Výstupní soubor III, pokračování Cycle 1 Pass 1 IDiag 1: E= -75.9710832672194 DUS: error= 4.94D-02 at cycle 1 NSaved= 1. NSaved= 1 IEnMin= 1 EnMin= -75.9710832672194 IErMin= ErrMax= 4.94D-02 EMaxC= 1.00D-01 BMatC= 1.08D-01 BMatP= IDIUse=3 WtCom= 5.06D-01 WtEn= 4.94D-01 : 1 ErrMin= 1 . 08D-01 4.94D-02 Coeff-Com Coeff-En: Coeff: Gap= GapD= RMSDP=6 0 . 100D+01 : 0.100D+01 0 . 100D+01 0.463 Goal= None Shift= 0.000 0.463 DampG=2.000 DampE=0.500 DampFc=l.0000 IDamp= 04D-03 MaxDP=l.13D-01 OVMax= 1.12D-01 kráceno Cycle 10—Pass 1 IDiag 1: E= -7 6^418*!>i64 8 07 68 Delta-E= 0.000000000000 Rises=F Damp=F DIIS: *or= 2ľS^D-0 8 at cycle 10 NSaved= 10. NSaved=10^EnMin=lN^nMin= -76.0418076480768 IErMin=10 ErrMin= 2.15D-08 ErrMax= 2.36D-08 EMaxcV^I.00D-01 BMatC= 8.09D-15 BMatP= 1.54D-13 IDIUse=l WtC Coeff-Com: -Coeff-Com: -Coeff: Coeff: Gap= 0.546 Goal= RMSDP=4.57D-09 MaxDP 1.0 0D+0 0 WS^n ^9D-06 0.154D 1X1D-01 0.669D-01 ■06 0.154D-05 0.18' 0.00D+00 0.181D-04-0.155D-03 0.397D-03 0.649D-03 41D+00 0.138D+01 04-0.155D-03 0.397D-03 0.649D-03 0.138D+01 1 0.669D-01-0.441D+ None Shift= 0 52D-08 DE=-3.84D-13 OVMa" SCF Done E (RHF) Convg -76. 0418076481 0 . 4573D-08 A.U. 92D-0Í 10 cycles 2.0008 C7800 P očítačová chemie a molekulové modelování I - cvičení [výstupní soubor IV 4. Výpočet gradientu sílv _ záporně vzatý gradient energie (v atomových jednotkách) \ Center Atomic Forces (Hartrees/Bohr) Number Numb e r X Y Z 1 8 0. 000000000 0. 000000000 0. 131133317 2 1 0. 000000000 -0. 085534245 -0. 065566658 3 1 0. 000000000 0. 085534245 -0. 065566658 Cartesian Forces: Max 0.131133317 RMS 0.067020838 očítačová chemie a molekulové modelování I - cvičení -34- Výstupní soubor V 5. Optimalizace geometrie GradGradGradGradGradGradGradGradGradGradGradGradGradGradGradGradGradGrad Berny optimization. Internal Forces: Max 0.107734851 RMS 0 . 088033019 kráceno Variable Old X Rl 2.07857 R2 2.07857 Al 1.78011 I tern Maximum Force RMS Force Maximum Displacement RMS Displacement Predicted change in GradGradGradGradGradG ■0 ■0 ■0 -DE/DX 10773 10773 00599 Delta X (Linear) 0 . 00000 0 . 00000 0 . 00000 Delta X (Quad) -0.21160 -0.21160 -0 . 02126 '1'hresholTt 0 .000450 0 .000300 0 .001800 0 .001200 Delta X (Total) -0.21160 -0.21160 -0 . 02126 ged? New X 86697 86697 75885 060008D-02 Tads dGradG GradGradGradGradGráklGrad konvergenční kritéria pro optimalizaci geometrie plánovaná změna geometrie C7800 P očítačová chemie a molekulové modelování I - cvičení Výstupní soubor VI 6. Optimalizovaná geometrie Item Value Threshold Converged? Maximum Force 0.000267 0.000450 YES RMS Force 0.000172 0.000300 YES Maximum Displacement 0.000999 0.001800 YES RMS Displacement 0.000967 0.001200 YES Predicted change in Energy=-2.315429D-07 Optimization completed. -- Stationary point found. ! Optimized Parameters ! ! (Angstroms and Degrees) ! ! Name Definition Value Derivative Info. ! ! Rl R(l,2) 0.9463 -DE/DX = -0.0001 ! ! R2 R(l,3) 0.9463 -DE/DX = -0.0001 ! ! Al A(2,l,3) 104.696 -DE/DX = -0.0003 ! GradGradGradGradGradGradGradGradGradGradGradGradGradGradGradGradGradGrad očítačová chemie a molekulové modelování I - cvičení -36- Výstupní soubor VII 7. Populační analýza a výpočet vlastností systému (kráceno) / obsazené Alpha oce. eigenvalues — -20.54839 -1.34218 -0.70556 -0.56822 -0.49389 Alpha virt. eigenvalues — 0.18742 0.25779 0.79793 0.86396 1.16257 y x energie molekulových orbitalů > neobsazené Mulliken atomic charges: 1 1 O -0.292462 2 h 0.14623i "— Mullikenovy atomové (bodové) náboje 3 H 0.146231 Sum of Mulliken atomic charges = 0.00000 Dipole moment (field-independent basis, Debye): X= 0.0000 Y= 0.0000 Z= Quadrupole moment (field-independent basis, Debye-Ang): XX= -7.0085 YY= -4.1369 ZZ= XY= 0.0000 XZ= 0.0000 YZ= 2.0430 Tot= 2.0430 -5.7327 0 . 0000 dipólový a kvadrupólový elektrostatický moment Typ populační analýzy lze měnit pomocí klíčového slova: pop očítačová chemie a molekulové modelování I - cvičení -37- Výstupní soubor VIII Differentiating once with respect to electric field. with respect to dipole field. Electric field/nuclear overlap derivatives assumed to be zero. Keep Rl ints in memory in canonical form, NReq=873499. There are 3 degrees of freedom in the 1st order CPHF. IDoFFX=0 3 vectors produced by pass 0 Testl2= 3. 17D- 15 3 . 33D- 08 XBigl2= 1 . 15D+00 6. 84D- 01 . 11 . form 3 AO Fock derivatives at one time. 3 vectors produced by pass 1 Testl2= 3. 17D- 15 3 . 33D- 08 XBigl2= 5. 63D- 02 1 . 22D- 01 . 3 vectors produced by pass 2 Testl2= 3. 17D- 15 3 . 33D- 08 XBigl2= 7 . 24D- 03 3. 79D- 02 . 3 vectors produced by pass 3 Testl2= 3. 17D- 15 3 . 33D- 08 XBigl2= 2 . 25D- 04 7 . 06D- 03 . 3 vectors produced by pass 4 Testl2= 3. 17D- 15 3 . 33D- 08 XBigl2= 4 . 04D- 06 8 . 17D- 04 . 3 vectors produced by pass 5 Testl2= 3. 17D- 15 3 . 33D- 08 XBigl2= 2 . 96D- 08 7 . 41D- 05. 3 vectors produced by pass 6 Testl2= 3. 17D- 15 3 . 33D- 08 XBigl2= 3. 37D- 10 9. 10D- 06. 2 vectors produced by pass 7 Testl2= 3. 17D- 15 3 . 33D- 08 XBigl2= 1 . 95D- 12 5. 99D- 07 . 1 vectors produced by pass 8 Testl2= 3. 17D- 15 3 . 33D- 08 XBigl2= 9. 82D- 15 4 . 57D- 08 . Inverted reduced A of dimension 24 with in-core refinement. End of Minotr Frequency-dependent properties file 721 does not exist. End of Minotr Frequency-dependent properties file 722 does not exist. Symmetrizing basis deriv contribution to polar: IMax=3 JMax=2 DiffMx= 0.00D+00 G2DrvN: will do 4 centers at a time, making 1 passes doing MaxLOS=2 Calling FoFCou, ICntrl= 3107 FMM=F IlCent= 0 AccDes= 0.00D+00. FoFDir/FoFCou used for L=0 through L=2. End of G2Drv Frequency-dependent properties file End of G2Drv Frequency-dependent properties file 721 does not exist 722 does not exist očítačová chemie a molekulové modelování I - cvičení -38- Výstupní soubor IX Full mass-weighted force constant matrix: Low frequencies --- -40.7995 -0.0019 -0.0015 Low frequencies --- 1774.9584 4112.7795 4211.8138 6(5) vibrací musí mít nízké frekvence (ideálně nulové) - 3 translační a 3 rotační stupně volnosti systému 0 . 0005 37 . 6815 55.235É Harmonie frequencies cnr ■1), IR intensities (KM/Mole), Raman scattering activities (A**4/AMU), depolarization ratios for plane and unpolarized incident light, reduced masses and normal coordinates: 1 (AMU), force constants (mDyne/A), frekvence normálních vibrací Al B2 Frequencies -- 1774 . 9584 4112 .7795 4211.8137 Red. masses -- 1.0818 1 . 0460 1 . 0821 Frc consts 2.0080 10 . 4246 11.3093 IR Inten ŕ 30.8414 21 .1644 60 .8089 Raman Activ -- 4.7861 68 . 9225 34.7495 Depolar (P) -- 0.5271 0 . 1703 0 .7500 Depolar (U) -- 0.6903 0 .2910 0.8571 Atom AN Y Y X Y Z X Y Z 1 8 0 . 00 0 . 00 0. 07 0 . 00 0 . 00 0. 05 0.00 0.07 0. 00 2 1 0 . 00 -0 . 43 -0. 56 0 . 00 0 .58 -0. 40 0.00 -0.56 0. 43 3 1 0 . 00 0 .43 -0. 56 0 . 00 -0 . 58 -0. 40 0.00 -0.56 -0. 43 směr pohybů atomů pro normální vibraci očítačová chemie a molekulové modelování I - cvičení -39- Výstupní soubor X E (RHF) -76 . 0270533118 Je použit model ideálního plynu! E = E k m (R- J + E opt ,m VRT ,1 - Thermochemistry - Temperature podmínky 298.150 Kelvin. Pressure .^0 000 Atm. Zero-point correction Thermal correction to Thermal correction to Thermal correction to Sum of electronic and Sum of electronic and Sum of electronic and Sum of electronic and Energy= Enthalpy= Gibbs Free Energy= zero-point Energies= thermal Energies= thermal Enthalpies= thermal Free Energies= aceno 0 . 02300 0 . 025843 0.026787 0 . 005411 -76.004045 -76.001211 -76.000267 -76.021642 energie vibrací pro základní stav Gibbsova energie Total Electronic Translational Rotational Vibrational E (Thermal) KCal/Mol 16.216 0 . 000 0 .889 0 .889 14.439 CV Cal/Mol-Kelvin 5. 989 0 . 000 2 . 981 2 . 981 0 . 028 Cal/Mol-Kelvin 44.988 0 .000 34.608 10.376 0 .004 http://gaussian.com/g_whitepap/thermo.htm očítačová chemie a molekulové modelování I - cvičení -40- Výstupní soubor XI 7. Výpočet s korekcí na BSSE Counterpoise: corrected energy = -152.060077655641 Counterpoise: BSSE energy = 0.002829142358 . , energie komplexu s korekcí na BSSE velikost chyby očítačová chemie a molekulové modelování I - cvičení -41- Extrakce dat - optimalizace 1) Aktivace modulu qmutil: $ module add qmutil pouze jednou v daném terminálu 2) Zobrazení průběhu optimalizace (energie): $ extract-gopt-ene soubor.log 3) Průběh optimalizace (všechny geometrie): $ extract-gopt-xyz soubor.log > soubor_opt.xyz 4) Získání optimalizované geometrie (poslední): / $ extract-xyz-str soubor_opt. xyz laŕst > soubor_last. xyz Je vhodné analyzovat průběh optimalizace, např. v programu vmd nebo Avogadro ) Počítačová chemie a molekulové modelování I - cvičení Extrakce dat - optimalizace, pokr. $ extract-gopt-ene soubor.log # Coordinate: # Step Energy [kcal/mol] Energy [au] #-------------------------------------- 1 2 3 4 číslo optimalizačního kroku 0.000 0.171 0.188 0.190 -0.028650961 -0.028922822 -0.028950914 -0.02895 relativní energie vůči absolutní energie v Hartree Energie optimalizované struktury, tedy geometrie obsažené v Geometrii si můžeme prohlédnout v programu v Hartree. 7800 Počítačová chemie a molekulové modelování I - cvičení Infinity https://lcc.ncbr.muni.cz/whitezone/development/infinity/ očítačová chemie a molekulové modelování I - cvičení -44- Přehled příkazů Správa software: * module aktivace/deaktivace software Správa úloh: • psubmit • pinfo • pgo • pjobs zadání úlohy do dávkového systému informace o úloze přihlásí uživatele na výpočetní uzel, kde se úloha vykonává přehled úloh uživatele zadaných do dávkového systému Další užitečné příkazy: pqueues, pnodes, pqstat očítačová chemie a molekulové modelování I - cvičení -45- Životní cyklus úlohy výpočetní uzel#l výpočetní uzel#2 np=l torque-mom torque-mom (pbs_mom) (pbs_mom) stav uzlů (nodes) 1 řazení úloh na výpočetní uzly dle požadovaných zdrojů torque-scheduler (pbs_sched) připravená úloha čekající úloha I V V X F F I I bezici úloha 7800 Počítačová chemie a molekulové modelování I - cvičení Úloha Úloha musí splňovat následující podmínky: • každá úloha se spouští v samostatném adresáři • všechny vstupní data úlohy musí být v adresáři úlohy • adresáře úloh nesmí být do sebe zanořené • průběh úlohy je řízen skriptem nebo vstupním souborem (u automaticky detekovaných úloh) • skript úlohy musí být v bashi • ve skriptu úlohy se nesmí používat absolutní cesty, všechny cesty musí být uvedeny relativně k adresáři úlohy Spuštění úlohy v gaussianu Úlohy v programu Gaussian budeme spouštět vždy přes rozhraní Infinity. Úlohy budeme spouštět na 1 CPU ve frontě normál, kde výpočetní uzel bude vaše pracovní stanice. aktivace modulu gaussian $ module add gaussian *--- pouze jednou v daném terminálu $ psubmit localhost@normal soubor.com ..... „ . _ , vstupní soubor pro program gaussian úloha bude zaražena do fronty normál (max imální doba běhu úlohy 12 hodin) běh úlohy je omezen pouze na vaší pracovní stanici Po doběhnutí úlohy bude výsledek výpočtu uložen v souboru soubor.log K monitorování průběhu úlohy lze použít příkaz pinfo, který se spouští v adresáři úlohy. Dalšími možnostmi jsou příkazy pjobs a pqstat. Pokud je úloha spuštěna, je možné použít příkaz pgo, který vás přihlásí na výpočetní uzel a změní aktuální adresář do pracovního adresáře úlohy. Lze tak kontrolovat průběh výpočtu za běhu úlohy. Servisní soubory V adresáři úlohy vznikají při zadání úlohy do dávkového systému a dále v průběhu života úlohy a po jejím ukončení servisní soubory Jedná se o textové soubory, které je možné analyzovat pomocí libovolného textového prohlížeče (příkazy less, more) či editoru (vi, kwrite, gedit). Pro účely cvičení je možné obsah souborů ignorovat, kromě případů, kdy úloha nedoběhne v pořádku. Význam souborů je následující: kontrolní soubor s informacemi o průběhu úlohy vlastní skript (wrapper), který se spouští dávkovým systémem standardní výstup z běhu *.infex skriptu, nutno analyzovat při nestandardním ukončení úlohy seznam uzlů vyhrazených pro úlohu seznam GPU karet vyhrazených pro úlohu unikátní identifikátor úlohy standardní výstup z běhu skriptu úlohy, nutno analyzovat při nestandardním ukončení úlohy *.info *.infex *.infout *.nodes *.gpus *.key *.stdout očítačová chemie a molekulové modelování I - cvičení -49- Přenos dat Typ synchronizace Význam syne Data úlohy jsou při spuštění úlohy zkopírována ze vstupního adresáře do pracovního adresáře. Pak je spuštěn vlastní výpočet. Po jeho doběhnutí jsou všechna data z pracovního adresáře zkopírována zpět do vstupního adresáře úlohy. Uživatelské rozhraní (Ul) (Frontend) /job/input/di rsync Výpočetní uzel #1 orker Node (WN) /scratch/job_id/ rsync Poznámky: • výchozí synchronizační mód • při použití localhost@normal je výpočetní uzel a uživatelské rozhraní identický stroj C7800 Počítačová chemie a molekulové modelování I - cvičení Nemesis očítačová chemie a molekulové modelování I - cvičení -51- Vizualizace optimalizace geometric 1) Projekt: Trajectory 2) File->lmport Trajectory as -> Gaussian Geometry Optimization Project 1: NEMESIS - Molecular Modelling Package ■ Trajectories ) IS Name Trajectory Al ► Trajectory 1 Structure 1 d voj k CT Number of trajectories: 1 í*^ Trajectories Structures lActive ProFile I Active profile: Profile 1 Profile objects Name TV " Light 1 Light | Background 1 Backgrc * Standard Model 1 Standar Freeze d Atoms 1 Freezec er ik ram U S T C C D P nl |i? í? El [š len m @ í t|, H>) 9:01 PM X Petr Kulhánek i|} -v * Trajectory Basic Play Segments Filters Referenced by SI Name Snapshosts 1 ethan 6 Gaussian Geometry ( Basic Energy info dvoj kli k Properties Manager Gaussian Geometry c\titaization ID Energy [a.u.] Relative Energy [kcc 1 0.028617175 0.00 2 -0.028893932 -0.17 3 -0.028929498 -0.20 4 -0.028931 S79 -0.20 5 -0.028934299 -0.20 ö -0.028935040 -0.20 ; ■■■ lvi no * 7 i !try m m tion □istance Angle Torsion -q (M) [4j [►] II a - M Restrain Property Label průběh optimalizace 7800 Počítačová chemie a molekulové modelování I - cvičení 1) Projekt: Trajectory 2) File->lmport Trajectory as -> Gaussian Vibrations Project 1: NEMESIS - Molecular Modelling Package Trajectories rBD t> Name Trajectory At ► Trajectory 1 Structure 1 Number of trajectories: 1 5*^ Trajectories | Structures Active Profile OH IS Active profile: Profile 1 Profile objects Name TV * Light 1 Light I Background 1 Backgrc * Standard Model 1 Standar Freezed Atoms 1 Freezed cm IE ZTXE) ILüLUniunii H1 i tj, H))J 10:08 PM X Petr Kulhánek Preperties Manager ID S) dvoj kli k UČ « Trajectory Basic Play Segments Filters Referenced by SI Name Snapshasts Type 1 ethanfreq 180 Gaussian Vibratk \ dvojl| ik Name Type Description ID Gaussian Vibrations Basic Vibrations setup ^volíme vibraci Frequency IR Intensity Sec □ □ □ □ □ 224.Ů 878.2 878.2 11 20.0 11 20.0 11 37.8 1359.3 1408.2 1408.2 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 ) 50 Nummer of vibrations: 24 ActiveVibrations: 0 , Active imaginary | Deactivate all tance Angle Torsion spustíme animaci nu'i |u s t c c d p &\ \x% ,1? q mi r Avogadro očítačová chemie a molekulové modelování I - cvičení -54- Vizualizace vibrací Do programu Avogadro načteme soubor.log, obsahující výsledky vibrační analýzy. Souhrn frekvencí jednotlivých normálních vibrací najdeme v menu Extensions->Vibrations. Avogadro Navigate Settings S Display visual cues Toolsettings... j| Display Settings.. View 1 Molecular Vibrations Frequency (cm-1) Intensity (km/mol) 1 223.7 0.0 2 877.8 1.5 3 878.0 1.5 4 1,119.8 0.0 5 ^L119.9 0.0 6 1,13*^ 0.0 7 1,359.2 0.3 8 1,408.1 0.3 Options Scale: 1 i Display force vectors □ Animation speed set by Frequency Export... Close vizualizace vibrací Messages les ~ m $ *l 4D 12:22PM X frekvence vibrací C7800 P očítačová chemie a molekulové modelování I - cvičení očítačová chemie a molekulové modelování I - cvičení -56- Vizualizace optimalizace geometric Do programu načteme xyz trajektorii (průběh optimalizace) vyextrahovanou skriptem extract-gopt-xyz z modulu qmutil. počet optimalizačních kroků L me<ž VMD Main ™ k/ Extension File Molecule Graphics Display Mouse Help ID T A D F Molecule Atoms y Frames Vol 0 TAD opt.xyz 3 4 0 l-il« om □ Loop step J 1 ]► speed É přesouvání mezi jednotlivými geometriemi očítačová chemie a molekulové modelování I - cvičení -57- [ Vizualizace volumetrických dat Volumetrické data (cube soubory vytvořené programem cubegen) lze načíst přímo do programu vmd. Ve výchozí vizualizaci se zobrazí molekula v čárovém modelu bez volumetrických dat. Volumetrická data lze zobrazit jako (Representations): • isoplochu (isosurface) • řez (volumeslice) Mapování elektrostatického potenciálu na isoplochu elektronové hustoty: 1) načteme hustotu a elektrostatický potenciál do jedné molekuly 2) zobrazíme isoplochu elektronové hustoty 3) pro isoplochu zvolíme vybarvení (Coloring Method) podle volumetrických dat (Volume) - zvolíme elektrostatický potenciál, škála barev se nastavuje na záložce Trajectory (Color Scale Data Range) očítačová chemie a molekulové modelování I - cvičení -58- Extrapolace na CBS očítačová chemie a molekulové modelování I - cvičení -59- Extrapolace E(x) = ECBS + Ae Protože máme čtyři vstupní energie (pro x=2, 3, 4, a 5) a jen tři neznámé (ECBS, A a B) musíme použít metodu nejmenších čtverců. Cílem metody je nalézt hodnotu parametrů Ecbs' A a B tak, aby účelová (chybová) funkce byla minimální. 5 fCE^ ,A,B) = ^[E(x,Ecbs ,A,B)-EHF(x)f = min! x=2 K hledání optimálních parametrů můžeme použít metodu fit z programu gnuplot. Viz originální dokumentace gnuplotu nebo: http://www.root.cz/clanky/gnuplot-prikaz-fit/ -60- očítačová chemie a molekulové modelování I - cvičení Postup ■ Připravíme textový soubor data.txt, který bude obsahovat dva sloupce: kardinální číslo báze (2, 3, 4,...) a vypočtenou energii metodou HF. ■ Spustíme program gnuplot a zobrazíme průběh energie ze souboru data.txt: gnuplot> plot './data.txt' using 1:2 with points ■ Definujeme funkci pro extrapolaci: gnuplot> E(x) = Ecbs + A*exp(-B*x) ■ Nastavíme výchozí hodnoty parametru pro optimalizaci: gnuplot> A=l gnuplot> B=l gnuplot> Ecbs=-80 # nizsi nez nejmensi vypočtena energie ■ Provedeme před-optimalizaci parametrů (ECBS a A a pak finální optimalizaci všech parametrů: gnuplot> fit E(x) "./data.txt" via Ecbs, A gnuplot> fit E(x) "./data.txt" via Ecbs, A, B očítačová chemie a molekulové modelování I - cvičení Postup, pokračování Zobrazíme vstupní data, funkci E(x) a hodnotu ECBS. Provedeme vizuální kontrolu získaných výsledů. Funkce E(x) musí procházet všemi body a limitně se blížit nalezené hodnotě E -CBS- gnuplot> set xrange[2:7] gnuplot> plot './data.txť using 1:2 with points, E(x), Ecbs Vypíšeme přesnou hodnotu E CBS gnuplot> print Ecbs Water Molecule - HF/CBS 3 .21 LU -76.020 i--76.025 --76.030 --76.035 --76.040 --76.045 --76.050 --76.055 --76.060 --76.065 --76.070 : -76.075 L T E(x) ECBS 1 3 4 5 Cardinal number x 6 očítačová chemie a molekulové modelování I - cvičení -62-