Zobrazení zdroje stránky Gramova–Schmidtova ortogonalizace

{{Upravit}}
'''Gramův-Schmidtův proces''' neboli '''Gramova-Schmidtova ortogonalizace''' (nesprávně<ref>{{Citace elektronické monografie
 | titul = Internetová jazyková příručka
 | url = http://prirucka.ujc.cas.cz/?id=400#nadpis1
 | datum vydání = 
 | datum aktualizace = 
 | datum přístupu = 2011-10-31
 | vydavatel = Ústav pro jazyk český Akademie věd ČR
}}</ref> ''Gram-Schmidtova ortogonalizace'') je metoda, která v daném [[unitární prostor|unitárním prostoru]] (neboli [[vektorový prostor|vektorovém prostoru]] se [[skalární součin|skalárním součinem]]) umožňuje pro zadanou [[konečná množina|konečnou množinu]] [[vektor]]ů nalézt [[Ortonormální báze|ortonormální bázi]] [[podprostor]]u jimi generovaného.

== Algoritmus ==
Uvažujme pro jednoduchost <math>\mathbb{R}^m</math> reálný lineární vektorový prostor sloupcových vektorů o <math>m</math> složkách (se standardním skalárním součinem). Nechť <math>a_1,\ldots,a_n\in\mathbb{R}^m</math> jsou, opět pro jednoduchost, lineárně nezávislé, tedy <math>n\leq m</math>. Úkolem je nalézt ortonormální bázi <math>q_1,\ldots,q_n</math> <math>n</math>-rozměrného podprostoru <math>\mathbb{R}^m</math>, který vektory <math>a_i</math> generují; má tedy platit

:<math>\mathrm{span}\{a_1,\ldots,a_n\}=\mathrm{span}\{q_1,\ldots,q_n\},\qquad \langle q_k,q_j\rangle=q_j^Tq_k = \delta_{k,j}</math>

kde <math>\mathrm{span}</math> značí [[lineární obal]] množiny v závorce.

Algoritmus danou sadu vektorů prochází postupně přičemž v každém kroku vygeneruje nový vektor hledané báze. Omezíme-li se pouze na první vektor, a protože požadujeme aby <math>\|q_1\|=1</math>, musí platit

:<math>a_1 = q_1r_{1,1},\qquad \text{kde}\quad r_{1,1}=\|a_1\|_2,</math>

a dostáváme vztah pro výpočet prvního vektoru ortonormální báze <math>q_1=a_1/\|a_1\|_2</math>. Protože <math>a_2</math> je lineárně nezávislý na <math>a_1</math> a tedy i na <math>q_1</math>, můžeme ho vyjádřit jako

:<math>a_2 = q_1r_{1,2} + q_2r_{2,2},</math>

kde <math>q_2</math> je nějaký nový vektor takový, že <math>q_1^Tq_2=0,\;\|q_2\|_2=1</math>. Po pronásobení předchozího vztahu <math>q_1^T</math> zleva,

:<math>q_1^Ta_2 = q_1^Tq_1r_{1,2} + q_1^Tq_2r_{2,2}=r_{1,2}</math>

(připomeňme, že <math>q_1^Tq_1=\|q_1\|_2^2=1</math>), dostaneme vztah pro výpočet <math>r_{1,2}</math> (ortogonalizační koeficient; tj. velikost projekce <math>a_2</math> do směru <math>q_1</math>). Protože známe <math>a_2,\,q_1,\,r_{1,2}</math> dostáváme

:<math>a_2-q_1r_{1,2} = q_2r_{2,2}, \qquad\text{kde}\quad r_{2,2}=\|a_2-q_1r_{1,2}\|_2</math>

je norma zbytku vektoru <math>a_2</math> po ortogonalizaci proti <math>q_1</math>. Všimněme si, že po dosazení za <math>r_{1,2}</math> dostáváme

:<math>a_2-q_1q_1^Ta_2 = (I_m-q_1q_1^T)a_2 = q_2r_{2,2},\qquad\text{kde matice}\quad (I_m-q_1q_1^T)</math>

není nic jiného, než [[ortogonální projektor]] do [[ortogonální doplněk|ortogonálního doplňku]] <math>\mathrm{span}\{q_1\}^\perp</math> lineárního obalu vektoru <math>q_1</math> v <math>\mathbb{R}^m</math>.

Tento postup lze zřejmě opakovat do vyčerpání všech vektorů <math>a_k</math>.

Algoritmicky zapsáno:

:00: '''vstup''' <math>a_1,\ldots,a_n</math>
:01: <math>r_{1,1} := \|a_1\|_2</math>
:02: <math>q_1 := a_1/r_{1,1}</math>
:03: '''for''' <math>k := 2,\ldots,n</math>
:04: &nbsp;&nbsp;&nbsp; <math>p := a_k</math>
:05: &nbsp;&nbsp;&nbsp; '''for''' <math>j := 1,\ldots,k-1</math>
:06: &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; <math>r_{j,k} := q_j^Tp = \langle p,q_j\rangle</math>
:07: &nbsp;&nbsp;&nbsp; '''end'''
:08: &nbsp;&nbsp;&nbsp; '''for''' <math>j := 1,\ldots,k-1</math>
:09: &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; <math>p := p - q_jr_{j,k}</math>
:10: &nbsp;&nbsp;&nbsp; '''end'''
:11: &nbsp;&nbsp;&nbsp; <math>r_{k,k} := \|p\|_2</math>
:12: &nbsp;&nbsp;&nbsp; <math>q_k := p/r_{k,k}</math>
:13: '''end'''

Tato varianta algoritmu se nazývá '''klasický Gramův-Schmidtův algoritmus (CGS)''' a je novější než varianta původní, dnes zvaná '''modifikovaný Gramův-Schmidtův algoritmus (MGS)'''. MGS získáme z výše popsaného CGS prostým vypuštěním řádků 07 a 08, tedy, spojením obou vnitřních cyklů.

== Varianty algoritmu a jejich chování ==
Oba algoritmy CGS i MGS jsou matematicky ekvivalentní, jejich reálné implementace mají výrazně odlišné chování.

CGS algoritmus je výrazně paralelní. Výpočet první vnitřní smyčky (tj. výpočet koeficientů <math>r_{j,k}</math>) lze provádět nezávisle pro jednotlivá <math>j</math>; tedy, jednotlivá <math>r_{j,k}</math> mohu počítat na různých procesorech, jejich výpočet se neovlivňuje, nezávisí na sobe a může probíhat paralelně. Zatímco MGS je z tohoto pohledu sekvenční.

Na druhou stranu výpočet pomocí MGS je numericky výrazně stabilnější než výpočet pomocí CGS, kde může, vlivem zaokrouhlovacích chyb, dojít k úplné ztrátě ortogonality mezi vektory <math>q_1,\ldots,q_n</math>.

Označíme-li <math>Q_{k}=[q_1,\ldots,q_{k}]\in\mathbb{R}^{m\times k},\;Q_{k}^TQ_{k}=I_{k}</math>, lze vztah pro ortogonalizaci vektoru <math>a_{k+1}</math> psát pomocí projektorů dvěma matematicky ekvivalentními způsoby

:<math>p := (I_m-Q_{k}Q_{k}^T)a_{k+1}=(I_m-q_kq_k^T)\ldots(I_m-q_2q_2^T)(I_m-q_1q_1^T)a_{k+1}.</math>

První projekce odpovídá výpočtu pomocí CGS, druhá postupná výpočtu pomocí MGS. Je zřejmé že CGS ortogonalizace (projekce) se počítá paralelně do všech směrů najednou, kdežto sekvenční ortogonalizace (projekce) MGS umožňuje v <math>j</math>-tém kroku částečně eliminovat chyby vzniklé zaokrouhlováním v předchozích krocích <math>(j-1),\ldots,2,1</math>.

Řešením v praxi používaným bývá tzv. '''klasický Gramův-Schmidtův algoritmus s iteračním zpřesněním (ICGS)''', který obsahuje dvě vnitřní smyčky jako CGS (je tedy paralelizovatelný), ale obě smyčky se provedou dvakrát (čímž se výrazně zlepší numerické vlastnosti, ztráta ortogonality mezi vektory <math>q_i</math> je pak dokonce menší než u MGS).

=== Ztráta ortogonality ===
Nechť <math>\hat{Q}_n</math> je matice vektorů spočtených pomocí některé varianty Gramova-Schmidtova algoritmu v počítači se standardní konečnou aritmetikou s plovoucí řádovou čárkou, tj. <math>\hat{Q}_n=Q_n+E_n\approx Q_n</math> a <math>\hat{Q}_n^T\hat{Q}_n\approx I_n</math>. Veličina

:<math>\|\hat{Q}_n^T\hat{Q}_n-I_n\|_2</math>

se nazývá '''ztráta ortogonality''' a je jednou z klíčových veličin sloužících k posouzení kvality spočtené ortonormální báze.

Uvažujme tzv. '''Lauchliho matici'''<ref>{{Citace elektronické monografie
 | titul = Lauchli matrix
 | url = http://math.nist.gov/MatrixMarket/deli/Lauchli/information.html
 | datum vydání = 
 | datum aktualizace = 
 | datum přístupu = 2011-11-03
 | vydavatel = Matrix Market
}}</ref>

:<math>A=\left[\begin{array}{ccc}1&\ldots&1\\\rho&&0\\&\ddots&\\0&&\rho\end{array}\right]\in\mathbb{R}^{(n+1)\times n}, \qquad n=20,\; \rho=10^{-7},\; \kappa_2(A)\approx 4.47\times10^7,</math>

kde <math>\kappa_2(A)</math> je [[Podmíněnost matice|podmíněnost]] matice <math>A</math>.
Uvažujeme-li standardní aritmetiku se [[Strojová přesnost|strojovou přesností]] <math>\epsilon_M\approx2.22\times10^{-16}</math> (double), pak ztráta ortogonality odpovídající jednotlivým výše zmíněným algoritmům aplikovaným na danou Lauchliho matici, je ve druhém sloupci následující tabulky. Ve třetím sloupci je obecný vztah platný pro libovolnou matici <math>A</math>

{| class="wikitable"
|+ 
! Algoritmus
! Ztráta ortogonality (Lauchliho matice)
! Ztráta ortogonality (obecně)
|-
|CGS
|<math>2.2\times10^{-2}</math>
|<math>\kappa_2^2(A)\epsilon_M</math>
|-
|MGS
|<math>2.2\times10^{-9}</math>
|<math>\kappa_2(A)\epsilon_M</math>
|-
|ICGS
|<math>2.4\times10^{-16}</math>
|<math>\epsilon_M</math>
|-
|}

== Vztah Gramova-Schmidtova algoritmu a QR rozkladu ==
Srovnáním sloupcových vektorů <math>a_k,\;q_j</math> a koeficientů <math>r_{j,k}</math> do matic,

:<math>A=[a_1,\ldots,a_n],\quad Q=[q_1,\ldots,q_n]\in\mathbb{R}^{m\times n},\quad R=\left[\begin{array}{cccc}r_{1,1}&r_{1,2}&\ldots&r_{1,n} \\ 0&r_{2,2}&\ldots&r_{2,n}\\\vdots&\ddots&\ddots&\vdots\\0&\ldots&0&r_{n,n}\end{array}\right]\in\mathbb{R}^{n\times n},</math>

kde <math>Q^TQ=I_n</math> a <math>R</math> je čtvercová regulární matice dostáváme

:<math>A=QR</math>

tedy [[QR rozklad]] matice <math>A</math> (pro ověření stačí porovnat <math>k</math>-tý sloupec rovnosti, tedy <math>a_k</math> s <math>k</math>-tým sloupcem součinu <math>QR</math>).

== Ortogonální polynomy ==
Gramův-Schmidtův algoritmus lze aplikovat na prvky libovolného prostoru se skalárním součinem. Uvažujeme například prostor polynomů <math>\mathcal{P}(a,b)</math> se skalárním součinem

:<math>\langle p(x),q(x)\rangle_w = \int_a^b p(x)q(x)w(x)dx</math>,

kde <math>w(x)</math> je nějaká [[váhová funkce]]. Aplikací Gramova-Schmidtova algoritmu na sadu vektorů (polynomů) <math>1,x,x^2,\ldots,x^{n-1}</math> (v tomto pořadí) dostáváme, pro vhodně volené <math>a,\,b</math> a váhovou funkci <math>w(x)</math> libovolnou sadu [[Ortogonální polynomy|ortogonálních polynomů]] (normalizovanou vzhledem k normě indukované daným skalárním součinem).

Pro <math>a=-1,\,b=1,\,w(x)=1</math> Gramův-Schmidtův algoritmus generuje [[Legenderovy polynomy]],
pro <math>a=-1,\,b=1,\,w(x)=(1-x^2)^{-1/2}</math> dostaneme [[Čebyševovy polynomy]] prvního druhu,
atd.

== Reference ==
<references/>

== Literatura ==
* [[Gene H. Golub|Gene Howard Golub]], Charles F. Van Loan: ''Matrix Computations'', Johns Hopkins University Press, 1996 (3rd Ed.). (Zejména kapitoly 5.2.7 CGS, 5.2.8 MGS a 5.2.9 Work and Accuracy.)
* J. Duintjer Tebbens, I. Hnětynková, M. Plešinger, [[Zdeněk Strakoš|Z. Strakoš]], P. Tichý: ''Analýza metod pro maticové výpočty, základní metody.'' Matfyzpress 2012. {{ISBN|978-80-7378-201-6}}. (Kapitola 3, Ortogonální transformace a QR rozklady, str. 53-88.)
{{Autoritní data}}

[[Kategorie:Lineární algebra]]