Zagadnienie początkowe dla równań różniczkowych zwyczajnych pierwszego rzędu

1. Zagadnienie początkowe dla równań różniczkowych zwyczajnych pierwszego rzędu Dany jest układ różniczkowych Należy znaleźć jego rozwiązanie x(t) dla warunków początkowych zadanych przez x(t0).

2. Przykłady problemów które można zapisać bezpośrednio jako układu równań różniczkowych rzędu pierwszego: • Kinetyka chemiczna (x jest stężeniem lub w ogólności postępem reakcji a t czasem). • Obliczanie trajektorii „wewnętrznej” współrzędnej reakcji (Intrinsic Reaction Coordinate; IRC); wtedy x jest wektorem współrzędnych układu reagującego a t współrzędną reakcji (z). Obliczenia startuje się z punktu siodłowego w dwóch kierunkach określonych przez wektor własny hesjanu odpowiadający ujemnej wartości własnej.

3. E x1 x2 H H3N Cl x2 NH4+...Cl- x1 x1* NH3...H...Cl E NH3...HCl x2 x2* x1

4. Dygresja: jeżeli prawa strona jest układem funkcji liniowych względem x (tak jak w kinetyce reakcji pierwszego rzędu) to rozwiązanie jest analityczne i ma postać kombinacji liniowej funkcji eksponencjalnych:

5. Często nie jest aż tak prosto (przykład: reakcje Biełousowa-Żabotyńskiego).

6. Numeryczne rozwiązywanie zagadnienia początkowego: ogólnie. Dzielimy przedział [t0,b] w którym szukamy rozwiązania na N odcinków t0<t1<t2<...<tN=b Definiujemy hi=ti+1-ti oraz przybliżone rozwiązanie i różnicę pomiędzy przybliżonym a dokładnym rozwiązaniem w punkcie ti X(ti)=Xi»x(ti)=yi ei=||xi-X||i Mamy:

7. Ogólny podział metod numerycznych rozwiązywania równań różniczkowych zwyczajnych pierwszego rzędu: • Metody jednokrokowe: do obliczenia x w kolejnym kroku t wykorzystują tylko wartości x z poprzedniego kroku (np. metoda Eulera-Cauchy’ego,metoda Rungego-Kutty). • Metody wielokrokowe: wykorzystują wartości x z kilku kroków wstecz (np. metoda Adamsa-Bashforda, metoda Adamsa-Stromera). • Metody ekstrapolacyjne: wykorzystują kwadraturę Romberga. • Metody predyktora i korektora: w danym kroku t najpierw oblicza się przewidywane wartości x (krok predykcyjny) a potem się je udokładnia (krok korekcyjny) (np. algorytm Geara). Tych metod można używać w połączeniu zarówno z metodami jedno- jak i wielokrokowymi.

8. Metoda Eulera-Cauchy’ego Błąd popełniany w i-tym kroku całkowania (lokalny) Całkowity błąd którym obarczone jest rozwiązanie w i-tym punkcie (globalny) rozwiązanie przybliżone rozwiązanie dokładne f x t1 t2 t1 t2 t5 t5 t3 t3 t3 t4 t3 t4 t t

9. Udoskonalona metoda Eulera-Cauchy’ego (drugiego rzędu)

10. Metoda predyktora-korektora Heuna

11. Metody Rungego-Kutty Sformułowanie ogólne: Rząd zbieżności metod Rungego-Kutty w zależności od m m=1: metoda Eulera-Cauchy’ego m=2: ulepszona metoda Eulera-Cauchy’ego lub metoda Heuna m=4: klasyczna metoda Rungego-Kutty.

12. Klasyczna metoda Rungego-Kutty

13. Porównanie metody Eulera-Cauchy’ego, zmodyfikowanej metody Eulera-Cauchy’ego oraz metody Rungego-Kutty dla zagadnienia z h=0.5 y x

14. Metody implicite Rungego-Kutty W najprostszej wersji (m=2 wychodzimy z ulepszonej metody Eulera-Cauchy’ego i doprowadzamy do samouzgodnienia wartości f na końcu przedziału: Ogólnie, dla rzędu m prowadzimy kwadraturę Gaussa-Legendre’a (z węzłami w miejscach zerowych wielomianu Legendre’a rzędu m). Globalny błąd rozwiązania jest wtedy rzędu O(h2m) (w porównaniu z O(hm) dla formuł explicite Rungego-Kutty.

15. Dla m=2.

16. Metody wielokrokowe W metodzie s-krokowej do całkowania wykorzystujemy wartości X i f obliczone w s poprzednich krokach. Przez te punkty prowadzimy wielomian interpolacyjny Fs(t)=Fs(t,X(t)), który następnie całkujemy. Przykład: metoda Adamsa-Bashforda z wykorzystaniem trzech kroków wstecz.

17. Metoda Geara dla układów sztywnych

18. Zagadnienie brzegowe dla równań drugiego rzędu. Całkowanie równań dynamiki molekularnej

19. Algorytm Verleta:

20. Prędkościowy algorytm Verleta (velocity Verlet) Krok 1: Krok 2:

21. Algorytm “zabiego skoku” (leapfrog): Wszystkie trzy algorytmy są algorytmami symplektycznymi, tj, całkowita energia układu oscyluje wokół pewnej stałej wartości bliskiej początkowej energii całkowitej (inaczej: zachowują “cień hamiltonianu” (shadow Hamiltonian). Takiej właściwości nie mają wszystkie algorytmy dynamiki molekularnej (np. algorytm Geara). Algorytmy symplektyczne zaprojektowano również do symulacji MD w warunkach izokinetycznych (stała temperatura) oraz izotermiczno-izobarycznych (stała temperatura i ciśnienie).

22. Zależność składowych energii i energii całkowitej od czasu dla symulacji MD Ac-Ala10-NHMe (Khalili et al., J. Phys. Chem. B, 2005, 109, 13785-13797) Energia kinetyczna Energia [kcal/mol] Energia całkowita Energia potencjalna Energia całkowita 0.0 1.0 2.0 3.0 4.0 5.0 Czas [ns]

23. Literatura dotycząca algorytmów całkowania równań dynamiki molekularnej: • Frenkel, D.; Smit, B. Understanding molecular simulations, Academic Press, 1996, rozdział 4. • Calvo, M. P.; Sanz-Serna, J. M. Numerical Hamiltonian Problems; Chapman & Hall: London, U. K., 1994. • Verlet, L. Phys. Rev. 1967, 159, 98. • Swope, W. C.; Andersen, H. C.; Berens, P. H.; Wilson, K. R. J. Chem. Phys. 1982, 76, 637. • Tuckerman, M.; Berne, B. J.; Martyna, G. J. J. Chem. Phys. 1992, 97, 1990. • Ciccotti, G.; Kalibaeva, G. Philos. Trans. R. Soc. London, Ser. A 2004, 362, 1583.