¨ zel´ıto ˝ e ´s szimbolikus sza ´ m´ıta ´ sok Ko 3. gyakorlat
Gauss elimin´ aci´ o, LU felbont´ as K´esz´ıtette: Gelle Kitti
Csendes Tibor Somogyi Viktor London Andr´as De´ak G´abor jegyzetei alapj´an
1 EGYENLETRENDSZEREK
1.
Egyenletrendszerek
A numerikus elj´ar´asok egyik legfontosabb r´esz´et adj´ak az Ax = b alak´ u egyenletrendszerek megold´asai. Line´aris egyenletrendszer nek nevezz¨ uk a k¨ovetkez˝o alak´ u egyenletrendszereket: a11 x1 + a12 x2 + · · · + a1n xn = b1 a21 x1 + a22 x2 + · · · + a2n xn = b2 .. . an1 x1 + an2 x2 + · · · + ann xn = bn ahol aij ∈ C, teh´at az egy¨ utthat´ok komplex sz´amok is lehetnek. Egy egyszer˝ u p´elda line´aris egyenletrendszerre: 2x1 +4x2 +8x3 = 5 6x1 +1x2 +5x3 = 7 2x1 +3x2 +3x3 = 1 M´atrix alakban fel´ırva: x1 5 2 4 8 A = 6 1 5, b = 7, x = x2 x3 1 2 3 3
Az Ax = b alak´ u line´aris egyenletrendszereknek pontosan akkor van egyetlen megold´asuk, ha az A m´atrix nem szingul´aris (ekkor regul´aris, azaz det(A) 6= 0, A rangja n), ´es ekkor ez a megold´as x = A−1 b. Defin´ıci´ o 1.1. Az A m´atrix rangja a line´arisan f¨ uggetlen oszlopvektorainak sz´ama. Defin´ıci´ o 1.2. Vektorok egy halmaz´at line´arisan f¨ uggetlennek nevezz¨ uk, ha egyik¨ uk sem fejezhet˝o ki a t¨obbi vektor line´aris kombin´aci´ojak´ent. Defin´ıci´ o 1.3. A szingul´aris m´atrix egy olyan n´egyzetes m´atrix, amely nem invert´alhat´ o. Egy A n´egyzetes m´atrix akkor ´es csak akkor szingul´aris, ha a determin´ansa 0, vagyis det(A) = 0. Defin´ıci´ o 1.4. Egy n´egyzetes A m´atrix akkor invert´alhat´o, ha l´etezik olyan n × n-es B m´atrix, melyre igaz, hogy: AB = BA = In , ahol In az n × n-es egys´egm´atrixot jel¨oli. Ekkor B-t az A m´atrix inverz´enek nevezz¨ uk (A−1 ).
1
´ O ´ 2 GAUSS-ELIMINACI
Feladat 1. Oldjuk meg a fenti egyenletet a Matlab seg´ıts´eg´evel! Megold´ as: >> A = [2 4 8; 6 1 5; 2 3 3]; >> b = [5 7 1]'; >> x = inv(A)*b %ez egyenerteku az x=A\b utasitassal
A fenti m´odszer h´atr´anya, hogy nagyon m˝ uveletig´enyes, illetve numerikusan instabil, emiatt ezt nem fogjuk haszn´alni (b´ar ebben a feladatban mindent kisz´amol u ´gy ahogy kell, n´ezz¨ uk meg egy 1000 × 1000 m´eret˝ u m´atrixra).
2.
Gauss-elimin´ aci´ o
Olyan m´atrixot viszonylag egyszer˝ uen meg tudunk adni, amellyel balr´ol val´o beszorz´as azt eredm´enyezi, hogy egy a vektor k-adik alatti egy¨ utthat´oi elt˝ unnek: a1 a1 1 a2 a2 . . . .. .. . . 1 ak = ak Mk a = −ak+1 /ak ak+1 0 .. .. . .. . .. . . −an /ak 1 an 0 Itt az Mk m´atrixban minden ki nem ´ırt egy¨ utthat´o nulla. Az ak egy¨ utthat´ot gener´al´oelemnek h´ıvjuk, mag´at az Mk m´atrixot pedig elimin´aci´os m´atrix nak, vagy Gausstranszform´aci´os m´atrix nak. A Gauss-elimin´aci´o l´enyege, hogy egy Ax = b alak´ u line´aris egyenletrendszerb˝ol elimin´aci´os m´atrixok seg´ıts´eg´evel olyat alak´ıtunk ki, amely baloldal´an egy fels˝o h´aromsz¨ogm´atrix van, ´ıgy az egyenletrendszerbe val´o behelyettes´ıt´est meg tudjuk oldani, ´es ´ıgy megkapjuk mag´anak az egyenletrendszernek a megold´as´at is. P´ elda. Hajtsuk v´egre az adott egyenletrendszeren a Gauss-elimin´aci´ot! 2x1 +6x2 +3x3 = 2 x1 +5x2 +15x3 = 4 3x1 +10x2 +4x3 = 6
2
´ O ´ 2 GAUSS-ELIMINACI
Megold´ as. ´Irjuk fel m´atrixos alakban az egyenletrendszert! (Ax = b)
2 6 3 x1 2 1 5 15 x2 = 4 3 10 4 x3 6 Egyenletek megold´asakor a bevett gyakorlat, hogy mindk´et oldalon azonos m˝ uveleteket v´egz¨ unk. Jelen esetben Gauss-transzform´aci´os m´atrixokat kell keresn¨ unk, ´es ezekkel beszoroznunk az egyenlet mindk´et oldal´at balr´ol, mivel a c´elunk az A m´atrix fels˝o h´aromsz¨ogm´atrixsz´a alak´ıt´asa. R¨oviden teh´at a k¨ovetkez˝o ´atalak´ıt´asokat kell v´egrehajtanunk: M2 M1 Ax = M2 M1 b El˝osz¨or hajtsuk v´egre az A m´atrixon a Gauss-elimin´aci´o l´ep´eseit: M1 A = A1 M2 A1 = A2 Azaz: 2 1 0 0 2 6 3 −1/2 1 0 1 5 15 = 0 0 −3/2 0 1 3 10 4 1 0 0 2 6 3 2 1 0 0 2 27/2 = 0 0 0 −1/2 1 0 1 −1/2 0
6 3 2 27/2 1 −1/2 6 3 2 27/2 0 −29/4
Majd a b vektoron is, ugyanazokkal a transzform´aci´os m´atrixokkal: M1 b = b1 M2 b1 = b2 Azaz:
1 0 −1/2 1 −3/2 0 1 0 1 0 0 −1/2
2 0 2 0 4 = 3 6 3 1 0 2 2 0 3 = 3 1 3 3/2
3
´ 3 LU FELBONTAS
Az egyenlet¨ unket ´atrendezt¨ uk, most m´ar u ´jra fel´ırhatjuk az a´talak´ıtott form´at ´es leolvashatjuk a megold´asokat:
2 6 3 x1 2 0 2 27/2 x2 = 3 0 0 −29/4 x3 3/2 ⇓ x3 = − x2 =
3−
6 29
27 x 2 3
=
84 29
2 2 − 6x2 − 3x3 214 x1 = =− 2 29 Ellen˝ orz´ es.
84 6 214 +6 −3 =2 29 29 29 214 84 6 − + 5 − 15 = 4 29 29 29 214 84 6 −3 + 10 − 4 = 6 29 29 29 −2
3.
LU felbont´ as
Az el˝obbi elj´ar´as sor´an ´eszrevehett¨ uk, hogy egy als´o ´es egy fels˝o h´aromsz¨og-m´atrix szorzat´ara bontottuk az A m´atrixot. Ez fogja motiv´alni az LU felbont´ast, amit a k¨ovetkez˝okben defini´alunk. Defin´ıci´ o 3.1. Az LU felbont´as az eredeti A m´atrixot k´et h´aromsz¨og-m´atrix szorzat´ara bontja: A = LU , ahol −1 L = M −1 = (Mn−1 . . . M1 )−1 = M1−1 . . . Mn−1 = L1 . . . Ln−1
De m´egis hogy j¨ott ez ki? Az el˝oz˝o konstrukci´ob´ol tudjuk, hogy Mn−1 Mn−2 . . . M1 A m´atrix egy fels˝otriangul´aris m´atrix, ez lesz most az U . Azt is tudjuk, hogy b´armilyen regul´aris X m´atrixra A = X −1 XA igaz, mert az X −1 ´es az X ,,ki¨ uti egym´ast”, azaz egys´egm´atrix lesz bel˝ole. Teh´at ha X = Mn−1 ...M1 , akkor az el˝obbi alapj´an kapjuk, hogy A = (Mn−1 ..M1 )−1 (Mn−1 ..M1 )A, ebb˝ol az elej´et, az (Mn−1 ..M1 )−1 -et nevezz¨ uk L-nek, a marad´ek pedig az el˝obbi U m´atrix. Foglalkozzunk egy kicsit az L = (Mn−1 Mn−2 ...M1 )−1 sorozattal. Ugyan ´ıgy is ki lehetne sz´amolni, de invert´alni lassan tudunk m´atrixot. Ehelyett azt haszn´aljuk ki, hogy −1 (Mn−1 Mn−2 ...M1 )−1 = M1−1 M2−1 ...Mn−1 (m´atrixok szorzat´anak inverze mindig egyenl˝o a 4
´ 3 LU FELBONTAS
m´atrixok inverzeinek ford´ıtott sorrendben vett szorzat´aval), mert ezeket az Mi Gausstranszform´aci´os m´atrixokat k¨onny˝ u invert´alni, ´es ´ıgy az Li = Mi−1 m´atrixokat k¨onnyen kisz´amolhatjuk, csak ¨ossze kell majd szorozni ´es megvan az L m´atrixunk. Mivel az Li k als´otriangul´aris m´atrixok, ez´ert az eg´esz L is az lesz, mivel als´otriangul´aris m´atrixok ´ ha ez a felbont´as megvan, akkor az egyenletrendszer¨ szorzata als´otriangul´aris. Es unk LU x = b alakba ker¨ ul. P´ elda. Adjuk meg az el˝oz˝o p´elda A m´atrix´anak LU felbont´as´at a most tanult m´odszerekkel! Megold´ as. Az inverz m´atrixot ebben az esetben k¨onnyen megkapjuk, mind¨ossze annyit kell tenn¨ unk, hogy a nem f˝o´atl´obeli elemek el˝ojel´et az ellent´etesre v´altoztatjuk. Teh´at: 1 0 0 1 0 0 M1 = −1/2 1 0 =⇒ L1 = 1/2 1 0 −3/2 0 1 3/2 0 1 0 1 0 0 M2 = 0 1 0 =⇒ L2 = 0 1 0 1/2 0 −1/2 1
1 0 0 1
⇓
1 0 0 L = L1 · L2 = 1/2 1 0 3/2 1/2 1 2 6 3 U = M2 · M1 · A = 0 2 27/2 0 0 −29/4 1 0 0 2 6 3 2 6 3 LU = 1/2 1 0 0 2 27/2 = 1 5 15 = A 3/2 1/2 1 0 0 −29/4 3 10 4
3.1
Line´ aris egyenletrendszer megold´ asa LU felbont´ assal
Ahhoz, hogy eljussunk egy line´aris egyenletrendszer megold´as´ahoz, a k¨ovetkez˝o l´ep´eseket kell megtenn¨ unk: 1. Bontsuk fel az A m´atrixot LU alakra (⇒ LU x = b). 2. Oldjuk meg az Ly = b egyenletet y-ra. 3. Oldjuk meg az U x = y egyenletet x-re. 5
´ 3 LU FELBONTAS
P´ elda. Az eddigi p´eld´ank megold´asa a fenti LU felbont´asos egyenlettel. Megold´ as. Ly = b
1 0 0 y1 2 2 1/2 1 0 y2 = 4 =⇒ 3 6 3/2 3/2 1/2 1 y3 Ux = y
2 6 3 x1 2 −214/29 0 2 27/2 x2 = 3 =⇒ 84/29 0 0 −29/4 x3 3/2 −6/29
3.2
LU felbont´ as Matlabban
A Matlabban a k¨ovetkez˝o f¨ uggv´ennyel tudunk egy LU felbont´ast v´egrehajtani: >> [L,U] = lu(A) – egy fels˝ o triangul´aris ´es egy permut´alt als´o triangul´aris
m´atrixot eredm´enyez u ´gy, hogy A = LU . >> [L,U,P] = lu(A) – egy fels˝ o triangul´aris, egy als´o triangul´aris ´es egy permut´aci´os m´atrixot eredm´enyez u ´gy, hogy P A = LU . >> A = [2 6 3; 1 5 15; 3 10 4]; >> b = [2 4 6]'; >> [L,U] = lu(A); >> y = L\b; >> x = U\y x = −7.3793 2.8966 −0.2069
6
4 FELADATOK
4.
Feladatok 1. Legyen A egy 100 × 100 m´eret˝ u, 1-t˝ol 10-ig v´eletlen elemeket tartalmaz´o m´atrix, valamint legyen b egy 100×1 m´eret˝ u, 1-t˝ol 10-ig v´eletlen elemeket tartalmaz´o vektor. Oldjuk meg az Ax = b egyenletet Matlabban (LU felbont´as n´elk¨ ul, az inverzzel val´o szorz´ast haszn´alva)! 2. Hasonl´o m´odon gener´aljunk egy A ´es b m´atrixot, majd oldjuk meg azt a Matlab seg´ıts´eg´evel, az LU felbont´ast haszn´alva! 3. Hasonl´ıtsuk o¨ssze a k´et m´odszer fut´asidej´et! Teh´at: gener´alni kell egy A-t ´es egy b-t, majd v´egrehajtani a hagyom´anyos m´odszert ´es az LU felbont´ast is ugyan azon a gener´alt A ´es b m´atrixon. A fut´asi id˝oket a tic() ´es toc() f¨ uggv´enyh´ıv´asokkal ellen˝orizhetj¨ uk.
7