Wpisany przez Tomasz Lubiński,
13 marca 2006 18:42
Metoda Gaussa - Seidela jest metodą iteracyjną i pozwala nam obliczyć układ n równań z n niewiadomymi Ax = b. Wektor x0 będący początkowym przybliżeniem rozwiązania układu będzie dany (zwykle przyjmuje się go jako wektor złożony z samych zer). By zastosować tą metodę należy najpierw tak zamienić kolejność równań układu, aby na głównej przekątnej były elementy różne od zera.
Na początku macierz współczynników A rozłożymy na sumę trzech macierzy A = L + D + U, gdzie L jest macierzą w której znajdują się elementy których numer wiersza jest większy od numeru kolumny, D to macierz diagonalna z elementami tylko na głównej przekątnej, a U to macierz, w której znajdują się elementy których numery wiersza są mniejsze od numerów kolumny. Można to zapisać następująco:
Obliczymy następujący układ równań:
4x1 - x2 - 0.2x3 + 2x4 = 30
-1x1 + 5x2 - 2x4 = 0
0.2x1 + x2 + 10x3 - x4 = -10
- 2x2 - x3 + 4x4 = 5
Zapiszmy go teraz w postaci Ax = b
Obliczmy pierwszą iterację metody, według przytoczonego na początku wzoru:
x11 = 7,5 + 0,25x20 + 0,05x30 - 0,5x40
x11 = 7,5
x21 = 0 + 0,2x11 + 0,4x40
x21 = 1,5
x31 = -1 - 0,02x11 - 0,1x21 + 0,1x40
x31 = -1,30
x41 = 1,25 + 0,5x21 + 0,25x31
x41 = 1,675
Kolejna iteracja
x12 = 7,5 + 0,25x21 + 0,05x31 - 0,5x41
x12 = 6,9725
x22 = 0 + 0,2x12 + 0,4x41
x22 = 2,0645
x32 = -1 - 0,02x12 - 0,1x22 + 0,1x41
x32 = -1,1784
x42 = 1,25 + 0,5x22 + 0,25x32
x42 = 1,98765
Można teraz obliczyć kolejną iterację. Każda z nich przybliża nas do dokładnego wyniku.
Na początku macierz współczynników A rozłożymy na sumę trzech macierzy A = L + D + U, gdzie L jest macierzą w której znajdują się elementy których numer wiersza jest większy od numeru kolumny, D to macierz diagonalna z elementami tylko na głównej przekątnej, a U to macierz, w której znajdują się elementy których numery wiersza są mniejsze od numerów kolumny. Można to zapisać następująco:
\begin{bmatrix}
a_{1,1} & a_{1,2} & ... & a_{1,n}\\
a_{2,1} & a_{2,2} & ... & a_{2,n}\\
... & ... & ... & ...\\
a_{n,1} & a_{n,2} & ... & a_{n,n}
\end{bmatrix}
=
\begin{bmatrix}
0 & 0 & ... & 0\\
a_{2,1} & 0 & ... & 0\\
... & ... & ... & ...\\
a_{n,1} & a_{n,2} & ... & 0
\end{bmatrix}
+
\begin{bmatrix}
a_{1,1} & 0 & ... & 0\\
0 & a_{2,2} & ... & 0\\
... & ... & ... & ...\\
0 & 0 & ... & a_{n,n}
\end{bmatrix}
+
\begin{bmatrix}
0 & a_{1,2} & ... & a_{1,n}\\
0 & 0 & ... & a_{2,n}\\
... & ... & ... & ...\\
0 & 0 & ... & 0
\end{bmatrix}
Następnie obliczymy macierz odwrotną do macierzy D, czyli D-1. Otrzymamy ją po podniesieniu do potęgi -1 wszystkich wartości na głównej przekątnej macierzy. Po tych operacjach możemy przystąpić już do iteracyjnego obliczania kolejnych przybliżeń rozwiązania według następującego wzoru:
x^{n+1} = D^{-1}b - D^{-1}Lx^{n+1} - D^{-1}Ux^n
(indeksy n oznaczają tutaj numer iteracji)Przykład:
Obliczymy następujący układ równań:
4x1 - x2 - 0.2x3 + 2x4 = 30
-1x1 + 5x2 - 2x4 = 0
0.2x1 + x2 + 10x3 - x4 = -10
- 2x2 - x3 + 4x4 = 5
Zapiszmy go teraz w postaci Ax = b
\begin{bmatrix}
4 & -1 & -0.2 & 2\\
-1 & 5 & 0 & -2\\
0.2 & 1 & 10 & -1\\
0 & -2 & -1 & 4
\end{bmatrix}
\begin{bmatrix}
x_1\\
x_2\\
x_3\\
x_4
\end{bmatrix}
=
\begin{bmatrix}
30\\
0\\
-10\\
5
\end{bmatrix}
Podzielmy teraz macierz A na sumę macierzy L + D + U
\begin{bmatrix}
4 & -1 & -0.2 & 2\\
-1 & 5 & 0 & -2\\
0.2 & 1 & 10 & -1\\
0 & -2 & -1 & 4
\end{bmatrix}
=
\begin{bmatrix}
0 & 0 & 0 & 0\\
-1 & 0 & 0 & 0\\
0.2 & 1 & 0 & 0\\
0 & -2 & -1 & 0
\end{bmatrix}
+
\begin{bmatrix}
4 & 0 & 0 & 0\\
0 & 5 & 0 & 0\\
0 & 0 & 10 & 0\\
0 & 0 & 0 & 4
\end{bmatrix}
+
\begin{bmatrix}
0 & -1 & -0.2 & 2\\
0 & 0 & 0 & -2\\
0 & 0 & 0 & -1\\
0 & 0 & 0 & 0
\end{bmatrix}
Obliczmy teraz macierz D-1
\begin{bmatrix}
4 & 0 & 0 & 0\\
0 & 5 & 0 & 0\\
0 & 0 & 10 & 0\\
0 & 0 & 0 & 4
\end{bmatrix}
=
\begin{bmatrix}
0.25 & 0 & 0 & 0\\
0 & 0.2 & 0 & 0\\
0 & 0 & 0.1 & 0\\
0 & 0 & 0 & 0.25
\end{bmatrix}^{-1}
A teraz kolejno D-1b, D-1L, D-1U
\begin{bmatrix}
7.5\\
0\\
-1\\
1.25
\end{bmatrix}
\begin{bmatrix}
0 & 0 & 0 & 0\\
-0.2 & 0 & 0 & 0\\
0.02 & 0.1 & 0 & 0\\
0 & -0.5 & -0.25 & 0
\end{bmatrix}
\begin{bmatrix}
0 & -0.25 & -0.05 & 0.5\\
0 & 0 & 0 & -0.4\\
0 & 0 & 0 & -0.1\\
0 & 0 & 0 & 0
\end{bmatrix}
Rozpoczynamy od zerowego przybliżenia czyli x10 = 0, x20 = 0, x30 = 0, x40 = 0Obliczmy pierwszą iterację metody, według przytoczonego na początku wzoru:
x11 = 7,5 + 0,25x20 + 0,05x30 - 0,5x40
x11 = 7,5
x21 = 0 + 0,2x11 + 0,4x40
x21 = 1,5
x31 = -1 - 0,02x11 - 0,1x21 + 0,1x40
x31 = -1,30
x41 = 1,25 + 0,5x21 + 0,25x31
x41 = 1,675
Kolejna iteracja
x12 = 7,5 + 0,25x21 + 0,05x31 - 0,5x41
x12 = 6,9725
x22 = 0 + 0,2x12 + 0,4x41
x22 = 2,0645
x32 = -1 - 0,02x12 - 0,1x22 + 0,1x41
x32 = -1,1784
x42 = 1,25 + 0,5x22 + 0,25x32
x42 = 1,98765
Można teraz obliczyć kolejną iterację. Każda z nich przybliża nas do dokładnego wyniku.
Implementacje
| Autor | Język programowania | Komentarz | Otwórz | Pobierz | Ocena |
| Tomasz Lubiński | C/C++ | .cpp | .cpp | ***** / 21 | |
| Marian | C/C++ | C++ | .cpp | .cpp | ***** / 10 |
| Codename | C/C++ | .cpp | .cpp | ***** / 4 | |
| Tomasz Lubiński | Delphi/Pascal | Borland Delphi 5 | .pas | .pas | ***** / 3 |
| Tomasz Lubiński | Java | .java | .java | ***** / 6 |
Poprawiony: 23 września 2012 14:11
x^(n+1) = D-1b - D-1Lx^(n+1) - D-1Ux^(n)
nie powinien wyglądać tak:
x^(n+1) = D-1b - D-1Lx^(n) - D-1Ux^(n)
Bo mam wrażenie, że jest tam błąd. Może mi ktoś wytłumaczyć dlaczego jest tak?
Bo z analizy kodów wynika, że powinno być tak jak na dole jak podałem.
Rozwiązaniem są liczby 5 i 2 co łatwo sprawdzić. A program niestety pokazuje inne rozwiązanie.
x^(n+1) = D-1b - D-1Lx^(n) - D-1Ux^(n) a nie tak jak tu podany:
x^(n+1) = D-1b - D-1Lx^(n+1) - D-1Ux^(n)
proszę o wytlumaczenie...
jest błąd w pierwszej iteracji 0^0=1 a tu autor przyjął 0^0=0