Katedra Chemii Fizycznej SFU (RSU)
METODY NUMERYCZNE I PROGRAMOWANIE
Materiały do ​​zajęć wykładowych
Wykładowca – art. Obrót silnika. Szczerbakow I.N.

ROZWIĄZANIE RÓWNAŃ RÓŻNICZKOWYCH ZWYKŁYCH

Opis problemu

Rozwiązując problemy naukowe i inżynieryjne, często konieczne jest matematyczne opisanie jakiegoś układu dynamicznego. Najlepiej zrobić to w formie równań różniczkowych ( DU) lub układy równań różniczkowych. Najczęściej problem ten pojawia się przy rozwiązywaniu problemów związanych z modelowaniem kinetyki reakcji chemicznych i różnych zjawisk przenoszenia (ciepła, masy, pędu) – przenoszenia ciepła, mieszania, suszenia, adsorpcji, przy opisie ruchu makro- i mikrocząstek.

Równanie różniczkowe zwyczajne(ODE) n-tego rzędu jest równaniem, które zawiera jedną lub więcej pochodnych żądanej funkcji y(x):

Tutaj y(n) oznacza pochodną rzędu n jakiejś funkcji y(x), x jest zmienną niezależną.

W niektórych przypadkach równanie różniczkowe można przekształcić do postaci, w której jawnie wyrażona jest najwyższa pochodna. Ta forma zapisu nazywa się równaniem, rozstrzygane w odniesieniu do najwyższej pochodnej(w tym przypadku po prawej stronie równania nie ma najwyższej pochodnej):

To właśnie ta forma nagrywania jest akceptowana standard przy rozważaniu metod numerycznych rozwiązywania ODE.

Liniowe równanie różniczkowe jest równaniem liniowym względem funkcji y(x) i wszystkich jej pochodnych.

Na przykład poniżej znajdują się liniowe ODE pierwszego i drugiego rzędu

Rozwiązywanie równania różniczkowego zwyczajnego jest funkcją y(x), która dla dowolnego x spełnia to równanie w pewnym skończonym lub nieskończonym przedziale. Nazywa się proces rozwiązywania równania różniczkowego poprzez całkowanie równania różniczkowego.

Ogólne rozwiązanie ODE N-ty rząd zawiera n dowolnych stałych C 1 , C 2 , …, C n

Wynika to oczywiście z faktu, że całka nieoznaczona jest równa funkcji pierwotnej całki plus stała całkowania

Ponieważ do rozwiązania równań różniczkowych n-tego rzędu potrzeba n całkowań, w rozwiązaniu ogólnym pojawia się n stałych całkowania.

Prywatne rozwiązanie ODE uzyskuje się z ogólnego, jeśli stałym całkowania nadano określone wartości, definiując dodatkowe warunki, których liczba pozwala nam obliczyć wszystkie niepewne stałe całkowania.

Dokładne (analityczne) rozwiązanie (ogólnego lub szczegółowego) równania różniczkowego implikuje uzyskanie pożądanego rozwiązania (funkcji y(x)) w postaci wyrażenia z funkcji elementarnych. Nie zawsze jest to możliwe nawet w przypadku równań pierwszego rzędu.

Rozwiązanie numeryczne DE (iloraz) polega na obliczeniu funkcji y(x) i jej pochodnych w pewnych punktach leżących na pewnym odcinku. Oznacza to, że w rzeczywistości rozwiązanie równania różniczkowego n-tego rzędu postaci uzyskuje się w postaci poniższej tabeli liczb (kolumnę wartości najwyższej pochodnej oblicza się przez podstawienie wartości do równanie):

Na przykład w przypadku równania różniczkowego pierwszego rzędu tabela rozwiązań będzie zawierać dwie kolumny - x i y.

Nazywa się zbiór wartości odciętych, w którym określana jest wartość funkcji oczko, na którym zdefiniowana jest funkcja y(x). Same współrzędne są wywoływane węzły siatki. Najczęściej dla wygody się je stosuje jednolite siatki, w którym różnica między sąsiednimi węzłami jest stała i nazywa się odstęp siatki Lub etap integracji równanie różniczkowe

Lub , I= 1, …, N

Aby określić rozwiązanie prywatne konieczne jest ustawienie dodatkowych warunków, które pozwolą na obliczenie stałych całkowania. Co więcej, powinno być dokładnie n takich warunków. Dla równań pierwszego rzędu - jeden, dla drugiego - 2 itd. W zależności od sposobu ich określenia przy rozwiązywaniu równań różniczkowych wyróżnia się trzy rodzaje problemów:

· Problem Cauchy'ego (problem początkowy): Musze znaleźć coś takiego rozwiązanie prywatne równanie różniczkowe spełniające pewne warunki początkowe określone w jednym punkcie:

czyli podana jest pewna wartość zmiennej niezależnej (x 0) oraz wartość funkcji i wszystkich jej pochodnych aż do rzędu (n-1) w tym punkcie. Ten punkt (x 0) nazywa się podstawowy. Na przykład, jeśli rozwiązywane jest DE pierwszego rzędu, wówczas warunki początkowe wyrażane są jako para liczb (x 0 , y 0)

Tego rodzaju problem pojawia się podczas rozwiązywania ODA, które opisują na przykład kinetykę reakcji chemicznych. W tym przypadku znane są stężenia substancji w początkowym momencie ( t = 0) i konieczne jest znalezienie stężeń substancji po pewnym czasie ( T) . Jako przykład można również przytoczyć problem wymiany ciepła lub przenoszenia masy (dyfuzji), równania ruchu punktu materialnego pod wpływem sił itp.

· Problem wartości brzegowych . W takim przypadku wartości funkcji i (lub) jej pochodnych są znane w więcej niż jednym punkcie, na przykład w początkowym i końcowym momencie czasu, i konieczne jest znalezienie konkretnego rozwiązania równania różniczkowego pomiędzy tymi punktami. W tym przypadku nazywane są same dodatkowe warunki regionalny (linia graniczna) warunki. Oczywiście problem wartości brzegowych można rozwiązać dla ODE co najmniej drugiego rzędu. Poniżej znajduje się przykład ODE drugiego rzędu z warunkami brzegowymi (podano wartości funkcji w dwóch różnych punktach):

· Problem Sturma-Liouville'a (problem wartości własnych). Problemy tego typu są podobne do problemów wartości brzegowych. Rozwiązując je, należy znaleźć rozwiązanie, przy jakich wartościach dowolnego parametru DU spełnia warunki brzegowe (wartości własne) i funkcje będące rozwiązaniem DE dla każdej wartości parametru (funkcje własne). Na przykład wiele problemów w mechanice kwantowej to problemy z wartościami własnymi.

Numeryczne metody rozwiązywania problemu Cauchy'ego ODE pierwszego rzędu

Rozważmy kilka metod numerycznych rozwiązywania Problemy Cauchy’ego(zadanie początkowe) równania różniczkowe zwyczajne pierwszego rzędu.

(6.2)

Zapiszmy to równanie w postaci ogólnej, rozwiązanej względem pochodnej (prawa strona równania nie zależy od pierwszej pochodnej):

Znalezienie wartości funkcji y w danych punktach siatki jest konieczne, jeżeli znane są wartości początkowe, gdzie w punkcie początkowym x 0 znajduje się wartość funkcji y(x).

Przekształćmy równanie, mnożąc przez d x

(6.3)

I integrujemy lewą i prawą stronę pomiędzy i-tym i i+ pierwszymi węzłami siatki.

Spośród wielu metod opracowanych do rozwiązywania ODE pierwszego rzędu, rozważamy metody i . Są dość proste i dają wstępne wyobrażenie o podejściu do rozwiązania tego problemu w ramach rozwiązania numerycznego.

Metoda Eulera

Historycznie rzecz biorąc, pierwszym i najprostszym sposobem numerycznego rozwiązania problemu Cauchy'ego dla ODE pierwszego rzędu jest metoda Eulera. Opiera się na przybliżeniu pochodnej przez stosunek skończonych przyrostów zależnej ( y) i niezależne ( X) zmienne pomiędzy węzłami siatki jednolitej:

gdzie y i+1 jest pożądaną wartością funkcji w punkcie x i+1.

Jeśli teraz przekształcimy to równanie i uwzględnimy jednorodność siatki całkowej, otrzymamy iteracyjny wzór, dzięki któremu możemy obliczyć tak, ja+1, jeśli y i jest znane w punkcie x i:

Porównując wzór Eulera z otrzymanym wcześniej wyrażeniem ogólnym, widać, że do przybliżonego obliczenia całki w metodzie Eulera stosuje się najprostszy wzór na całkowanie - wzór na prostokąty wzdłuż lewej krawędzi odcinka.

Graficzna interpretacja metody Eulera jest również łatwa (patrz rysunek poniżej). Rzeczywiście, z postaci rozwiązywanego równania () wynika, że ​​wartość jest wartością pochodnej funkcji y(x) w punkcie x=x i - , a zatem jest równa tangensowi kąt styczny poprowadzony do wykresu funkcji y(x) w punkcie x =x i .

Z prawego trójkąta na rysunku możesz znaleźć

Stąd pochodzi wzór Eulera. Zatem istotą metody Eulera jest zastąpienie funkcji y(x) na odcinku całkującym linią prostą styczną do wykresu w punkcie x=x i. Jeśli pożądana funkcja różni się znacznie od liniowej na segmencie całkowania, wówczas błąd obliczeniowy będzie znaczny. Błąd metody Eulera jest wprost proporcjonalny do kroku całkowania:

Błąd~o

Struktura procesu obliczeniowego jest następująca. Dla danych warunków początkowych x 0 I y 0 można obliczyć

W ten sposób tworzona jest tabela wartości funkcji y(x) z pewnym krokiem ( H) Przez X na segmencie. Błąd w określeniu wartości y(x ja) H w tym przypadku im mniejsza wybrana długość kroku, tym będzie ona mniejsza

(o czym decyduje dokładność wzoru na całkowanie).

W przypadku dużego h metoda Eulera jest bardzo niedokładna. Daje coraz dokładniejsze przybliżenie w miarę zmniejszania się kroku całkowania. Jeżeli segment jest zbyt duży, wówczas każdy odcinek dzieli się na N segmentów całkowych i do każdego z nich stosuje się wzór Eulera z krokiem, czyli krok całkowania h jest przyjmowany mniejszy niż krok siatki, na której znajduje się rozwiązanie jest zdeterminowany.

Korzystając z metody Eulera, skonstruuj przybliżone rozwiązanie następującego problemu Cauchy'ego:

Na siatce z krokiem 0,1 w przedziale (6,5)

Rozwiązanie:

Równanie to zostało już zapisane w postaci standardowej, rozwiązanej względem pochodnej pożądanej funkcji.

Dlatego dla rozwiązywanego równania mamy

Przyjmijmy krok całkowania równy krokowi siatki h = 0,1. W takim przypadku dla każdego węzła sieci zostanie obliczona tylko jedna wartość (N=1). Dla pierwszych czterech węzłów siatki obliczenia będą następujące:

Pełne wyniki (z dokładnością do piątego miejsca po przecinku) podano w trzeciej kolumnie – h =0,1 (N =1). Dla porównania druga kolumna tabeli pokazuje wartości obliczone z analitycznego rozwiązania tego równania .

Druga część tabeli pokazuje błąd względny otrzymanych rozwiązań. Można zauważyć, że przy h = 0,1 błąd jest bardzo duży, sięgający 100% dla pierwszego węzła x = 0,1.

Tabela 1 Rozwiązanie równania metodą Eulera (dla słupów wskazano stopień całkowania oraz liczbę odcinków całkowania N pomiędzy węzłami siatki)

Zmniejszmy krok całkowania o połowę, h = 0,05, w tym przypadku dla każdego węzła sieci obliczenia zostaną przeprowadzone w dwóch krokach (N = 2). Zatem dla pierwszego węzła x =0,1 otrzymujemy:

(6.6)

Wzór ten okazuje się ukryty w odniesieniu do y i+1 (wartość ta znajduje się zarówno po lewej, jak i prawej stronie wyrażenia), czyli jest równaniem względem y i+1, które można rozwiązać, np. numerycznie, stosując jakąś metodę iteracyjną (w takiej postaci można to uznać za iteracyjny wzór prostej metody iteracyjnej). Można to jednak zrobić inaczej i około obliczyć wartość funkcji w węźle ja+1 używając zwykłego wzoru:

,

które można następnie wykorzystać w obliczeniach zgodnie z (6.6).

To daje metodę Guna lub metoda Eulera z przeliczeniem. Dla każdego węzła integracji wykonywany jest następujący ciąg obliczeń

(6.7)

Dzięki dokładniejszemu wzorowi całkowania błąd metody Hün jest proporcjonalny do kwadratu kroku całkowania.

Błąd~ godz. 2

Podejście zastosowane w metodzie Guna służy do konstruowania tzw. metod prognoza i korekta, co zostanie omówione później.

W przypadku dużego h metoda Eulera jest bardzo niedokładna. Daje coraz dokładniejsze przybliżenie w miarę zmniejszania się kroku całkowania. Jeżeli segment jest zbyt duży, wówczas każdy odcinek dzieli się na N segmentów całkowych i do każdego z nich stosuje się wzór Eulera z krokiem, czyli krok całkowania h jest przyjmowany mniejszy niż krok siatki, na której znajduje się rozwiązanie jest zdeterminowany.

Przeprowadźmy obliczenia równania () metodą Huna.

Przy kroku całkowania h = 0,1 w pierwszym węźle siatki x 1 otrzymujemy:

Co jest znacznie dokładniejsze niż wartości uzyskane metodą Eulera przy tym samym kroku całkowania. W tabeli 2 poniżej przedstawiono porównawcze wyniki obliczeń dla h = 0,1 metod Eulera i Gün.

Tabela 2 Rozwiązanie równania metodami Eulera i Guna

Zauważmy znaczny wzrost dokładności obliczeń metody Hün w porównaniu z metodą Eulera. Zatem dla węzła x = 0,1 odchylenie względne wartości funkcji wyznaczonej metodą Huyna okazuje się 30 (!) razy mniejsze. Tę samą dokładność obliczeń przy użyciu wzoru Eulera uzyskuje się, gdy liczba segmentów całkujących N wynosi około 30. W związku z tym przy zastosowaniu metody Hün z tą samą dokładnością obliczeń zajmie to około 15 razy mniej czasu komputerowego niż przy zastosowaniu metody Eulera .

Sprawdzenie stabilności rozwiązania

Rozwiązanie ODE w pewnym punkcie x i nazywa się stabilnym, jeżeli wartość funkcji znaleziona w tym punkcie tak, ja zmienia się nieznacznie wraz ze zmniejszaniem się stopnia integracji. Aby sprawdzić stabilność, należy zatem przeprowadzić dwa obliczenia wartości ( tak, ja) – z krokiem całkowania h i ze zmniejszoną (np. o dwa) wielkością kroku

Jako kryterium stabilności można przyjąć małość względnej zmiany otrzymanego rozwiązania przy zmniejszeniu kroku całkowania (ε jest z góry określoną małą wartością)

Kontrolę tę można przeprowadzić dla wszystkich rozwiązań w całym zakresie wartości X. Jeśli warunek nie jest spełniony, krok jest ponownie dzielony na pół i znajdowane jest nowe rozwiązanie itp. aż do uzyskania stabilnego roztworu.

Metody Runge-Kutty

Dalsza poprawa dokładności rozwiązywania ODE pierwszego rzędu jest możliwa poprzez zwiększenie dokładności przybliżonego obliczenia całki w wyrażeniu.

Widzieliśmy już zaletę przejścia od całkowania przy użyciu wzoru na prostokąt () do stosowania wzoru na trapez () podczas przybliżania tej całki.

Korzystając ze sprawdzonego wzoru Simpsona, można uzyskać jeszcze dokładniejszy wzór na rozwiązanie problemu Cauchy'ego dla ODE pierwszego rzędu - szeroko stosowanej w praktyce obliczeniowej metody Runge-Kutty.

Zaletą wieloetapowych metod Adamsa rozwiązywania ODE jest to, że w każdym węźle obliczana jest tylko jedna wartość prawej strony ODE - funkcja F(x,y). Wady obejmują niemożność rozpoczęcia metody wieloetapowej od jednego punktu początkowego, ponieważ obliczenia z wykorzystaniem wzoru na k-etap wymagają znajomości wartości funkcji w k węzłach. Dlatego konieczne jest uzyskanie rozwiązania (k-1) w pierwszych węzłach x 1 , x 2 , …, x k-1 jakąś metodą jednoetapową, np. metodą

Rozważamy jedynie rozwiązanie problemu Cauchy’ego. Układ równań różniczkowych lub jedno równanie należy przekształcić do postaci

Gdzie ,
–N-wektory wymiarowe; y– nieznana funkcja wektorowa; X– niezależny argument,
. W szczególności, jeśli N= 1, wówczas układ zamienia się w jedno równanie różniczkowe. Warunki początkowe ustala się następująco:
, Gdzie
.

Jeśli
w pobliżu punktu
jest ciągły i ma ciągłe pochodne cząstkowe względem y, to twierdzenie o istnieniu i niepowtarzalności gwarantuje, że istnieje tylko jedna ciągła funkcja wektorowa
, zdefiniowany w Niektóre sąsiedztwo punktu , spełniający równanie (7) i warunek
.

Zwróćmy uwagę na fakt, że sąsiedztwo punktu , w którym określono rozwiązanie, może być bardzo mała. Zbliżając się do granicy tego sąsiedztwa, rozwiązanie może sięgać w nieskończoność, oscylować z nieskończenie rosnącą częstotliwością, generalnie zachowuje się tak źle, że nie można go kontynuować poza granicę sąsiedztwa. W związku z tym takiego rozwiązania nie można prześledzić metodami numerycznymi na większym segmencie, jeśli jest to określone w opisie problemu.

Rozwiązywanie problemu Cauchy'ego na [ A; B] jest funkcją. W metodach numerycznych funkcję zastępuje się tabelą (tabela 1).

Tabela 1

Tutaj
,
. Odległość między sąsiednimi węzłami tabeli przyjmuje się zwykle jako stałą:
,
.

Istnieją tabele ze zmiennymi krokami. Krok tabeli jest określony przez wymagania problemu inżynierskiego i nie podłączony z dokładnością znalezienia rozwiązania.

Jeśli y jest wektorem, wówczas tabela wartości rozwiązań przyjmie postać tabeli. 2.

Tabela 2

W systemie MATHCAD zamiast tabeli używana jest macierz, która jest transponowana względem określonej tabeli.

Rozwiąż problem Cauchy'ego z dokładnością ε oznacza pobranie wartości z określonej tabeli (liczby lub wektory),
, takie że
, Gdzie
- dokładne rozwiązanie. Możliwe jest, że rozwiązanie segmentu określonego w zadaniu nie będzie kontynuowane. Następnie musisz odpowiedzieć, że problemu nie można rozwiązać na całym segmencie i musisz znaleźć rozwiązanie na segmencie, w którym on istnieje, starając się, aby ten segment był jak największy.

Należy pamiętać, że jest to dokładne rozwiązanie
nie wiemy (w przeciwnym razie po co stosować metodę numeryczną?). Stopień
należy uzasadnić innymi względami. Co do zasady nie ma możliwości uzyskania 100% gwarancji, że ocena zostanie przeprowadzona. Dlatego do oszacowania wartości stosuje się algorytmy
, które okazują się skuteczne w większości problemów inżynierskich.

Ogólna zasada rozwiązywania problemu Cauchy'ego jest następująca. Odcinek [ A; B] jest podzielony na kilka segmentów przez węzły integracji. Liczba węzłów k nie musi odpowiadać liczbie węzłów M tabela końcowa wartości decyzyjnych (tabele 1, 2). Zwykle, k > M. Dla uproszczenia założymy, że odległość między węzłami jest stała,
;H zwany etapem integracji. Następnie, zgodnie z pewnymi algorytmami, znając wartości Na I < S, oblicz wartość . Im mniejszy krok H, tym niższa wartość będzie się różnić od wartości rozwiązania dokładnego
. Krok H w tym podziale jest już zdeterminowany nie wymaganiami problemu inżynierskiego, ale wymaganą dokładnością rozwiązania problemu Cauchy'ego. Ponadto należy go tak dobrać, aby w jednym kroku znajdował się stół. 1, 2 pasują do całkowitej liczby kroków H. W tym przypadku wartości y, uzyskany w wyniku obliczeń z krokami H w punktach
, stosuje się odpowiednio w tabeli. 1 lub 2.

Najprostszym algorytmem rozwiązania problemu Cauchy'ego dla równania (7) jest metoda Eulera. Wzór obliczeniowy to:

(8)

Zobaczmy, jak oceniana jest dokładność znalezionego rozwiązania. Załóżmy, że
jest dokładnym rozwiązaniem problemu Cauchy’ego i także tym
, chociaż prawie zawsze tak nie jest. W takim razie gdzie jest stała C zależy od funkcji
w pobliżu punktu
. Zatem na jednym etapie integracji (znajdowaniu rozwiązania) pojawia się błąd kolejności . Ponieważ trzeba podjąć kroki
, to naturalne jest oczekiwanie, że całkowity błąd w ostatnim punkcie
wszystko będzie dobrze
, tj. zamówienie H. Dlatego metodę Eulera nazywa się metodą pierwszego rzędu, tj. błąd ma rząd pierwszej potęgi kroku H. W rzeczywistości na jednym etapie integracji uzasadnione może być następujące oszacowanie. Pozwalać
– dokładne rozwiązanie problemu Cauchy’ego z warunkiem początkowym
. To jasne
nie pokrywa się z wymaganym dokładnym rozwiązaniem
pierwotny problem Cauchy'ego równania (7). Jednakże w małych H i „dobrą” funkcję
te dwa dokładne rozwiązania będą się nieznacznie różnić. Zapewnia to wzór na resztę Taylora
, daje to błąd kroku całkowania. Na błąd końcowy składają się nie tylko błędy na każdym etapie całkowania, ale także odchylenia od pożądanego rozwiązania dokładnego
od dokładnych rozwiązań
,
, a odchylenia te mogą być bardzo duże. Jednak ostateczne oszacowanie błędu metody Eulera dla „dobrej” funkcji
nadal wygląda
,
.

Stosując metodę Eulera, obliczenia przebiegają w następujący sposób. Według określonej dokładności ε określić przybliżony krok
. Określanie liczby kroków
i ponownie w przybliżeniu wybierz krok
. Następnie ponownie regulujemy go w dół tak, aby na każdym kroku stół. 1 lub 2 odpowiadają całkowitej liczbie kroków całkowania. Dostajemy krok H. Zgodnie ze wzorem (8), wiedząc I , znajdujemy. Według znalezionej wartości I
znajdujemy tak dalej.

Wynikowy wynik może nie mieć i generalnie nie będzie miał pożądanej dokładności. Dlatego zmniejszamy krok o połowę i ponownie stosujemy metodę Eulera. Porównujemy wyniki pierwszego zastosowania metody i drugiego zastosowania identyczny zwrotnica . Jeśli wszystkie rozbieżności są mniejsze niż określona dokładność, ostatni wynik obliczeń można uznać za odpowiedź na problem. Jeśli nie, to ponownie zmniejszamy krok o połowę i ponownie stosujemy metodę Eulera. Teraz porównujemy wyniki ostatniego i przedostatniego zastosowania metody itp.

Metodę Eulera stosuje się stosunkowo rzadko ze względu na to, że dla osiągnięcia zadanej dokładności ε wymagana jest duża liczba kroków w kolejności
. Jednakże, jeśli
ma nieciągłości lub nieciągłe pochodne, wówczas metody wyższego rzędu dadzą ten sam błąd, co metoda Eulera. Oznacza to, że wymagana będzie taka sama ilość obliczeń, jak w metodzie Eulera.

Spośród metod wyższego rzędu najczęściej stosowana jest metoda Runge-Kutty czwartego rzędu. W nim obliczenia przeprowadza się według wzorów

Metoda ta w obecności ciągłych czwartych pochodnych funkcji
powoduje błąd w jednym kroku zamówienia , tj. w notacji wprowadzonej powyżej,
. Ogólnie rzecz biorąc, na przedziale całkowania, jeśli na tym przedziale zostanie określone rozwiązanie dokładne, błąd całkowania będzie rzędu .

Wybór stopnia całkowania odbywa się analogicznie jak w metodzie Eulera, z tą różnicą, że początkową przybliżoną wartość kroku wybiera się z zależności
, tj.
.

Większość programów używanych do rozwiązywania równań różniczkowych wykorzystuje automatyczny wybór kroków. Istota tego jest taka. Niech wartość zostanie już obliczona . Wartość jest obliczana
w przyrostach H, wybrane podczas obliczeń . Następnie wykonywane są dwa etapy całkowania za pomocą step , tj. dodany zostanie dodatkowy węzeł
pośrodku między węzłami I
. Obliczane są dwie wartości
I
w węzłach
I
. Wartość jest obliczana
, Gdzie P– kolejność metod. Jeśli δ jest mniejsza od dokładności określonej przez użytkownika, wówczas przyjmuje się, że
. Jeśli nie, wybierz nowy krok H równe i powtórz kontrolę dokładności. Jeśli podczas pierwszej kontroli δ jest znacznie mniejsza niż określona dokładność, wówczas podejmowana jest próba zwiększenia kroku. W tym celu jest obliczany
w węźle
w przyrostach H z węzła
i jest obliczany
w krokach 2 H z węzła . Wartość jest obliczana
. Jeśli mniejsza niż określona dokładność, następnie krok 2 H uznane za akceptowalne. W takim przypadku przypisywany jest nowy krok
,
,
. Jeśli większa dokładność, wówczas krok pozostaje taki sam.

Należy wziąć pod uwagę, że programy z automatycznym wyborem kroku całkowania osiągają zadaną dokładność tylko przy wykonaniu jednego kroku. Dzieje się tak ze względu na dokładność przybliżenia rozwiązania przechodzącego przez punkt
, tj. przybliżenie rozwiązania
. Takie programy nie uwzględniają ile rozwiązania
różni się od pożądanego rozwiązania
. Dlatego nie ma gwarancji, że określona dokładność zostanie osiągnięta w całym przedziale całkowania.

Opisane metody Eulera i Runge-Kutty należą do grupy metod jednoetapowych. Oznacza to, że należy obliczyć
w tym punkcie
wystarczy znać znaczenie w węźle . Naturalnym jest oczekiwanie, że w przypadku wykorzystania większej ilości informacji na temat decyzji, uwzględnionych zostanie kilka jej poprzednich wartości
,
itd., a następnie nową wartość
będzie można znaleźć dokładniej. Strategię tę stosuje się w metodach wieloetapowych. Aby je opisać, wprowadzamy notację
.

Przedstawicielami metod wieloetapowych są metody Adamsa – Bashfortha:

Metoda k-ta kolejność powoduje błąd kolejności lokalnej
lub globalny – porządek .

Metody te należą do grupy metod ekstrapolacyjnych, tj. nowe znaczenie jest wyraźnie wyrażone poprzez poprzednie. Innym rodzajem są metody interpolacyjne. W nich na każdym kroku trzeba rozwiązać równanie nieliniowe dla nowej wartości . Weźmy jako przykład metody Adamsa – Moultona:

Aby skorzystać z tych metod, musisz znać kilka wartości na początku liczenia
(ich liczba zależy od kolejności metody). Wartości te należy uzyskać innymi metodami, na przykład metodą Runge-Kutty z małym krokiem (w celu zwiększenia dokładności). Metody interpolacyjne w wielu przypadkach okazują się bardziej stabilne i pozwalają na podjęcie większych kroków niż metody ekstrapolacyjne.

Aby nie rozwiązywać równania nieliniowego na każdym etapie metod interpolacji, stosuje się metody korekcji predyktorów Adamsa. Najważniejsze jest to, że metodę ekstrapolacji stosuje się najpierw na etapie i na wynikowej wartości
jest podstawiony po prawej stronie metody interpolacji. Na przykład w metodzie drugiego rzędu

Główne zagadnienia poruszane na wykładzie:

1. Opis problemu

2. Metoda Eulera

3. Metody Runge-Kutty

4. Metody wieloetapowe

5. Rozwiązanie problemu brzegowego liniowego równania różniczkowego II rzędu

6. Numeryczne rozwiązanie równań różniczkowych cząstkowych

1. Opis problemu

Najprostsze równanie różniczkowe zwyczajne (ODE) to równanie pierwszego rzędu rozwiązywane ze względu na pochodną: y” = f (x, y) (1). Główny problem związany z tym równaniem znany jest jako problem Cauchy’ego: znajdź a rozwiązanie równania (1) w postaci funkcji y (x), spełniającej warunek początkowy: y (x0) = y0 (2).
DE n-tego rzędu y (n) = f (x, y, y",:, y(n-1)), dla którego problemem Cauchy'ego jest znalezienie rozwiązania y = y(x) spełniającego warunki początkowe:
y (x0) = y0 , y" (x0) = y"0 , :, y(n-1)(x0) = y(n-1)0 , gdzie y0 , y"0 , :, y(n- 1)0 - podane liczby można sprowadzić do układu DE pierwszego rzędu.

· Metoda Eulera

Metoda Eulera opiera się na idei graficznego skonstruowania rozwiązania równania różniczkowego, ale ta sama metoda zapewnia również postać numeryczną pożądanej funkcji. Niech będzie dane równanie (1) z warunkiem początkowym (2).
Uzyskanie tablicy wartości pożądanej funkcji y(x) metodą Eulera polega na cyklicznym zastosowaniu wzoru: , i = 0, 1, :, n. Aby geometrycznie skonstruować łamaną Eulera (patrz rysunek), wybieramy biegun A(-1,0) i na osi rzędnych nanosimy odcinek PL=f(x0, y0) (punkt P jest początkiem współrzędnych). Oczywiście współczynnik kątowy promienia AL będzie równy f(x0, y0), zatem aby otrzymać pierwsze ogniwo łamanej Eulera wystarczy poprowadzić prostą MM1 z punktu M równoległego do półprostej AL aż przetnie się z linią prostą x = x1 w pewnym punkcie M1(x1, y1). Przyjmując punkt M1(x1, y1) jako początkowy, nanosimy odcinek PN = f (x1, y1) na oś Oy i rysujemy prostą przez punkt M1 M1M2 | | AN aż do przecięcia w punkcie M2(x2, y2) z prostą x = x2 itd.

Wady metody: niska dokładność, systematyczne kumulowanie się błędów.

· Metody Runge-Kutty

Główna idea metody: zamiast używać pochodnych cząstkowych funkcji f (x, y) we wzorach roboczych, używaj tylko samej tej funkcji, ale na każdym kroku obliczaj jej wartości w kilku punktach. W tym celu będziemy szukać rozwiązania równania (1) w postaci:


Zmieniając α, β, r, q, otrzymamy różne wersje metod Runge-Kutty.
Dla q=1 otrzymujemy wzór Eulera.
Przy q=2 i r1=r2=½ otrzymujemy, że α, β= 1 i stąd mamy wzór: , który nazywamy ulepszoną metodą Eulera-Cauchy’ego.
Dla q=2 i r1=0, r2=1 otrzymujemy, że α, β = ½ i stąd mamy wzór: - druga ulepszona metoda Eulera-Cauchy'ego.
Dla q=3 i q=4 istnieją także całe rodziny formuł Runge-Kutty. W praktyce stosuje się je najczęściej, ponieważ nie zwiększaj błędów.
Rozważmy schemat rozwiązywania równania różniczkowego za pomocą metody Runge-Kutty o czwartym rzędzie dokładności. Obliczenia przy zastosowaniu tej metody przeprowadza się według wzorów:

Wygodnie jest uwzględnić je w poniższej tabeli:

· Metody wieloetapowe

Omówione powyżej metody to tzw. metody stopniowego całkowania równania różniczkowego. Charakteryzują się tym, że wartości rozwiązania w kolejnym kroku szuka się korzystając z rozwiązania otrzymanego tylko w jednym poprzednim kroku. Są to tak zwane metody jednoetapowe.
Główną ideą metod wieloetapowych jest wykorzystanie kilku poprzednich wartości rozwiązania przy obliczaniu wartości rozwiązania w kolejnym kroku. Metody te nazywane są również metodami m-etapowymi w oparciu o liczbę m używaną do obliczenia wartości poprzedniego rozwiązania.
W ogólnym przypadku, aby określić przybliżone rozwiązanie yi+1, m-krokowe schematy różnicowe zapisuje się w następujący sposób (m 1):
Rozważmy konkretne formuły, które implementują najprostsze jawne i ukryte metody Adamsa.

Jawna metoda Adamsa drugiego rzędu (2-etapowa jawna metoda Adamsa)

Mamy a0 = 0, m = 2.
Są to zatem wzory obliczeniowe jawnej metody Adamsa II rzędu.
Dla i = 1 mamy niewiadomą y1, którą znajdziemy metodą Runge-Kutty dla q = 2 lub q = 4.
Dla i = 2, 3: znane są wszystkie niezbędne wartości.

Ukryta metoda Adamsa pierwszego rzędu

Mamy: a0 0, m = 1.
Są to zatem wzory obliczeniowe ukrytej metody Adamsa pierwszego rzędu.
Główny problem ze schematami ukrytymi jest następujący: yi+1 jest zawarte zarówno po prawej, jak i lewej stronie przedstawionej równości, więc mamy równanie na znalezienie wartości yi+1. Równanie to jest nieliniowe i zapisane w postaci odpowiedniej dla rozwiązania iteracyjnego, dlatego do jego rozwiązania zastosujemy prostą metodę iteracyjną:
Jeśli krok h zostanie wybrany dobrze, proces iteracyjny szybko osiągnie zbieżność.
Ta metoda również nie uruchamia się automatycznie. Aby obliczyć y1, musisz znać y1(0). Można to znaleźć za pomocą metody Eulera.