Metoda różnic skończonych

W analizie numerycznej The metody różnic skończonych jest powszechną techniką w celu znalezienia przybliżeniu roztworów o równań różniczkowych , które składają się na rozwiązaniu układu stosunków (schemat cyfrowego) łączący nieznanych wartości funkcji w określonych punktach wystarczająco blisko siebie. Inną .

Metoda ta wydaje się być najłatwiejsza do wdrożenia, ponieważ przebiega w dwóch etapach: z jednej strony dyskretyzacja przez skończone różnice operatorów wyprowadzenia / różniczkowania, z drugiej strony zbieżność tak otrzymanego diagramu numerycznego, gdy odległość między punkty maleją.

Aproksymacja operatorów za pomocą wzorów Taylora

Dyskretyzację operatorów różniczkowych (pierwsze pochodne, sekundy itp., Częściowe lub nie) można uzyskać za pomocą wzorów Taylora .

Sformułowanie Taylora-Younga jest preferowane w jego prostym użyciu, wyrażenie Taylora z resztą całki Laplace'a umożliwia pomiar błędów (por. Poniżej).

Przykłady przybliżeń operatorów

Przybliżenia wyśrodkowane

W punkcie $x$ i dla takiej wartości $h$ kroku dyskretyzacji, że $u$ jest trzykrotnie różniczkowalne w przedziale $[ x - h , x + h ]$ , wzór Taylora-Younga prowadzi do dwóch relacji:

{\ Displaystyle u (x + h) = u (x) + \ suma _ {n = 1} ^ {3} {{\ Frac {h ^ {n}} {n!}} u ^ {(n)} (x)} + h ^ {3} \ varepsilon _ {1} (x, h)}

{\ Displaystyle u (xh) = u (x) + \ suma _ {n = 1} ^ {3} {{\ Frac {(-h) ^ {n}} {n!}} u ^ {(n) } (x)} + h ^ {3} \ varepsilon _ {2} (x, h)}

gdzie dwie funkcje $ε i ( x , h )$ zbiegają się do 0 z $h$ . W związku z tym

{\ Displaystyle {\ Frac {u (x + h) -u (x)} {h}} = u '(x) + {\ Frac {h} {2}} u' '(x) + {\ Frac {h ^ {2}} {3!}} u ^ {(3)} (x) + h ^ {2} \ varepsilon _ {1} (x, h)}

{\ Displaystyle {\ Frac {u (xh) -u (x)} {h}} = - u '(x) + {\ Frac {h} {2}} u' '(x) - {\ Frac { h ^ {2}} {3!}} u ^ {(3)} (x) + h ^ {2} \ varepsilon _ {2} (x, h)}

odpowiadają dwóm przybliżeń $U ( x )$ w 1 ul rzędu $godzin$ .

Odejmując poprzednie zmiany, co sprowadza się do uśrednienia dwóch skończonych różnic przed i za $u ( x )$ , otrzymujemy

{\ Displaystyle {\ Frac {u (x + h) -u (xh)} {2h}} = u '(x) + {\ Frac {h ^ {2}} {6}} u ^ {(3) } (x) + h ^ {2} \ varepsilon _ {3} (x, h)}

co jest przybliżeniem $U ( x )$ od 2 e rzędu $godzin$ .

Zdecentralizowane przybliżenia Offset w górę strumienia

W punkcie $x$ i dla takiej wartości $h$ kroku dyskretyzacji, że $u$ jest trzykrotnie różniczkowalne w przedziale $[ x , x + 2 h ]$ , wzór Taylora-Younga prowadzi do zależności:

{\ Displaystyle u (x + h) = u (x) + \ suma _ {n = 1} ^ {3} {{\ Frac {h ^ {n}} {n!}} u ^ {(n)} (x)} + h ^ {3} \ varepsilon _ {1} (x, h)}

gdzie funkcja zbiega się do 0 z $h$ . W związku z tym ${\ Displaystyle \ varepsilon _ {i} (x, h)}$

{\ Displaystyle {\ Frac {u (x + h) -u (x)} {h}} = u '(x) + {\ Frac {h} {2}} u' '(x) + {\ Frac {h ^ {2}} {3!}} u ^ {(3)} (x) + h ^ {2} \ varepsilon _ {1} (x, h)}

Odpowiada to w przybliżeniu $litery U ( x )$ z 1 ul rzędu $godzin$ .

Powtarzając operację dla dalszego poza centrum, pisząc, że:

{\ Displaystyle u '' (x) = {\ Frac {{\ Frac {u (x + 2h) -u (x + h)} {h}} - {\ Frac {u (x + h) -u ( x)} {h}}} {h}} + \ varepsilon _ {2} (x, h)}

Skąd

{\ Displaystyle {\ Frac {u (x + 2h) -2u (x + h) + u (x)} {h ^ {2}}} = u '' (x) + {\ Frac {h ^ {2 }} {6}} u ^ {(3)} (x) + h ^ {2} \ varepsilon _ {3} (x, h)}

co jest przybliżeniem $U „” ( x )$ od 2 e rzędu $godzin$ .

Formuły rozszerzone na kolejne zamówienia

Rozszerzając rozmiar szablonu , można wyznaczyć skończone różnice wyższych rzędów podobnymi metodami (zwiększając rząd we wzorze Taylora i określając odpowiednią kombinację liniową, aby zlikwidować zbędne wyrażenia).

Na przykład w punkcie $x$ i dla wartości $h$ kroku dyskretyzacji takiej, że $u$ jest czterokrotnie różniczkowalne w przedziale $[ x - 2 h , x + 2 h ]$ , przez rozszerzenie wzoru Taylora możemy pokazać niż pięć diagramy punktowe

{\ Displaystyle u '(x) = {\ Frac {u (x-2h) -8u (xh) + 8u (x + h) -u (x + 2h)} {12h}} + \ varepsilon _ {1} '(x, h)}

{\ Displaystyle u '' (x) = {\ Frac {-u (x-2h) + 16u (xh) -30u (x) + 16u (x + h) -u (x + 2h)} {12h ^ { 2}}} + \ varepsilon _ {2} '(x, h)}

są przybliżeniami pierwszej i drugiej pochodnej rzędu 4.

Rozszerzenie do funkcji wielowymiarowych Offset w górę strumienia

W punkcie $( x , y )$ i dla wartości $h$ kroku dyskretyzacji (takiej samej w dwóch wymiarach) takiej, że $u ( x , y )$ jest 4 razy różniczkowalna na prostokącie $[0, x + 2 h ] \times [ 0, y + 2 h ]$ , możemy napisać

{\ Displaystyle {\ Frac {(u (x + 2h, 0) -2u (x + h, 0) + u (x, 0)) + (u (0, y + 2h) -2u (0, y + h) + u (0, y))} {h ^ {2}}} = \ Delta u (x, y) + {\ frac {h ^ {2}} {12}} (\ częściowe _ {x} ^ {4} u (x, y) + \ częściowe _ {y} ^ {4} u (x, y)) + h ^ {2} \ varepsilon _ {5} (x, y, h)}

co jest przybliżeniem Laplace'a $hemibursztynianu u ( x , y ),$ z 2 e rzędu $godzin$ (por równania Laplace'a i równania Poissona ).

Szablon wyśrodkowany

W punkcie $( x , y )$ i dla wartości $h$ kroku dyskretyzacji (taka sama w dwóch wymiarach), że $u ( x , y )$ jest 4 razy różniczkowalna na prostokącie $[ x - h , x + h ] \times [ y - h , y + h ]$ , możemy napisać

{\ Displaystyle {\ Frac {u (x + h, y) + u (xh, y) + u (x, y + h) + u (x, yh) -4u (x, y)} {h ^ { 2}}} = \ Delta u (x, y) + {\ frac {h ^ {2}} {12}} (\ częściowe _ {x} ^ {4} u (x, y) + \ częściowe _ { y} ^ {4} u (x, y)) + h ^ {2} \ varepsilon _ {5} (x, y, h)}

co jest przybliżeniem Laplace'a $hemibursztynianu u ( x , y ),$ z 2 e rzędu $godzin$ (por równania Laplace'a i równania Poissona ).

Wspomniane wyżej pojęcie „porządku” odpowiada koncepcji lokalnej zbieżności operatora dyskretyzowanego. Globalna konwergencja rozwiązania dyskretnego to zupełnie inna koncepcja, mimo że istnieje między nimi związek rodzinny.

Siatka

W przypadku metody różnic skończonych siatka jest zbiorem izolowanych punktów (zwanych węzłami ) znajdujących się w domenie definicji funkcji podlegających równaniom różniczkowym cząstkowym, których siatką na jedynych węzłach są zdefiniowane niewiadome odpowiadające przybliżonej wartości tych funkcji.

Siatka obejmuje również węzły zlokalizowane na granicy pola (lub przynajmniej „blisko” tej granicy), aby móc narzucić warunki brzegowe i / lub warunek początkowy z wystarczającą precyzją.

A priori, pierwszą cechą siatki jest jak najlepsze pokrycie pola, w którym się rozwija, ograniczenie odległości między każdym węzłem a jego najbliższym sąsiadem. Jednak siatka musi także umożliwiać wyrażenie dyskretnego sformułowania operatorów różnicowania: z tego powodu węzły siatki najczęściej znajdują się na siatce, której głównymi kierunkami są osie zmiennych.

Krok siatki nazywa się odległością między dwoma sąsiednimi węzłami położonymi na linii równoległej do jednej z osi. W tym sensie krok jest zarówno pojęciem lokalnym, jak i kierunkowym. Będziemy mówić o globalnym skoku, aby wyznaczyć największy lokalny skok , pojęcie, które pozostaje kierunkowe.

Chociaż najczęściej zachowuje się stały ton (bez teoretycznego problemu z rozdzielczością), czasami rozsądne jest wprowadzenie zmiennej wysokości, która zostanie wybrana dokładniej w obszarach, w których dokładne rozwiązanie ma większe różnice: ta sztuczka umożliwia zmniejszenie liczby niewiadomych bez uszczerbku dla dokładności wyników. Z drugiej strony sformułowanie jest nieco bardziej złożone, ponieważ dyskretyzacja operatorów różniczkowych musi to uwzględniać.

Przykłady siatek

Dla równania różniczkowego dotyczącego funkcji zmiennej, której dziedziną (in ) jest przedział $[0;$ $1]$ , siatka pak stałe charakteryzuje $M$ $+ 1$ węzłów $x$ $I$ $=$ $H$ $, 0 \leq$ $i$ $\leq$ $M$ z etapem $h$ $= 1 /$ $M$ . Ta siatka obejmuje dwa punkty graniczne $x$ $0$ i $x$ $M,$ na które nakładane są możliwe warunki brzegowe. $\ mathbb {R}$

Rozważmy równanie różniczkowe cząstkowe dotyczące funkcji dwóch zmiennych (dziedzina ): ${\ Displaystyle \ Omega \ subset \ mathbb {R} ^ {2}}$

Jeśli $Ω$ jest prostokątem $[0; 1] \times [0; 1]$ (którego boki są równoległe do osi), siatka wynikająca z siatki $( x i , y j ) = ( ih , jk ), 0 \leq i \leq M , 0 \leq j \leq N$ z krokami $h = 1 / M$ i $k = 1 / N$ to proste uogólnienie poprzedniego przypadku.
Jeśli Ω jest dyskiem wyśrodkowanym na początku io promieniu 1, rozważamy siatkę utworzoną przez węzły siatki, które znajdują się w dysku, to znaczy ( x i , y j ) ∈ Ω gdzie ( x i , r j ) = ( IH , jk ) o skoku h = 1 / M . Aby narzucić możliwe warunki brzegowe (na przykład te z Dirichleta, które ustalają wartość funkcji na ∂Ω ), rzadkie węzły znajdujące się dokładnie na granicy są zbyt mało reprezentatywne. Wskazane jest wówczas rozszerzenie właściwości „bycia na granicy” na inne węzły, które są blisko niej, np. Poprzez uwzględnienie wszystkich węzłów siatki, które nie mają czterech bezpośrednich sąsiadów. Wartości graniczne, które mają być ustawione w tych nowych węzłach granicznych, można zdefiniować na różne sposoby:
- Przyjmując wartość dokładnego problemu, który jest nałożony na najbliższy punkt $\partialΩ$ : w tym przypadku nieregularności geograficzne węzłów granicznych siatki (obserwowane wraz ze spadkiem $h$ ) generują zakłócenia rozwiązania dyskretnego, w najlepszym przypadku anomalie lokalne niezwiązane z dokładnym rozwiązaniem.
- Biorąc pod uwagę, że wartości nowych węzłów granicznych są nieznane, ale przez dodanie dodatkowych dyskretyzowanych relacji różnicowych, łączących „naturalnie” te niewiadome z wartościami węzłów sąsiednich oraz z określonymi punktami $\partialΩ$ . Jeśli podejście jest nieco bardziej złożone w realizacji, znacznie zmniejsza wadę poprzedniego.

Schemat cyfrowy

Schemat numeryczny można zdefiniować jako algebraiczne sformułowanie problemu dyskretnego zaprojektowanego przy użyciu metody różnic skończonych. Proces obejmuje następujące kroki:

Wybierz operatory dyskretne, które są przybliżeniami operatorów różniczkowych dokładnego sformułowania.
Wygeneruj siatkę domeny definicji, zwracając uwagę na węzły graniczne i sposób tłumaczenia warunków brzegowych.
Opierając się na wyrażeniach wynikających z operatorów dyskretnych, ustalić relacje wiążące wartości funkcji z węzłami siatki (nieznane czynniki).
Upewnij się, że zbiór niewiadomych i relacji, które je łączą, stanowi problem numeryczny, który nie jest przesadzony lub niedostatecznie określony. To sprawdzenie jest minimalnym warunkiem nadziei na znalezienie rozwiązania, ale nie daje żadnej gwarancji globalnej konwergencji.

Po ustaleniu diagramu numerycznego i sformułowaniu problemu dyskretnego, należy nie tylko go rozwiązać, ale także upewnić się, że rozwiązanie dyskretne jest zbieżne w kierunku rozwiązania dokładnego, gdy kroki siatki zmierzają do 0.

W przypadku niektórych tzw. Diagramów jawnych można uporządkować niewiadome w taki sposób, aby każdą z nich można było określić rekurencyjnie z poprzednich, które mają być już obliczone ( macierz trójkątna ). W przypadku schematów niejawnych czasami można uniknąć rozwiązania całego układu wszystkich równań. Dotyczy to w szczególności systemu ewoluującego, którego stan, charakteryzujący się zmiennymi przestrzennymi, jest określany przez warunki początkowe (t = 0), a następnie stopniowo ewoluuje w czasie: diagram numeryczny pozostaje jawny w zmiennej czasowej, a jej ukryty charakter dotyczy tylko zmienne przestrzenne.

We wszystkich przypadkach każde równanie na diagramie numerycznym odnosi się tylko do niewielkiej liczby niewiadomych. W środowisku liniowym właściwość ta prowadzi do sformułowania problemu dyskretnego za pomocą rzadkich macierzy i wykorzystania go do rozwiązania go odpowiednimi metodami . Ta zaleta jest niezaprzeczalna, gdy rozmiar siatki przekracza ramy studiów dydaktycznych.

Rozdzielczość diagramów numerycznych jest generalnie oparta na klasycznych metodach algebraicznych. Jednak inne równoważne sformułowania mogą wymagać metod optymalizacji .

Przykład schematu cyfrowego

Rozważ następujący problem:

{\ Displaystyle u '(x) - \ tau u (x) = 0 \, \ forall x \ w [0,1], \ u (0) = u_ {0}.}

Ten problem pozostaje akademicki, o ile znane jest dokładne rozwiązanie:

{\ Displaystyle u (x) = u_ {0} \, e ^ {\ tau x}.}

Przy jawnym schemacie Eulera rzędu 1 zastosowanym do regularnej siatki o skoku $h = 1 / M$ , niewiadome $u n$ odzwierciedlające $u ( nh )$ są powiązane relacjami

{\ Displaystyle {\ Frac {u_ {n + 1} -u_ {n}} {h}} - \ tau u_ {n} = 0, \; \ forall \, 0 \ leqslant n <M.}

Ten diagram prowadzi do relacji rekurencji

{\ Displaystyle u_ {n + 1} = (1 + h \ tau) u_ {n}}

którego jawnym rozwiązaniem jest

{\ Displaystyle u_ {n} = \ lewo (1 + {\ Frac {\ tau} {M}} \ prawo) ^ {n} u_ {0}.}

Inne sformułowanie otrzymane za pomocą wykresu rzędu 2 (z wyjątkiem węzła $n = 1,$ dla którego zachowuje się wykres rzędu 1) daje

{\ Displaystyle {\ Frac {u_ {1} -u_ {0}} {h}} - \ tau u_ {0} = 0,}

{\ Displaystyle {\ Frac {u_ {n + 1} -u_ {n-1}} {2h}} - \ tau u_ {n} = 0, \; \ forall \, 0 <n <M.}

Podobnie jak pierwszy, ten drugi diagram jest wyraźny .

Bardzo łatwo jest numerycznie określić rozwiązania tych dwóch diagramów, aby porównać je z dokładnym rozwiązaniem. Wydaje się zasadne oczekiwanie lepszych wyników na drugim diagramie, ponieważ jego kolejność jest wyższa niż na pierwszym (można wykazać, że oba diagramy cyfrowe są jednakowo zbieżne):

Ta intuicja zachowuje się, gdy $τ > 0$ .
Z drugiej strony okazuje się błędny, gdy $τ$ jest dostatecznie ujemne (np. $Τ <-3$ ). Rzeczywiście, diagram rzędu 1, który jest zastosowany do pierwszego węzła, generuje w tym przypadku początkową różnicę $u 1 - u ( h ),$ która będzie generować oscylacje na różnicach obserwowanych w kolejnych węzłach, oscylacje, które są wzmacniane aż do ostatnich węzłów: Eulera Schemat niewątpliwie generuje znacznie lepsze rozwiązanie.

Porównanie to wyraźnie pokazuje, że dobra reprezentacja operatorów różniczkowych nie jest warunkiem wystarczającym do uzyskania dobrego diagramu numerycznego.

Konwergencja

Zbieżność diagramu cyfrowego właściwość teoretycznej globalny zapewnienia, że różnica (w rozumieniu z normą ) pomiędzy przybliżeniu roztworu i dokładnie roztworem dąży do 0, gdy etap dyskretyzacji dąży do 0 (lub po każdym z etapów globalnej powiązany z różnymi kierunkami zmierza w kierunku 0).

Przybliżone rozwiązanie schematu cyfrowego pozostaje mało wiarygodne, o ile nie wykazano jego zbieżności. Dowód ten jest niewątpliwie najdelikatniejszym punktem metody drobnych różnic, w każdym razie wymagającej użycia narzędzi analitycznych .

Nie wystarczy sprawdzić na konkretnych przykładach liczbowych, że zachowanie rozwiązania dyskretnego jest zgodne z oczekiwaniami zapewniającymi zbieżność. Z drugiej strony takie przykłady mogą pomóc w udowodnieniu czegoś przeciwnego.

Koncepcyjnie różnice między rozwiązaniem przybliżonym a rozwiązaniem dokładnym przejawiają się połączeniem dwóch zjawisk:

Zmienność wywołana lokalnie przez przybliżenia właściwe operatorowi dyskretyzowanemu. Jest to pojęcie spójności lub odporności diagramu. Na przykład wysoki stopień aproksymacji (uzyskany za pomocą twierdzenia Taylora ) nie ma dużej wartości, gdy dokładne rozwiązanie nie ma wymaganej regularności.
W przypadku wyraźnego wykresu, propagacja zmian, które na etapach obliczeń łączą się z poprzednimi zmianami i mogą się stopniowo zwiększać. Jest to pojęcie stabilności diagramu, które pojawia się również wtedy, gdy jest domniemane.

Koncepcje te nie uwzględniają błędów zaokrągleń, które mogą dodatkowo komplikować sprawę, jak pokazano na poniższym rysunku, otrzymanym na konkretnym przykładzie :

Norma, dla której badana jest konwergencja, musi pozostać niezależna od etapów dyskretyzacji. Jednak często stosuje się standardy związane z przestrzeniami L p . Dla funkcji zmiennej:

Prosta zbieżność: zbieżność w dowolnym punkcie przestrzeni definicji.
Jednolita zbieżność (w L ∞ ): ${\ Displaystyle \ | u_ {h} -u \ | _ {\ infty} = h \ max _ {i} | u_ {i} -u (ih) |.}$
Zbieżność w L 2 : ${\ displaystyle \ | u_ {h} -u \ | _ {2} = h {\ sqrt {\ sum _ {i} | u_ {i} -u (ih) | ^ {2}}}.}$
Zbieżność w L 1 : ${\ Displaystyle \ | u_ {h} -u \ | _ {1} = h \ suma _ {i} | u_ {i} -u (ih) |.}$

W kontekście problemu ewolucyjnego z warunkiem początkowym , twierdzenie Laxa rygorystycznie precyzuje pojęcia spójności i stabilności , przy czym drugi jest warunkiem koniecznym i wystarczającym do zapewnienia konwergencji .

Przykład konwergencji

W ostatnim przedstawionym powyżej przykładzie, dla którego znamy jednocześnie rozwiązanie dokładne i przybliżone (diagram Eulera) , raport spełnia ${\ Displaystyle \ Displaystyle u (x) = u_ {0} e ^ {\ tau x}}$ ${\ Displaystyle \ Displaystyle u_ {n} = u_ {0} (1+ \ tau h) ^ {n}}$

{\ Displaystyle \ ln \ lewo (\ lewo | {\ Frac {u_ {n}} {u (nh)}} \ prawo | \ prawo) = n | \ ln (1+ \ tau h) - \ tau h | = n \, O (h ^ {2}) = O (h)}

która dąży do 0, gdy zmierza do 0, to jednakowo dla ${\ Displaystyle \ Displaystyle h}$ ${\ Displaystyle n \ leqslant 1 / h.}$

W ten sposób dąży równomiernie do 0, co świadczy o zbieżności tego schematu Eulera w normie ${\ Displaystyle \ Displaystyle u (nh) -u_ {n}}$ ${\ Displaystyle \ |. \ | _ {\ infty}.}$

Uwagi i odniesienia

Obserwuje się tutaj różnicę w stosunku do metody elementów skończonych, w której nieznane funkcje są definiowane w dowolnym punkcie elementów siatki (utworzonej z trójkątów, prostokątów itp.). Ponieważ jednak możliwe jest interpolowanie wyników uzyskanych za pomocą różnic skończonych, ta odmienność nie jest fundamentalna.
Nawet jeśli sformułowanie pierwotnego problemu obejmuje całki, często korzystne jest wprowadzenie dodatkowych funkcji do ich reprezentacji.

Zobacz też

Powiązane artykuły

Linki zewnętrzne