Postupná nadměrná relaxace - Successive over-relaxation

v numerická lineární algebra, metoda postupná nadměrná relaxace (SOR) je variantou Gauss – Seidelova metoda pro řešení a lineární soustava rovnic, což vede k rychlejší konvergenci. Podobnou metodu lze použít pro jakékoli pomalu konvergující iterativní proces.

To bylo navrženo současně David M. Young Jr. a tím Stanley P. Frankel v roce 1950 za účelem automatického řešení lineárních systémů na digitálních počítačích. Před prací Younga a Frankela byly použity metody přílišné relaxace. Příkladem je metoda Lewis Fry Richardson a metody vyvinuté R. V. Southwell. Tyto metody však byly navrženy pro výpočet pomocí lidské kalkulačky, vyžadující určité odborné znalosti k zajištění konvergence k řešení, které je činí nepoužitelnými pro programování na digitálních počítačích. Tyto aspekty jsou diskutovány v práci Davida M. Younga Jr.^[1]

Formulace

Vzhledem k čtvercovému systému n lineární rovnice s neznámým X:

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

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, a přísně spodní a horní trojúhelníkový komponenty L a U:

${displaystyle A = D + L + U,}$

kde

${displaystyle D = {egin {bmatrix} a_ {11} & 0 & cdots & 0 0 & a_ {22} & cdots & 0 vdots & vdots & ddots & vdots 0 & 0 & cdots & a_ {nn} end {bmatrix}}, quad L = {egin {bmatrix} 0 & 0 & cdots & 0 a_ {21} & 0 & cdots & 0 vdots & vdots & ddots & vdots a_ {n1} & a_ {n2} & cdots & 0end {bmatrix}}, quad U = {egin {bmatrix} 0 & a_ {12} & cdots & a_ {1n} 0 & 0 & cdots & a_ {2n} vdots & vdots & ddots & vdots 0 & 0 & cdots & 0end {bmatrix}}.}$

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

${displaystyle (D + omega L) mathbf {x} = omega mathbf {b} - [omega U + (omega -1) D] mathbf {x}}$

pro konstantu ω > 1, nazvaný relaxační faktor.

Metoda postupné nadměrné relaxace je iterační technika který řeší levou stranu tohoto výrazu X, s použitím předchozí hodnoty pro X na pravé straně. Analyticky to lze napsat jako:

${displaystyle mathbf {x} ^ {(k + 1)} = (D + omega L) ^ {- 1} {ig (} omega mathbf {b} - [omega U + (omega -1) D] mathbf {x} ^ {(k)} {ig)} = L_ {w} mathbf {x} ^ {(k)} + mathbf {c},}$

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}}$Využitím trojúhelníkové formy (D+ωL), prvky X^(k+1) lze vypočítat postupně pomocí dopředné střídání:

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

Konvergence

Spektrální poloměr ${displaystyle ho (C_ {omega})}$ iterační matice pro metodu SOR ${displaystyle C_ {omega}}$. Graf ukazuje závislost na spektrálním poloměru matice Jacobiho iterace ${displaystyle mu: = ho (C_ {ext {Jac}})}$.

Volba relaxačního faktoru ω není nutně snadné a závisí na vlastnostech matice koeficientu. V roce 1947 Ostrowski dokázal, že pokud ${displaystyle A}$ je symetrický a pozitivní-definitivní pak ${displaystyle ho (L_ {omega}) <1}$ pro ${displaystyle 0$. Následuje tedy konvergence iteračního procesu, ale obecně nás zajímá spíše rychlejší konvergence než jen konvergence.

Míra konvergence

Konvergenční poměr pro metodu SOR lze analyticky odvodit. Je třeba předpokládat následující

• relaxační parametr je vhodný: ${displaystyle omega in (0,2)}$
• Jacobiho iterační matice ${displaystyle C_ {ext {Jac}}: = I-D ^ {- 1} A}$ má pouze skutečná vlastní čísla
• Jacobiho metoda je konvergentní: ${displaystyle mu: = ho (C_ {ext {Jac}}) <1}$
• existuje jedinečné řešení: ${displaystyle det Aeq 0}$.

Pak lze míru konvergence vyjádřit jako^[2]

${displaystyle ho (C_ {omega}) = {egin {cases} {frac {1} {4}} left (omega mu + {sqrt {omega ^ {2} mu ^ {2} -4 (omega -1)} } ight) ^ {2} ,, & 0$

kde je optimální relaxační parametr dán vztahem

${displaystyle omega _ {ext {opt}}: = 1 + left ({frac {mu} {1+ {sqrt {1-mu ^ {2}}}}} ight) ^ {2} ,.}$

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é indexování je vynecháno. Algoritmus je následující:

Vstupy: A, b, ωVý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 (1-omega) phi _ {i} + {frac {omega} {a_ {ii}}} (b_ {i} -sigma)}$    konec (i-loop) zkontrolujte, zda je dosaženo konvergencekonec (opakovat)

Poznámka

${displaystyle (1-omega) phi _ {i} + {frac {omega} {a_ {ii}}} (b_ {i} -sigma)}$ lze také psát ${displaystyle phi _ {i} + zbývající omega ({frac {b_ {i} -sigma} {a_ {ii}}} - phi _ {i} ight)}$, čímž šetří jednu násobení v každé iteraci vnější pro-smyčka.

Příklad

Představujeme lineární systém

${displaystyle {egin {aligned} 4x_ {1} -x_ {2} -6x_ {3} + 0x_ {4} = 2, - 5x_ {1} -4x_ {2} + 10x_ {3} + 8x_ {4} = 21, 0x_ {1} + 9x_ {2} + 4x_ {3} -2x_ {4} = - 12, 1x_ {1} + 0x_ {2} -7x_ {3} + 5x_ {4} = - 6 .end {zarovnáno}}}$

Pro řešení rovnic zvolíme relaxační faktor ${displaystyle omega = 0,5}$ a počáteční vektor odhadu ${displaystyle phi = (0,0,0)}$. Podle postupného algoritmu nadrelaxace je získána následující tabulka představující příkladnou iteraci s aproximacemi, která v ideálním případě, ale ne nutně, najde přesné řešení, (3, −2, 2, 1), v 38 krocích.

Opakování	${displaystyle x_ {1}}$	${displaystyle x_ {2}}$	${displaystyle x_ {3}}$	${displaystyle x_ {4}}$
1	0.25	-2.78125	1.6289062	0.5152344
2	1.2490234	-2.2448974	1.9687712	0.9108547
3	2.070478	-1.6696789	1.5904881	0.76172125
...	...	...	...	...
37	2.9999998	-2.0	2.0	1.0
38	3.0	-2.0	2.0	1.0

Níže je uvedena jednoduchá implementace algoritmu v programu Common Lisp. V obecném případě si dejte pozor na jeho sklon k přetečení s plovoucí desetinnou čárkou.

(defparametr + MAXIMÁLNÍ POČET ITERACÍ + 100  "Počet iterací, po jejichž překročení by algoritmus měl být."   ukončit svou činnost bez ohledu na její aktuální řešení. “)(defun get-chyby (vypočítané řešení přesné řešení)  "Pro každou složku vektoru COMPUTED-SOLUTION načte její   chyba s ohledem na očekávané PŘESNÉ ŘEŠENÍ, vrácení a   vektor chybových hodnot. "  (mapa 'vektor #'- přesné řešení vypočítané řešení))(defun je-konvergentní (chyby &klíč (tolerance chyb 0.001))  "Zkontroluje, zda je dosaženo konvergence s ohledem na   ERRORS vektor, který registruje nesrovnalost mezi vypočítanou   a přesný vektor řešení.   ---   Konvergence je splněna, pokud je každá složka absolutní chyby   menší nebo rovno CHYBOVÉ TOLERANCI. “  (fllet ((chyba je přijatelná (chyba)          (<= (břišní svaly chyba) tolerance chyb)))    (každý #'chyba je přijatelná chyby)))(defun make-zero-vector (velikost)  "Vytvoří a vrátí vektor SIZE se všemi prvky nastavenými na 0."  (make-pole velikost : počáteční prvek 0.0 : typ prvku 'číslo))(defun sor (A b omega            &klíč (phi (make-zero-vector (délka b)))                 (kontrola konvergence                   #'(lambda (opakování phi)                       (prohlásit (ignorovat phi))                       (>= opakování + MAXIMÁLNÍ POČET ITERACÍ +))))  „Implementuje postupnou metodu nadměrné relaxace (SOR), na kterou se aplikuje   lineární rovnice definované maticí A a pravou stranou   vektor B, využívající relaxační faktor OMEGA, vrací   vypočítaný vektor řešení.   ---   První krok algoritmu, volba počátečního odhadu PHI, je   představuje volitelný parametr klíčového slova PHI, který je výchozí   na nulový vektor stejné struktury jako B. Je-li zadán, je tento   vektor bude destruktivně upraven. V každém případě vektor PHI   představuje výslednou hodnotu funkce.   ---   Koncová podmínka je implementována pomocí CONVERGENCE-CHECK,   volitelný predikát     lambda (iterace phi) => zobecněný-boolean   který vrací T, což znamená okamžité ukončení, po dosažení   konvergence neboli NIL, signalizující kontinuální provoz, jinak. “  (nechat ((n (dimenze pole A 0)))    (smyčka pro opakování z 1 podle 1 dělat      (smyčka pro i z 0 níže n podle 1 dělat        (nechat ((rho 0))          (smyčka pro j z 0 níže n podle 1 dělat            (když (/= j i)              (nechat ((a [ij]  (aref A i j))                    (phi [j] (aref phi j)))                (příloha rho (* a [ij] phi [j])))))          (setf (aref phi i)                (+ (* (- 1 omega)                      (aref phi i))                   (* (/ omega (aref A i i))                      (- (aref b i) rho))))))      (formát T „~ & ~ d. solution = ~ a“ opakování phi)      ;; Zkontrolujte, zda je dosaženo konvergence.      (když (funcall kontrola konvergence opakování phi)        (vrátit se))))  phi);; Vyvolejte funkci s příkladnými parametry.(nechat ((A              (make-pole (seznam 4 4)                        : počáteční obsah                        '((  4  -1  -6   0 )                          ( -5  -4  10   8 )                          (  0   9   4  -2 )                          (  1   0  -7   5 ))))      (b              (vektor 2 21 -12 -6))      (omega          0.5)      (přesné řešení (vektor 3 -2 2 1)))  (sor    A b omega    : kontrola konvergence    #'(lambda (opakování phi)        (nechat ((chyby (get-chyby phi přesné řešení)))          (formát T "~ & ~ d. chyby = ~ a" opakování chyby)          (nebo (je-konvergentní chyby : tolerance chyb 0.0)              (>= opakování + MAXIMÁLNÍ POČET ITERACÍ +))))))

Jednoduchá implementace výše uvedeného pseudokódu v Pythonu.

import numpy tak jako npdef sor_solver(A, b, omega, initial_guess, konvergence_criteria):    """    Toto je implementace pseudokódu poskytovaného v článku Wikipedie.    Argumenty:        A: nxn numpy matrix.        b: n rozměrný numpy vektor.        omega: relaxační faktor.        initial_guess: Počáteční odhad řešení pro řešení, který má začít.        Convergence_criteria: Maximální odchylka přijatelná, pokud jde o současné řešení jako vhodné.    Vrací:        phi: vektor řešení dimenze n.    """    krok = 0    phi = initial_guess[:]    reziduální = np.linalg.norma(np.matmul(A, phi) - b)  # Počáteční zbytek    zatímco reziduální > konvergence_criteria:        pro i v rozsah(A.tvar[0]):            sigma = 0            pro j v rozsah(A.tvar[1]):                -li j != i:                    sigma += A[i, j] * phi[j]            phi[i] = (1 - omega) * phi[i] + (omega / A[i, i]) * (b[i] - sigma)        reziduální = np.linalg.norma(np.matmul(A, phi) - b)        krok += 1        tisk("Krok {} Reziduální: {: 10,6 g}".formát(krok, reziduální))    vrátit se phi# Příklad případu, který odráží případ uvedený v článku na Wikipediiresidual_convergence = 1e-8omega = 0.5  # Relaxační faktorA = np.pole([[4, -1, -6, 0],              [-5, -4, 10, 8],              [0, 9, 4, -2],              [1, 0, -7, 5]])b = np.pole([2, 21, -12, -6])initial_guess = np.nuly(4)phi = sor_solver(A, b, omega, initial_guess, residual_convergence)tisk(phi)

Symetrické postupné uvolnění

Verze pro symetrické matice A, ve kterém

${displaystyle U = L ^ {T} ,,}$

se označuje jako Symetrická postupná nadměrná relaxace, nebo (SSOR), ve kterém

${displaystyle P = left ({frac {D} {omega}} + Light) {frac {omega} {2-omega}} D ^ {- 1} left ({frac {D} {omega}} + Uight), }$

a iterační metoda je

${displaystyle mathbf {x} ^ {k + 1} = mathbf {x} ^ {k} -gamma ^ {k} P ^ {- 1} (Amathbf {x} ^ {k} -mathbf {b}), kgeq 0.}$

Metody SOR a SSOR jsou připisovány David M. Young Jr..

Další aplikace metody

Podobnou techniku ​​lze použít pro jakoukoli iterační metodu. Pokud původní iterace měla formu

${displaystyle x_ {n + 1} = f (x_ {n})}$

pak by se použila upravená verze

${displaystyle x_ {n + 1} ^ {mathrm {SOR}} = (1-omega) x_ {n} ^ {mathrm {SOR}} + omega f (x_ {n} ^ {mathrm {SOR}}).}$

Formulace uvedená výše, používaná k řešení soustav lineárních rovnic, však není zvláštním případem této formulace, pokud X je považován za úplný vektor. Pokud se místo toho použije tato formulace, bude vypadat rovnice pro výpočet dalšího vektoru

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

kde ${displaystyle L _ {*} = L + D}$. Hodnoty ${displaystyle omega> 1}$ se používají k urychlení konvergence procesu pomalé konvergence, zatímco hodnoty ${displaystyle omega <1}$ se často používají jako pomůcka pro vytvoření konvergence odlišného iteračního procesu nebo pro urychlení konvergence přestřelení proces.

Existují různé metody, které adaptivně nastavují parametr relaxace ${displaystyle omega}$ na základě pozorovaného chování konvergujícího procesu. Obvykle pomáhají dosáhnout u některých problémů superlineární konvergence, u jiných však selhávají.

Viz také

• Jacobiho metoda
• Gaussian Belief Propagation
• Rozdělení matice

Poznámky

Reference

• Tento článek včlení text z článku Metoda po sobě jdoucí relaxace _-_ SOR na CFD-Wiki to je pod GFDL licence.

• Abraham Berman, Robert J. Plemmons, Nezáporné matice v matematických vědách, 1994, SIAM. ISBN  0-89871-321-8.
• Black, Noel & Moore, Shirley. „Postupná metoda overrelaxation“. MathWorld.
• A. Hadjidimos, Postupná overrelaxace (SOR) a související metody, Journal of Computational and Applied Mathematics 123 (2000), 177–199.
• Yousef Saad, Iterační metody pro řídké lineární systémy, 1. vydání, PWS, 1996.
• Netlib kopie "Šablony pro řešení lineárních systémů", Barrett et al.
• Richard S.Varga 2002 Iterativní analýza matice, Druhé vydání. (vydání Prentice Hall z roku 1962), Springer-Verlag.
• David M. Young Jr. Iterativní řešení velkých lineárních systémů, Academic Press, 1971. (dotisk Dover, 2003)

externí odkazy

• Modul pro metodu SOR
• Tridiagonální lineární systémový řešič založený na SOR v C ++