Jacobiho metoda - Jacobi method

v numerická lineární algebra, Jacobiho metoda je iterační algoritmus pro určování řešení a striktně diagonálně dominantní soustava lineárních rovnic. Každý diagonální prvek je vyřešen pro a je připojena přibližná hodnota. Proces je poté iterován, dokud konverguje. Tento algoritmus je ořezanou verzí Jacobiho transformační metoda diagonalizace matic. Metoda je pojmenována po Carl Gustav Jacob Jacobi.

Popis

Nechat

${displaystyle Amathbf {x} = mathbf {b}}$

být čtvercový systém n lineární rovnice, kde:

${displaystyle A = {egin {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} = {egin {bmatrix} x_ {1} x_ {2} vdots x_ {n} end {bmatrix}}, qquad mathbf {b} = {egin {bmatrix} b_ {1} b_ {2} vdots b_ {n} konec {bmatrix}}.}$

Pak A lze rozložit na a úhlopříčka komponent D, spodní trojúhelníková část L a horní trojúhelníkovou část U:

${displaystyle A = D + L + Uqquad {ext {where}} qquad D = {egin {bmatrix} a_ {11} & 0 & cdots & 0 0 & a_ {22} & cdots & 0 vdots & vdots & ddots & vdots 0 & 0 & cdots & a_ {nn} end {bmatrix }} {ext {and}} L + U = {egin {bmatrix} 0 & a_ {12} & cdots & a_ {1n} a_ {21} & 0 & cdots & a_ {2n} vdots & vdots & ddots & vdots a_ {n1} & a_ {n2} & cdots & 0end {bmatrix}}.}$

Roztok se poté získá iteračně pomocí

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

kde ${displaystyle mathbf {x} ^ {(k)}}$ je kth aproximace nebo iterace ${displaystyle mathbf {x}}$ a ${displaystyle mathbf {x} ^ {(k + 1)}}$ je další nebo k + 1 iterace ${displaystyle mathbf {x}}$. Elementární vzorec je tedy:

${displaystyle x_ {i} ^ {(k + 1)} = {frac {1} {a_ {ii}}} vlevo (b_ {i} -sum _ {jeq i} a_ {ij} x_ {j} ^ { (k)} ight), quad i = 1,2, ldots, n.}$

Výpočet ${displaystyle x_ {i} ^ {(k + 1)}}$ vyžaduje každý prvek v X^(k) kromě sebe. Na rozdíl od Gauss – Seidelova metoda, nemůžeme přepsat ${displaystyle x_ {i} ^ {(k)}}$ s ${displaystyle x_ {i} ^ {(k + 1)}}$, protože tato hodnota bude potřebná po zbytek výpočtu. Minimální velikost úložiště jsou dva vektory velikosti n.

Algoritmus

Vstup: počáteční odhad ${displaystyle x ^ {(0)}}$ k řešení„(diagonální dominantní) matice ${displaystyle A}$, pravý vektor ${displaystyle b}$, kritérium konvergenceVýstup: po dosažení konvergenceKomentáře: pseudokód založený na výše uvedeném vzorci založeném na prvcích${displaystyle k = 0}$zatímco nebylo dosaženo konvergence dělat    pro i: = 1 krok do n dělat        ${displaystyle sigma = 0}$        pro j: = 1 krok do n dělat            -li j ≠ i pak                ${displaystyle sigma = sigma + a_ {ij} x_ {j} ^ {(k)}}$            konec        konec        ${displaystyle x_ {i} ^ {(k + 1)} = {{frac {1} {a_ {ii}}} vlevo ({b_ {i} -sigma} ight)}}$    konec    ${displaystyle k = k + 1}$konec

Konvergence

Standardní podmínkou konvergence (pro jakoukoli iterativní metodu) je, když spektrální poloměr iterační matice je menší než 1:

${displaystyle ho (D ^ {- 1} (L + U)) <1.}$

Postačující (ale není nutnou) podmínkou pro konvergování metody je matice A je přísně nebo neredukovatelně diagonálně dominantní. Přísná dominance úhlopříčky řádků znamená, že pro každý řádek je absolutní hodnota úhlopříčného členu větší než součet absolutních hodnot ostatních členů:

${displaystyle left | a_ {ii} ight |> součet _ {jeq i} {left | a_ {ij} ight |}.}$

Jacobiho metoda někdy konverguje, i když tyto podmínky nejsou splněny.

Všimněte si, že Jacobiho metoda nekonverguje pro každou symetrii pozitivně definitivní matice.Například

${displaystyle A = {egin {pmatrix} 29 & 2 & 1 2 & 6 & 1 1 & 1 & {frac {1} {5}} end {pmatrix}} quad Rightarrow quad D ^ {- 1} (L + U) = {egin {pmatrix} 0 & { frac {2} {29}} & {frac {1} {29}} {frac {1} {3}} & 0 & {frac {1} {6}} 5 & 5 & 0end {pmatrix}} quad Rightarrow quad ho (D ^ {- 1} (L + U)) přibližně 1,0661 ,.}$

Příklady

Příklad 1

Lineární systém formuláře ${displaystyle Ax = b}$ s počátečním odhadem ${displaystyle x ^ {(0)}}$ je dána

${displaystyle A = {egin {bmatrix} 2 & 1 5 & 7 end {bmatrix}}, b = {egin {bmatrix} 11 13 end {bmatrix}} quad {ext {and}} quad x ^ {(0)} = {egin {bmatrix} 1 1 end {bmatrix}}.}$

Používáme rovnici ${displaystyle x ^ {(k + 1)} = D ^ {- 1} (b- (L + U) x ^ {(k)})}$, popsáno výše, odhadnout ${displaystyle x}$. Nejprve přepíšeme rovnici do pohodlnější formy ${displaystyle D ^ {- 1} (b- (L + U) x ^ {(k)}) = Tx ^ {(k)} + C}$, kde ${displaystyle T = -D ^ {- 1} (L + U)}$ a ${displaystyle C = D ^ {- 1} b}$. Ze známých hodnot

${displaystyle D ^ {- 1} = {egin {bmatrix} 1/2 & 0 0 & 1/7 end {bmatrix}}, L = {egin {bmatrix} 0 & 0 5 & 0 end {bmatrix}} quad {ext {and} } quad U = {egin {bmatrix} 0 & 1 0 & 0 end {bmatrix}}.}$

určujeme ${displaystyle T = -D ^ {- 1} (L + U)}$ tak jako

${displaystyle T = {egin {bmatrix} 1/2 & 0 0 & 1/7 end {bmatrix}} vlevo {{egin {bmatrix} 0 & 0 -5 & 0 end {bmatrix}} + {egin {bmatrix} 0 & -1 0 a 0 end {bmatrix}} ight} = {egin {bmatrix} 0 & -1 / 2 -5 / 7 & 0 end {bmatrix}}.}$

Dále, ${displaystyle C}$ je nalezen jako

${displaystyle C = {egin {bmatrix} 1/2 & 0 0 & 1/7 end {bmatrix}} {egin {bmatrix} 11 13 end {bmatrix}} = {egin {bmatrix} 11/2 13/7 konec {bmatrix}}.}$

S ${displaystyle T}$ a ${displaystyle C}$ vypočítané, odhadujeme ${displaystyle x}$ tak jako ${displaystyle x ^ {(1)} = Tx ^ {(0)} + C}$:

${displaystyle x ^ {(1)} = {egin {bmatrix} 0 & -1 / 2 -5 / 7 & 0 end {bmatrix}} {egin {bmatrix} 1 1 end {bmatrix}} + {egin {bmatrix } 11/2 13/7 end {bmatrix}} = {egin {bmatrix} 5,0 8/7 end {bmatrix}} přibližně {egin {bmatrix} 5 1,143 end {bmatrix}}.}$

Další výnosy z iterace

${displaystyle x ^ {(2)} = {egin {bmatrix} 0 & -1 / 2 -5 / 7 & 0 end {bmatrix}} {egin {bmatrix} 5.0 8/7 end {bmatrix}} + {egin {bmatrix} 11/2 13/7 end {bmatrix}} = {egin {bmatrix} 69/14 -12/7 end {bmatrix}} přibližně {egin {bmatrix} 4,929 -1,714 end {bmatrix }}.}$

Tento proces se opakuje až do konvergence (tj. Do ${displaystyle | Ax ^ {(n)} - b |}$ je malý). Řešení po 25 iteracích je

${displaystyle x = {egin {bmatrix} 7.111 -3.222end {bmatrix}}.}$

Příklad 2

Předpokládejme, že máme následující lineární systém:

${displaystyle {egin {aligned} 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 {zarovnáno}}}$

Pokud se rozhodneme (0, 0, 0, 0) jako počáteční aproximace je první přibližné řešení dáno vztahem

${displaystyle {egin {aligned} x_ {1} & = (6 + 0- (2 * 0)) / 10 = 0,6, x_ {2} & = (25 + 0 + 0- (3 * 0)) / 11 = 25/11 = 2,2727, x_ {3} & = (- 11- (2 * 0) + 0 + 0) /10=-1,1,x_ {4} & = (15- (3 * 0) +0) /8=1,875.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 pěti iteracích.

${displaystyle x_ {1}}$	${displaystyle x_ {2}}$	${displaystyle x_ {3}}$	${displaystyle x_ {4}}$
0.6	2.27272	-1.1	1.875
1.04727	1.7159	-0.80522	0.88522
0.93263	2.05330	-1.0493	1.13088
1.01519	1.95369	-0.9681	0.97384
0.98899	2.0114	-1.0102	1.02135

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

Příklad 3 pomocí Pythonu a NumPy

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

def jacobi(A, b, x_init, epsilon=1e-10, max_iterations=500):    D = np.diag(np.diag(A))    LU = A - D    X = x_init    D_inv = np.diag(1 / np.diag(D))    pro i v rozsah(max_iterations):        x_new = np.tečka(D_inv, b - np.tečka(LU, X))        -li np.linalg.norma(x_new - X) < epsilon:            vrátit se x_new        X = x_new    vrátit se X# problémová dataA = np.pole([    [5, 2, 1, 1],    [2, 6, 2, 1],    [1, 2, 7, 1],    [1, 1, 2, 8]])b = np.pole([29, 31, 26, 19])# můžete zvolit libovolný počáteční vektorx_init = np.nuly(len(b))X = jacobi(A, b, x_init)tisk("X:", X)tisk("vypočteno b:", np.tečka(A, X))tisk("real b:", b)

Produkuje výstup:

x: [3,99275362 2,95410628 2,16183575 0,96618357] vypočítané b: [29. 31. 26. 19.] skutečné b: [29 31 26 19]

Vážená Jacobiho metoda

Vážená Jacobiho iterace používá parametr ${displaystyle omega}$ vypočítat iteraci jako

${displaystyle mathbf {x} ^ {(k + 1)} = omega D ^ {- 1} (mathbf {b} - (L + U) mathbf {x} ^ {(k)}) + vlevo (1-omega ight) mathbf {x} ^ {(k)}}$

s ${displaystyle omega = 2/3}$ být obvyklou volbou.^[1]Ze vztahu ${displaystyle L + U = A-D}$, lze to vyjádřit také jako

${displaystyle mathbf {x} ^ {(k + 1)} = omega D ^ {- 1} mathbf {b} + vlevo (I-omega D ^ {- 1} Aight) mathbf {x} ^ {(k)} }$.

Konvergence v symetrickém pozitivním určitém případě

V případě, že systémová matice ${displaystyle A}$ je symetrický pozitivní-definitivní typ jeden může ukázat konvergenci.

Nechat ${displaystyle C = C_ {omega} = I-omega D ^ {- 1} A}$ iterační matice. Pak je zaručena konvergence

${displaystyle ho (C_ {omega}) <1quad Longleftrightarrow quad 0$

kde ${displaystyle lambda _ {ext {max}}}$ je maximální vlastní číslo.

Spektrální poloměr lze minimalizovat pro konkrétní výběr ${displaystyle omega = omega _ {ext {opt}}}$ jak následuje

${displaystyle min _ {omega} ho (C_ {omega}) = ho (C_ {omega _ {ext {opt}}}) = 1- {frac {2} {kappa (D ^ {- 1} A) +1 }} quad {ext {for}} quad omega _ {ext {opt}}: = {frac {2} {lambda _ {ext {min}} (D ^ {- 1} A) + lambda _ {ext {max }} (D ^ {- 1} A)}} ,,}$

kde ${displaystyle kappa}$ je číslo podmínky matice.

Viz také

• Gauss – Seidelova metoda
• Postupná nadměrná relaxace
• Iterační metoda § Lineární systémy
• Gaussian Belief Propagation
• Rozdělení matice

Reference

externí odkazy

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

• Černá, Noel; Moore, Shirley & Weisstein, Eric W. „Jacobiho metoda“. MathWorld.
• Jacobi Method z www.math-linux.com