Gauss-Seidel迭代矩陣求法的思考

1、Gauss-Seide迭代矩陣求法的思考在迭代法收斂性的判別中，我們有充分條件：若迭代矩陣B的某種范數(shù)怛1=q<1,貝U迭代法x(k*)=Bx的+d,k=0,1，對(duì)任意的初始向量x(0)都收斂于方程組Ax=b的精確解x*。從這個(gè)條件中我們可以看出，想要知道迭代法是否收斂，就要知道迭代矩陣(當(dāng)然如果系數(shù)矩陣是正定的或嚴(yán)格對(duì)角占優(yōu)的，那就不用知道其迭代矩陣，因?yàn)檫@時(shí)它的Jacobi迭代和Gauss-Seidel迭代一定收斂)，Jacobi迭代矩陣為Bj=D,(L+U)=I-D/A,Gauss-Seidel迭代矩陣Bg=-(D+L),U，這兩個(gè)矩陣中都涉及到了矩陣的逆從上高等代數(shù)時(shí)學(xué)到矩陣的逆

2、開(kāi)始，就一直懼怕有關(guān)矩陣逆的題目，因?yàn)榍缶仃嘇的逆Aa1=A*=.一.*A,這就必須求出A的行列式A與A的伴隨矩陣A,對(duì)于求矩陣A的行列式，就是一個(gè)繁瑣的過(guò)程，計(jì)算量大且易出錯(cuò)，而這兒還不僅如此,這兒還要求出矩陣A的伴隨矩陣Ao如果矩陣a11a2a1na21a22a2naaa,則t3n1an2annA=A11A21.*A12A22A=AinA2na11An1aAn2,而具中的Aj=ai,1ai書(shū),1Annan1&,j1ma1,j1man9a1ai,j1-ai,nai1,j1-ai1,j1-ai1,naan,j1an,j1-ann因此求A*的計(jì)算量比求A的行列式的計(jì)算量還要大的多，所以A

3、，很難求。因此數(shù)學(xué)家便開(kāi)始尋找求A，的相對(duì)容易的方法，其中有一種初等變換的方法，即對(duì)(AE世行初等行變換，當(dāng)把A變成E時(shí)，E便變成了A,，此方法要簡(jiǎn)單的多，但在變換過(guò)程中要消耗大量空間。在用迭代法解線(xiàn)性方程組的方法中，都涉及到了一個(gè)矩陣的逆，而且其涉及到的還不僅僅是一個(gè)矩陣的逆那么簡(jiǎn)單，其涉及到的是用一個(gè)矩陣的逆去乘另一個(gè)矩陣，如果一步一步算，想要算出矩陣的逆，再算兩個(gè)矩陣相乘，沒(méi)有一步是簡(jiǎn)單的，兩步計(jì)算過(guò)程都很繁瑣，極易出錯(cuò)。仔細(xì)觀(guān)察后，我發(fā)現(xiàn)正是因?yàn)榫仃嚨哪媾c另一矩陣相乘，從而在整體上出現(xiàn)了相對(duì)簡(jiǎn)單的計(jì)算，其過(guò)程是略去矩陣逆的計(jì)算，從而簡(jiǎn)化計(jì)算。aiiXi,812X2'"

4、-amXn=b1對(duì)于n介線(xiàn)性方程組：，即Ax=b,其系數(shù)矩陣RniXi+an2X2+8nnXn=。A=3ijK非奇異且a.#0(i=i,2,3,，n),對(duì)k=0,i,2,則可建立Jacobi迭代格式：Xi(ki)ki)(ki)Xni/(k)(k)(k)(-ai2X2ai3X3-ainXnbl),aiii/(k)(k)(k)(-a2iXia23X3-a2nXnb2),a22(-aniXi-an2X2-_an,nJXn1bn)-ann我們知道Jacobi迭代矩陣為Bj=-D/(L+U)=I-D/A(2),其中aiia220ai2a13ain0a23a2nU0工：*an,nI00a2i0L=a3ia

5、320aa+,+.anian2an,n,0(3)。由式可看出，計(jì)算Bj,首先需求出D。然后再作矩陣乘法。當(dāng)然這兒由"an于D的特殊性，Da很好求，D-*=a22ia22ann10a12-a11ai3aiiain1a11a2ia23a2n0iaiia22a22Bj=I-DA=a3ia320_S3n_a33a33a33-sama*，Banian2an301annannann就會(huì)發(fā)現(xiàn)，其實(shí)Bj可以直接寫(xiě)出來(lái)，無(wú)需0ai2-a11a2i0aii計(jì)算，由可得Bj=a31a32a33aa33aanian2如果我們拋開(kāi)式，直接看,_ain_1aiia11a23a2na22a220一.a3n,其直接

6、從線(xiàn)a33a+.a_an30ann_annann性方程組中得來(lái)，顯然快于一步步的計(jì)算，而且第二中算法不僅簡(jiǎn)單還不容易出錯(cuò)，提高了求迭代矩陣的效率。當(dāng)然，第一種算法的D可以直接寫(xiě)出很好求,從而效率也沒(méi)提高多少，但對(duì)于Gauss-Seidel迭代，就不然了。Xi(k1)(-ai2x2k)-ai3x3k)k1)an1-a1nxnk)-bi),a22(-a2ixi(k1)-a23x3k)-a2nxnk)-b2),(k1)xn1z(ki)(ki)一(-anixi-an2x2ann(k1)-an,n4Xn/bn).我們知道Gauss-Seidel迭代矩陣Bg=(D+L),U(5),其中矩陣D,L,U與上述

7、中一樣。但此處（D+L），就不是太好求了，即使它是個(gè)下三角矩陣。然而求出（D+L）”后，還要進(jìn)行矩陣的乘法，因?yàn)閰^(qū)=-（D+L）“U即:Gauss-Seidel迭代格式:a2ia22Bg=(D+L)U=a3ia32a33janian2an3I0&2ai30a23I0annaina2naan,n0計(jì)算有點(diǎn)繁瑣，然而，我產(chǎn)生一種想法，其是否也可與Jacobi迭代矩陣那樣，直接寫(xiě)出來(lái)了？通過(guò)一番計(jì)算，再加上實(shí)例的體會(huì)，我找出了一種相對(duì)簡(jiǎn)潔的關(guān)系。把式寫(xiě)成xi(ki)=(-ai2x2k)-ai3x3k)-sinxnk)-h)iaii(ki)(ki)(k)(k)a22x2-a2ixia23x3-

8、a2nxnb22（k+）（k*（k書(shū)）.（Ui）,、fnnxn=-anixi-an2x2-an,n,Xn*bnn）把i）代入2）并整理得（由于我們的目的是得到矩陣，所以在此就不考慮以心,，bn了）a22x2ki)=(毋*3)x2k)aii'("a2i,-阻拓。,aii(-a2i*R-a2n)xnk)aii2')令bi2=一呢.=一加,，bin=一'1n，則i)變?yōu)閄i(k+)=bi2x2k)+h3x3k)+b1n姆aiiaiiaii(-a2i*hn-a2n)/a22xnk)。此時(shí)2'）式變?yōu)閤2k"=(-a2i*bi2)/a22x2k)(fi

9、3-a23)/a22x3k)2A)。令b22=-a21bl2,b23=-a21bi3-a23，,b2n=-a21b1n-a2n,則2A)式變?yōu)閍22a22a22(k+)u(k)x,(k).(k)。把、式代入3)式整理得x2=>2乂2飛23乂3飛2/門(mén)(k1).(k)(k)a33x3-(a31bi2-a32b22)x2(-a31bi3-a32b23)x33')a31b12-a32b22.,b33-a31b13a32b23,b34a31b14-a32b24-a34a33a33a33-a31bin-a32b2n-a3na33則3')式變?yōu)楸氐?b32x2k)+b33x3&quo

10、t;+b34x4k)+b3nXnk)如此一直在前一步的基礎(chǔ)上求后一步矩陣中的元素的值，一直進(jìn)行下去，則n-1)式變?yōu)閤nk=0,2$2+4/21+a-nxnk)則第n個(gè)式子變?yōu)閍nnxn)=(-an,1b12-an,2b22-_an,nJbnd,2)x2)，(-an,1b13-an,2b23-an,nJbnJ,3)x3)(-an,1b1,n-an,2b2,n-an,n二bi,n)xn=bn,2x2k)bn,3x3k)bn,4x4k)bn,nx(k)n從而得到Gauss-Seidel迭代矩陣b12b22bl3b23bb2一0b12b1nl10-a12-a1nBg03b22b2n9=01a11b2

11、2sa11asb2nm0bn2bnn_10bn2bnn.0b12b13-bina21b12-a21b13一a23a21bln-a2n00Ia22a22a220bn2bn3.bnn一a31bl2-a32b22一a31b1n一a32b2n-a3na33a33，%+b,naaa-a,1b1,iJL_ai,iJ，bJJ-a,1b,T"一小力山書(shū)一小斗-ai,1b,n=,'一2,ijb_|,n-"na,iai,iai,iiaabn,ibnJ+bn,n1b12b13b1nb22b23b2nb32b33b3naaa-an,1b12an,nbn2-an,1b13an,n'b

12、n,3-an,1b1nan,nbnnan,nan,nan,n0-0010-000接下來(lái)我用書(shū)上一個(gè)例子來(lái)展現(xiàn)上述方法求迭代矩陣的優(yōu)越性:8x1x2-2x3=9,例設(shè)方程組為43x1-10x2+x3=19,試分別寫(xiě)出Jacobi迭代和5x1-2x2+20x3=72,Gauss-Seidel迭代格式以及相應(yīng)的迭代矩陣。解：原方程的Jacobi迭代格式和Gauss-Seidel迭代格式分別為(k1)Xi(k1)X2(k1)X3(-x2k)2x3k)9),8(-3x1(k)-x3k)+19),和10(-5x1(k)-2x2k),72),20x2k1)=k1)x1(k1)(-x2k)2x3k)8-(-3

13、x1(k1)-x3k1)19),101(-5娟2x2k1)72),20由可直接得Jacobi迭代矩陣為Bj=0.3140.140.1而相應(yīng)的Gauss-Seidel迭代矩陣可由（*）式得:-3Bg二1811)I<8J10b321143-1410b33-0.125-0.0375-5(-0.125)2(-0.0375)0.250.175-50.2520.1752020-00-0-0.125-0.03750.02750.250.175-0.045與書(shū)上用公式算所得結(jié)果相同，但這種計(jì)算顯然很簡(jiǎn)潔。對(duì)于3階以上的迭代矩陣的計(jì)算，我的方法將會(huì)節(jié)約大量時(shí)間，而且還不容易出錯(cuò)。以上我們討論的是方陣，但從

14、（*）式可以看出，我們也可以求出不是方陣的Bg,這便給人一種想法，是否不是方陣時(shí)也可用迭代法求其解，但有一點(diǎn)是肯定的，當(dāng)方程個(gè)數(shù)少于未知數(shù)個(gè)數(shù)時(shí)，線(xiàn)性方程組有無(wú)窮多解，因此這個(gè)問(wèn)題有可能是否定的，即無(wú)法用迭代法求系數(shù)矩陣不是方陣的解，此問(wèn)題還有待研究。以下是我寫(xiě)的有關(guān)我的新方法的matlab代碼，其中也包含書(shū)上的方法求迭代矩陣的代碼，我輸入一個(gè)四階的系數(shù)矩陣，由兩種方法所得出的Gauss-Seidel迭代矩陣完全相同。Matlab代碼:%<Gauss-Seidel矩陣functionG_SB(A)m,n=size(A);ifm=ndisp('系數(shù)矩陣不是方陣)returnend%

15、用矩陣運(yùn)算求Gauss-Seidel迭代矩陣D=diag(diag(A);U=-triu(A,1);L=tril(A,-1);disp('矩陣運(yùn)算求出的迭代矩陣')B1=inv(D+L)*U,%B=zeros(m,n);forj=2:nB(1,j)=-A(1,j)/A(1,1);endfori=2:mforj=2:nk=1:i-1;ifi>=jB(i,j)=sum(-A(i,k)*B(k,j)/A(i,i);continueendB(i,j)=(sum(-A(i,k)*B(k,j)-A(i,j)/A(i,i);endenddisp(直接法求出的迭代矩陣)B,end給定系數(shù)矩陣，并得結(jié)果：>>A=51-1