Gauss – Seidelova metoda - Gauss–Seidel method

v numerická lineární algebra, Gauss – Seidelova metoda, také známý jako Liebmannova metoda nebo metoda postupného přemístění, je iterační metoda slouží k řešení a soustava lineárních rovnic. Je pojmenována po Němec matematici Carl Friedrich Gauss a Philipp Ludwig von Seidel, a je podobný Jacobiho metoda. Ačkoli to může být aplikováno na jakoukoli matici s nenulovými prvky na úhlopříčkách, konvergence je zaručena pouze v případě, že matice je buď striktně diagonálně dominantní,^[1] nebo symetrický a pozitivní určitý. Bylo to zmíněno pouze v soukromém dopise Gausse jeho studentovi Gerling v roce 1823.^[2] Publikace nebyla doručena před rokem 1874 Seidelem.

Popis

Metoda Gauss – Seidel je iterační technika pro řešení čtvercového systému n lineární rovnice s neznámým X:

${ displaystyle A mathbf {x} = mathbf {b}}$.

Je definována iterací

${ displaystyle L _ {*} mathbf {x} ^ {(k + 1)} = mathbf {b} -U mathbf {x} ^ {(k)},}$

kde ${ displaystyle mathbf {x} ^ {(k)}}$ je kth aproximace nebo iterace ${ displaystyle mathbf {x}, , mathbf {x} ^ {(k + 1)}}$ je další nebo k + 1 iterace ${ displaystyle mathbf {x}}$a matice A je rozložen na a spodní trojúhelníkový součástka ${ displaystyle L _ {*}}$a přísně horní trojúhelníkový součástka U: ${ displaystyle A = L _ {*} + U}$.^[3]

Podrobněji napište A, X a b v jejich součástech:

${ displaystyle A = { begin {bmatrix} a_ {11} & a_ {12} & cdots & a_ {1n} a_ {21} & a_ {22} & cdots & a_ {2n} vdots & vdots & ddots & vdots a_ {n1} & a_ {n2} & cdots & a_ {nn} end {bmatrix}}, qquad mathbf {x} = { begin {bmatrix} x_ {1} x_ {2} vdots x_ {n} end {bmatrix}}, qquad mathbf {b} = { begin {bmatrix} b_ {1} b_ {2} vdots b_ {n} end {bmatrix}}.}$

Pak rozklad A do své spodní trojúhelníkové složky a její přísně horní trojúhelníkové složky je dána vztahem:

${ displaystyle A = L _ {*} + U qquad { text {kde}} qquad L _ {*} = { begin {bmatrix} a_ {11} & 0 & cdots & 0 a_ {21} & a_ {22 } & cdots & 0 vdots & vdots & ddots & vdots a_ {n1} & a_ {n2} & cdots & a_ {nn} end {bmatrix}}, quad U = { begin { bmatrix} 0 & a_ {12} & cdots & a_ {1n} 0 & 0 & cdots & a_ {2n} vdots & vdots & ddots & vdots 0 & 0 & cdots & 0 end {bmatrix}}.}$

Systém lineárních rovnic lze přepsat jako:

${ displaystyle L _ {*} mathbf {x} = mathbf {b} -U mathbf {x}}$

Metoda Gauss – Seidel nyní řeší levou stranu tohoto výrazu pro X, s použitím předchozí hodnoty pro X na pravé straně. Analyticky to lze napsat jako:

${ displaystyle mathbf {x} ^ {(k + 1)} = L _ {*} ^ {- 1} ( mathbf {b} -U mathbf {x} ^ {(k)}).}$

Avšak využitím trojúhelníkového tvaru ${ displaystyle L _ {*}}$, prvky X^(k+1) lze vypočítat postupně pomocí dopředné střídání:

${ displaystyle x_ {i} ^ {(k + 1)} = { frac {1} {a_ {ii}}} vlevo (b_ {i} - součet _ {j = 1} ^ {i-1 } a_ {ij} x_ {j} ^ {(k + 1)} - sum _ {j = i + 1} ^ {n} a_ {ij} x_ {j} ^ {(k)} vpravo), quad i = 1,2, dots, n.}$ ^[4]

Postup obvykle pokračuje, dokud změny provedené iterací nebudou pod určitou tolerancí, například dostatečně malou reziduální.

Diskuse

Elementární vzorec pro metodu Gauss-Seidel je extrémně podobný jako u metody Jacobiho metoda.

Výpočet X^(k+1) používá prvky X^(k+1) , které již byly vypočítány, a pouze prvky X^(k) které nebyly vypočteny v iteraci k + 1. To znamená, že na rozdíl od Jacobiho metody je vyžadován pouze jeden úložný vektor, protože prvky mohou být při výpočtu přepsány, což může být výhodné pro velmi velké problémy.

Na rozdíl od Jacobiho metody však výpočty pro každý prvek nelze provést v paralelní. Kromě toho jsou hodnoty při každé iteraci závislé na pořadí původních rovnic.

Gauss-Seidel je stejný jako SOR (postupná nadměrná relaxace) s ${ displaystyle omega = 1}$.

Konvergence

Konvergenční vlastnosti metody Gauss – Seidel jsou závislé na matici A. Je známo, že postup konverguje, pokud:

• A je symetrický pozitivní-definitivní,^[5] nebo
• A je přísně nebo neredukovatelně diagonálně dominantní.^[6]

Metoda Gauss-Seidel někdy konverguje, i když tyto podmínky nejsou splněny.

Algoritmus

Vzhledem k tomu, že prvky lze přepsat při jejich výpočtu v tomto algoritmu, je potřeba pouze jeden úložný vektor a vektorová indexace je vynechána. Algoritmus je následující:

Vstupy: A, bVýstup: ${ displaystyle  phi}$Vyberte počáteční odhad ${ displaystyle  phi}$ k řešeníopakovat až do konvergence pro i z 1 dokud n dělat        ${ displaystyle  sigma  leftarrow 0}$        pro j z 1 dokud n dělat            -li j ≠ i pak                ${ displaystyle  sigma  leftarrow  sigma + a_ {ij}  phi _ {j}}$            skončit, pokud        konec (j-smyčka) ${ displaystyle  phi _ {i}  leftarrow { frac {1} {a_ {ii}}} (b_ {i} -  sigma)}$    konec (i-loop) zkontrolujte, zda je dosaženo konvergencekonec (opakovat)

Příklady

Příklad pro maticovou verzi

Lineární systém zobrazený jako ${ displaystyle A mathbf {x} = mathbf {b}}$ darováno:

${ displaystyle A = { begin {bmatrix} 16 & 3 7 & -11 end {bmatrix}}}$ a ${ displaystyle b = { begin {bmatrix} 11 13 end {bmatrix}}.}$

Chceme použít rovnici

${ displaystyle mathbf {x} ^ {(k + 1)} = L _ {*} ^ {- 1} ( mathbf {b} -U mathbf {x} ^ {(k)})}$

ve formě

${ displaystyle mathbf {x} ^ {(k + 1)} = T mathbf {x} ^ {(k)} + C}$

kde:

${ displaystyle T = -L _ {*} ^ {- 1} U}$ a ${ displaystyle C = L _ {*} ^ {- 1} mathbf {b}.}$

Musíme se rozložit ${ displaystyle A _ {} ^ {}}$ do součtu spodní trojúhelníkové složky ${ displaystyle L _ {*} ^ {}}$ a přísná horní trojúhelníková složka ${ displaystyle U _ {} ^ {}}$:

${ displaystyle L _ {*} = { začátek {bmatrix} 16 a 0 7 & -11 konec {bmatrix}}}$ a ${ displaystyle U = { begin {bmatrix} 0 & 3 0 & 0 end {bmatrix}}.}$

Inverzní z ${ displaystyle L _ {*} ^ {}}$ je:

${ displaystyle L _ {*} ^ {- 1} = { begin {bmatrix} 16 & 0 7 & -11 end {bmatrix}} ^ {- 1} = { begin {bmatrix} 0,0625 & 0,0000 0,0398 & -0.0909 end {bmatrix}}}$.

Nyní můžeme najít:

${ displaystyle T = - { begin {bmatrix} 0,0625 & 0,0000 0,0398 & -0,0909 end {bmatrix}} times { begin {bmatrix} 0 & 3 0 & 0 end {bmatrix}} = { begin {bmatrix} 0,000 & -0,1875 0,000 & -0,1194 end {bmatrix}},}$

${ displaystyle C = { begin {bmatrix} 0,0625 & 0,0000 0,0398 & -0,0909 end {bmatrix}} times { begin {bmatrix} 11 13 end {bmatrix}} = { begin { bmatrix} 0,6875 - 0,7439 end {bmatrix}}.}$

Nyní máme ${ displaystyle T _ {} ^ {}}$ a ${ displaystyle C _ {} ^ {}}$ a můžeme je použít k získání vektorů ${ displaystyle mathbf {x}}$ iterativně.

Nejprve si musíme vybrat ${ displaystyle mathbf {x} ^ {(0)}}$: můžeme jen hádat. Čím lepší bude odhad, tím rychleji bude algoritmus fungovat.

Předpokládáme:

${ displaystyle x ^ {(0)} = { begin {bmatrix} 1,0 1,0 end {bmatrix}}.}$

Poté můžeme vypočítat:

${ displaystyle x ^ {(1)} = { begin {bmatrix} 0,000 & -0,1875 0,000 & -0,1193 end {bmatrix}} times { begin {bmatrix} 1,0 1,0 end {bmatrix} } + { begin {bmatrix} 0,6875 - 0,7443 end {bmatrix}} = { begin {bmatrix} 0,5000 - 0,8636 end {bmatrix}}.}$

${ displaystyle x ^ {(2)} = { begin {bmatrix} 0,000 & -0,1875 0,000 & -0,1193 end {bmatrix}} times { begin {bmatrix} 0,5000 - 0,8636 end {bmatrix }} + { begin {bmatrix} 0,6875 - 0,7443 end {bmatrix}} = { begin {bmatrix} 0,8494 - 0,6413 end {bmatrix}}.}$

${ displaystyle x ^ {(3)} = { begin {bmatrix} 0,000 & -0,1875 0,000 & -0,1193 end {bmatrix}} times { begin {bmatrix} 0,8494 - 0,6413 end {bmatrix}} + { begin {bmatrix} 0,6875 - 0,7443 end {bmatrix}} = { begin {bmatrix} 0,8077 - 0,6678 end {bmatrix}}.}$

${ displaystyle x ^ {(4)} = { begin {bmatrix} 0,000 & -0,1875 0,000 & -0,1193 end {bmatrix}} times { begin {bmatrix} 0,8077 - 0,6678 end {bmatrix }} + { begin {bmatrix} 0,6875 - 0,7443 end {bmatrix}} = { begin {bmatrix} 0,8127 - 0,6646 end {bmatrix}}.}$

${ displaystyle x ^ {(5)} = { begin {bmatrix} 0,000 & -0,1875 0,000 & -0,1193 end {bmatrix}} times { begin {bmatrix} 0,8127 - 0,6646 end {bmatrix }} + { begin {bmatrix} 0,6875 - 0,7443 end {bmatrix}} = { begin {bmatrix} 0,8121 - 0,6650 end {bmatrix}}.}$

${ displaystyle x ^ {(6)} = { begin {bmatrix} 0,000 & -0,1875 0,000 & -0,1193 end {bmatrix}} times { begin {bmatrix} 0,8121 - 0,6650 end {bmatrix }} + { begin {bmatrix} 0,6875 - 0,7443 end {bmatrix}} = { begin {bmatrix} 0,8122 - 0,6650 end {bmatrix}}.}$

${ displaystyle x ^ {(7)} = { begin {bmatrix} 0,000 & -0,1875 0,000 & -0,1193 end {bmatrix}} times { begin {bmatrix} 0,8122 - 0,6650 end {bmatrix }} + { begin {bmatrix} 0,6875 - 0,7443 end {bmatrix}} = { begin {bmatrix} 0,8122 - 0,6650 end {bmatrix}}.}$

Podle očekávání algoritmus konverguje k přesnému řešení:

${ displaystyle mathbf {x} = A ^ {- 1} mathbf {b} přibližně { begin {bmatrix} 0,8122 - 0,6650 end {bmatrix}}.}$

Ve skutečnosti je matice A striktně úhlopříčně dominantní (ale ne kladně definitivní).

Další příklad pro maticovou verzi

Další lineární systém zobrazený jako ${ displaystyle A mathbf {x} = mathbf {b}}$ darováno:

${ displaystyle A = { begin {bmatrix} 2 & 3 5 & 7 end {bmatrix}}}$ a ${ displaystyle b = { begin {bmatrix} 11 13 konec {bmatrix}}.}$

Chceme použít rovnici

${ displaystyle mathbf {x} ^ {(k + 1)} = L _ {*} ^ {- 1} ( mathbf {b} -U mathbf {x} ^ {(k)})}$

ve formě

${ displaystyle mathbf {x} ^ {(k + 1)} = T mathbf {x} ^ {(k)} + C}$

kde:

${ displaystyle T = -L _ {*} ^ {- 1} U}$ a ${ displaystyle C = L _ {*} ^ {- 1} mathbf {b}.}$

Musíme se rozložit ${ displaystyle A _ {} ^ {}}$ do součtu spodní trojúhelníkové složky ${ displaystyle L _ {*} ^ {}}$ a přísná horní trojúhelníková složka ${ displaystyle U _ {} ^ {}}$:

${ displaystyle L _ {*} = { začátek {bmatrix} 2 a 0 5 a 7 konec {bmatrix}}}$ a ${ displaystyle U = { begin {bmatrix} 0 & 3 0 & 0 end {bmatrix}}.}$

Inverzní z ${ displaystyle L _ {*} ^ {}}$ je:

${ displaystyle L _ {*} ^ {- 1} = { begin {bmatrix} 2 & 0 5 & 7 end {bmatrix}} ^ {- 1} = { begin {bmatrix} 0,500 & 0,000 - 0,357 a 0,143 konec {bmatrix}}}$.

Nyní můžeme najít:

${ displaystyle T = - { begin {bmatrix} 0,500 & 0,000 - 0,357 & 0,143 end {bmatrix}} times { begin {bmatrix} 0 & 3 0 & 0 end {bmatrix} } = { begin {bmatrix} 0,000 & -1,500 0,000 & 1,071 end {bmatrix}},}$

${ displaystyle C = { begin {bmatrix} 0,500 & 0,000 - 0,357 & 0,143 konec {bmatrix}} krát { begin {bmatrix} 11 13 konec {bmatrix}} = { begin {bmatrix} 5 500 - 2,071 end {bmatrix}}.}$

Nyní máme ${ displaystyle T _ {} ^ {}}$ a ${ displaystyle C _ {} ^ {}}$ a můžeme je použít k získání vektorů ${ displaystyle mathbf {x}}$ iterativně.

Nejprve si musíme vybrat ${ displaystyle mathbf {x} ^ {(0)}}$: můžeme jen hádat. Čím lepší bude odhad, tím rychleji bude algoritmus proveden.

Předpokládáme:

${ displaystyle x ^ {(0)} = { začátek {bmatrix} 1,1 2,3 konec {bmatrix}}.}$

Poté můžeme vypočítat:

${ displaystyle x ^ {(1)} = { begin {bmatrix} 0 & -1,500 0 & 1,071 end {bmatrix}} times { begin {bmatrix} 1,1 2,3 konec { bmatrix}} + { begin {bmatrix} 5 500 - 2,071 end {bmatrix}} = { begin {bmatrix} 2,050 0,393 end {bmatrix}}.}$

${ displaystyle x ^ {(2)} = { begin {bmatrix} 0 & -1,500 0 & 1,071 end {bmatrix}} times { begin {bmatrix} 2,050 0,393 end { bmatrix}} + { begin {bmatrix} 5 500 - 2,071 end {bmatrix}} = { begin {bmatrix} 4,911 - 1,651 end {bmatrix}}.}$

${ displaystyle x ^ {(3)} = cdots. ,}$

Pokud otestujeme konvergenci, zjistíme, že se algoritmus odchyluje. Matice A ve skutečnosti není ani diagonálně dominantní, ani kladně definitivní. Pak konvergence k přesnému řešení

${ displaystyle mathbf {x} = A ^ {- 1} mathbf {b} = { begin {bmatrix} -38 29 end {bmatrix}}}$

není zaručeno a v tomto případě k němu nedojde.

Příklad pro verzi rovnice

Předpokládejme k rovnice kde X_n jsou vektory těchto rovnic a výchozí bod X₀.Od první rovnice řešení pro X₁ ve smyslu ${ displaystyle x_ {n + 1}, x_ {n + 2}, dots, x_ {n}.}$ Pro další rovnice nahraďte předchozí hodnotyXs.

Aby bylo jasné, zvažte příklad.

${ displaystyle { begin {array} {rrrrl} 10x_ {1} & - x_ {2} & + 2x_ {3} && = 6, - x_ {1} & + 11x_ {2} & - x_ {3 } & + 3x_ {4} & = 25, 2x_ {1} & - x_ {2} & + 10x_ {3} & - x_ {4} & = - 11, & 3x_ {2} & - x_ { 3} & + 8x_ {4} & = 15. end {pole}}}$

Řešení pro ${ displaystyle x_ {1}, x_ {2}, x_ {3}}$ a ${ displaystyle x_ {4}}$ dává:

${ displaystyle { begin {aligned} x_ {1} & = x_ {2} / 10-x_ {3} / 5 + 3/5, x_ {2} & = x_ {1} / 11 + x_ { 3} / 11-3x_ {4} / 11 + 25/11, x_ {3} & = - x_ {1} / 5 + x_ {2} / 10 + x_ {4} / 10-11 / 10, x_ {4} & = - 3x_ {2} / 8 + x_ {3} /8+15/8.end {zarovnáno}}}$

Předpokládejme, že jsme si vybrali (0, 0, 0, 0) jako počáteční aproximace je první přibližné řešení dáno vztahem

${ displaystyle { begin {aligned} x_ {1} & = 3/5 = 0,6, x_ {2} & = (3/5) /11+25/11=3/55+25/11=2,3272 , x_ {3} & = - (3/5) / 5 + (2,3272) / 10-11/10=-3/25+0,23272-1.1=-0.9873, x_ {4} & = - 3 (2.3272) / 8 + (- 0,9873) /8+15/8=0,8789.end {zarovnáno}}}$

Pomocí získaných aproximací se iterační postup opakuje, dokud není dosaženo požadované přesnosti. Následují přibližná řešení po čtyřech iteracích.

${ displaystyle { begin {pole} {llll} x_ {1} & x_ {2} & x_ {3} & x_ {4} hline 0,6 & 2,32727 & -0,987273 & 0,878864 1,03018 & 2,03694 & -1,01446 & 0 0,984341 1,00659 & 2,00356 & -1,00253 & 0,998351 1,00086 & 2,0003 & -1,00031 & 0,99985 end {pole}}}$

Přesné řešení systému je (1, 2, −1, 1).

Příklad použití Pythonu a NumPy

Následující numerický postup jednoduše iteruje, aby vytvořil vektor řešení.

import numpy tak jako npITERATION_LIMIT = 1000# inicializovat maticiA = np.pole([[10., -1., 2., 0.],              [-1., 11., -1., 3.],              [2., -1., 10., -1.],              [0., 3., -1., 8.]])# inicializovat vektor RHSb = np.pole([6.0, 25.0, -11.0, 15.0])tisk("Systém rovnic:")pro i v rozsah(A.tvar[0]):    řádek = ["{0: 3 g}*X{1}".formát(A[i, j], j + 1) pro j v rozsah(A.tvar[1])]    tisk("[{0}] = [{1: 3 g}]".formát(" + ".připojit se(řádek), b[i]))X = np.zeros_like(b)pro it_count v rozsah(1, ITERATION_LIMIT):    x_new = np.zeros_like(X)    tisk("Opakování {0}: {1}".formát(it_count, X))    pro i v rozsah(A.tvar[0]):        s1 = np.tečka(A[i, :i], x_new[:i])        s2 = np.tečka(A[i, i + 1 :], X[i + 1 :])        x_new[i] = (b[i] - s1 - s2) / A[i, i]    -li np.úplně zavřít(X, x_new, rtol=1e-8):        přestávka    X = x_newtisk("Řešení: {0}".formát(X))chyba = np.tečka(A, X) - btisk("Chyba: {0}".formát(chyba))

Produkuje výstup:

Systém z rovnice:[ 10*x1 +  -1*x2 +   2*x3 +   0*x4] = [  6][ -1*x1 +  11*x2 +  -1*x3 +   3*x4] = [ 25][  2*x1 +  -1*x2 +  10*x3 +  -1*x4] = [-11][  0*x1 +   3*x2 +  -1*x3 +   8*x4] = [ 15]Opakování 1: [ 0.  0.  0.  0.]Opakování 2: [ 0.6         2.32727273 -0.98727273  0.87886364]Opakování 3: [ 1.03018182  2.03693802 -1.0144562   0.98434122]Opakování 4: [ 1.00658504  2.00355502 -1.00252738  0.99835095]Opakování 5: [ 1.00086098  2.00029825 -1.00030728  0.99984975]Opakování 6: [ 1.00009128  2.00002134 -1.00003115  0.9999881 ]Opakování 7: [ 1.00000836  2.00000117 -1.00000275  0.99999922]Opakování 8: [ 1.00000067  2.00000002 -1.00000021  0.99999996]Opakování 9: [ 1.00000004  1.99999999 -1.00000001  1.        ]Opakování 10: [ 1.  2. -1.  1.]Řešení: [ 1.  2. -1.  1.]Chyba: [  2.06480930e-08  -1.25551054e-08   3.61417563e-11   0,00000000e + 00]

Program řešení libovolného č. rovnic pomocí Matlabu

Následující kód používá vzorec${ displaystyle x_ {i} ^ {(k + 1)} = { frac {1} {a_ {ii}}} vlevo (b_ {i} - součet _ {j i} a_ {ij} x_ {j} ^ {(k)} vpravo), quad i, j = 1,2, ldots , n}$

funkceX =gauss_seidel(A, b, x, iters)pro i = 1:iters        pro j = 1:velikost(A,1)            X(j) = (1/A(j,j)) * (b(j) - A(j,:)*X + A(j,j)*X(j));        konec    koneckonec

Viz také

• Postupná nadměrná relaxace
• Kaczmarzova metoda („řádkově orientovaná“ metoda, zatímco Gauss-Seidel je „sloupově orientovaná“. Viz např. tento papír.)
• Iterační metoda. Lineární systémy
• Šíření gaussovské víry
• Rozdělení matice
• Richardsonova iterace

Poznámky

Reference

• Gauss, Carl Friedrich (1903), Werke (v němčině), 9, Göttingen: Köninglichen Gesellschaft der Wissenschaften.
• Golub, Gene H.; Van Loan, Charles F. (1996), Maticové výpočty (3. vyd.), Baltimore: Johns Hopkins, ISBN  978-0-8018-5414-9.
• Black, Noel & Moore, Shirley. „Gauss-Seidelova metoda“. MathWorld.

Tento článek včlení text z článku Gauss-Seidel_method na CFD-Wiki to je pod GFDL licence.

externí odkazy

• "Seidelova metoda", Encyclopedia of Mathematics, Stiskněte EMS, 2001 [1994]
• Gauss – Seidel z www.math-linux.com
• Gauss – Seidel Z institutu holistických numerických metod
• Iterace Gauss Siedel z www.geocities.com
• Gauss-Seidelova metoda
• Bickson
• Kód Matlabu
• Příklad kódu C.