《列主元高斯消去法和列主元三角分解法解线性方程》由会员分享，可在线阅读，更多相关《列主元高斯消去法和列主元三角分解法解线性方程（16页珍藏版）》请在金锄头文库上搜索。
1、计算措施实验报告1 【课题名称】用列主元高斯消去法和列主元三角分解法解线性方程【目的和意义】高斯消去法是一种古老的求解线性方程组的措施，但由它改善得到的选主元的高斯消去法则是目前计算机上常用的解低阶稠密矩阵方程组的有效措施。用高斯消去法解线性方程组的基本思想时用矩阵行的初等变换将系数矩阵A约化为具有简朴形式的矩阵(上三角矩阵、单位矩阵等），而三角形方程组则可以直接回带求解用高斯消去法解线性方程组(其中An）的计算量为:乘除法运算环节为,加减运算环节为。相比之下，老式的克莱姆法则则较为繁琐,如求解20阶线性方程组，克莱姆法则大概要次乘法，而用高斯消去法只需要60次乘除法。在高斯消去法运算的过程中
2、,如果浮现abs（A(i,i）)等于零或过小的状况，则会导致矩阵元素数量级严重增长和舍入误差的扩散，使得最后的计算成果不可靠,因此目前计算机上常用的解低阶稠密矩阵方程的迅速有效的措施时列主元高斯消去法,从而使计算成果更加精确。2、列主元三角分解法高斯消去法的消去过程，实质上是将A分解为两个三角矩阵的乘积A=LU,并求解Ly=b的过程。回带过程就是求解上三角方程组Ux=。因此在实际的运算中,矩阵L和U可以直接计算出,而不需要任何中间环节,从而在计算过程中将高斯消去法的环节进行了进一步的简略,大大提高了运算速度，这就是三角分解法采用选主元的方式与列主元高斯消去法同样，也是为了避免除数过小，从而保证
3、了计算的精确度【计算公式】1、 列主元高斯消去法设有线性方程组Ax=,其中设为非奇异矩阵。方程组的增广矩阵为 第1步（=）：一方面在A的第一列中选用绝对值最大的元素,作为第一步的主元素:然后互换（A，b)的第1行与第l行元素，再进行消元计算。 设列主元素消去法已经完毕第1步到第k-1步的按列选主元，互换两行,消元计算得到与原方程组等价的方程组A(k）xb(k) 第k步计算如下: 对于k1,2，,n-1 (1)按列选主元：即拟定使 (2）如果t,则互换，b第t行与第k行元素。 (3）消元计算消元乘数mik满足: (4）回代求解2、 列主元三角分解法对方程组的增广矩阵 通过k-步分解后，可变成如下
4、形式：第k步分解，为了避免用绝对值很小的数作除数，引进量,于是有=。如果 ，则将矩阵的第t行与第行元素互换，将(i，j)位置的新元素仍记为或,然后再做第步分解，这时【列主元高斯消去法程序流程图】【列主元高斯消去法Mtlb主程序】ucion x=gauss1(,，c) %列主元法高斯消去法解线性方程Ax=bi（legt()=lengt(b)） %判断输入的方程组与否有误 dis(输入方程有误!) turn;edis(原方程为AXb：) 显示方程组bdip(-)nleth（A);fork=1:n- 找列主元 p,q=a(bs(（k:n,k); 找出第k列中的最大值,其下标为,q q=q+-1； q
5、在(:,k)中的行号转换为在A中的行号 as(p)ktpA(k,：）; %列主元所在行不是目前行，将目前行与列主A(k,:）=A（q,:); 元所在行互换(涉及) A(q,：）=temp1; temp2=（,:); b（k,：)=b(q,:); b(q，:)=tep2;nd %消元 for ik+1:n m(i，k)A(i,）/A(k,k); %A（k,k)将A（i，k)消为0所乘系数 A(i,k:n）=A(i,:)-m(i,k）A(k,：n）; %第i行消元解决 b(）(）-m(i,k)b(k); b消元解决 endendis(消元后所得到的上三角阵是) A %显示消元后的系数矩阵b(n)=
6、b(n)(n,n); %回代求解o i=n-1:1：1 (i）=（(）-su(（i,i:n)*b(i+：n)）A（,);edclearx;dsp(AXb的解x是)x=b;【调用函数解题】【列主元三角分解法程序流程图】【列主元三角分解法alb主程序】自己编的程序:funcion xPL(A,b,eps) %定义函数列主元三角分解法函数f (ength(A)=lenh(b) 判断输入的方程组与否有误 dsp(输入方程有误!) reurn;endp(原方程为AXb：) 显示方程组Abdisp(-)n=egth(A);A=A b; 将A与b合并，得到增广矩阵fo=1: if r=1 for =：n c
7、 d=max(as(A(:，1); 选用最大列向量，并做行互换 if cps %最大值不不小于e,主元太小,程序结束 beak; else end d+1-; p(1,:）; A(1，:)=A(,:); A(d,:)=p; A(1,i)=A(1,i); end A(1,2:)A（1,:n）; A(2:,1)=(2:n,）A（1，1); %求u（,i) lse u(r，r)=A（r,r)-A（r,1：r-1）*A(1:r-1，r）； %按照方程求取(r,) f bs(u(r，r)）=s %如果（,r)不不小于e,则互换行 p=A(,:); A(r,:)=A(r+1,：)； A(r+1,：)p；
8、else nd or :n (r,)=A(r，)(,1：-1)*A(1：r-1,）; %根据公式求解，并把成果存在矩阵A中 en fr ir+1:A(,)(A（,r）-A(i,1：r-1）*A(1:-1，r)）/A(,); %根据公式求解,并把成果存在矩阵A中 end edend y(1）=A(,n1); for =2: h=0； fr ：1 h=h（,k）*y(k); end y()A(i，n+1)-; %根据公式求解(i) d x（n)=y(n)/(,n); ori=-1:：1 h0; fok=i+1：n h=h+A(i，k)*x(k）; end x()=（(i)-h)(i,); %根据公
9、式求解x() ed Adis(A=b的解x是）x； 输出方程的解可直接得到,L,U并解出方程解的的程序（查阅资料得子函数1,其作用是将矩阵A分解成L乘以U的形式。PLU2为调用PLU1解题的程序,是自己编的)()fuctio l,u,=PLU1(A) %定义子函数,其功能为列主元三角分解系数矩阵A ，n=ize（)； %判断系数矩阵与否为方阵 f m= err(矩阵不是方阵) eurn ed i det（)=0 %判断系数矩阵能否被三角分解 rror(矩阵不能被三角分解) end uA;p=eye（m）；le(); 将系数矩阵三角分解，分别求出P,L, or i=1: for j=i:m t(j)=(j,i)； for k:i-1 t(j)=t（）(j，k)*(k,i）; en ed a=i；=as(t(i)）; foj=i+1: if babs(
