求解线性方程组的直接解法5.2LU分解① Gauss消去法实现了LU分解顺序消元结束时的上三角矩阵U和所用的乘数,严格下三角矩阵。
将下三角矩阵的对角元改成1,记为L,则有A=LU,这事实是一般的,我们不难从消去的第k个元素时的矩阵k行及k列元素的历史得到这一点.因为从消元的历史有u kj=a kj-m k1u1j- m k2u2j -…- m k,k-1u k-1,j, j=k,k+1,…,nm ik=(a ik-m i1u1k- m i2u2k -…-m i,k-1u k-1,k>/u kk i=k+1,k+2,…,n于是a kj=m k1u1j+m k2u2j+…+m k,k-1u k-1,j+u kj, j=k,k+1,…,na ik=m i1u1k+m i2u2k+…+m i,k-1u k-1,k+m ik u kk i=k+1,k+2,…,n从前面两个式子我们可以直接计算L和U(见下段>.将矩阵分解为单位下三角矩阵和上三角矩阵之积称为矩阵的LU分解.顺序消元实现了LU分解,同时还求出了g, Lg=b的解.②直接LU分解上段我们得到(l ij=m ij>u kj=a kj-l k1u1j-l k2u2j -…- l k,k-1u k-1,j, j=k,k+1,…,nl ik=(a ik-l i1u1k-l i2u2k -…-l i,k-1u k-1,k>/u kk i=k+1,k+2,…,n2诸元素对应乘积,只不过算L的元素时还要除以同列对角元.这一规律很容易记住.可写成算法(L和U可存放于A>:for k=1:n-1for j=k:nu kj=a kj-l k1u1j-l k2u2j -…- l k,k-1u k-1,jendfor i=k+1:nl ik=(a ik-l i1u1k-l i2u2k -…-l i,k-1u k-1,k>/u kkendend这一算法也叫Gauss消去法的紧凑格式,可一次算得L,U的元素,不需逐步计算存储.考察上面的表格会发现还可安排其它计算次序,只要在这一次序下每个元素左边的L的元素与上方的U的元素已计算在先。
例如,逐行自左而右的次序, 逐列自上而下的次序,…易知g的计算规律同U.利用LU分解解Ax=b分三步:1．分解A=LU2．解Lg=b 求g3．解Ux=y 求x例3.用直接LU分解法解解用分解公式计算得求解③其它分解我们用顺序消元和直接分解两种方法实现了LU分解.还有更一般的三角分解,比如,下三角矩阵和单位上三角矩阵之积,又如单位下三角矩阵,对角矩阵,单位上三角矩阵之积,等等.下面给出第二种分解形式的算法LDR分解法。
A=LDR,L是单位下三角矩阵,D是对角矩阵,R是单位上三角矩阵.逐列计算(逐列作LU分解,再用U的对角元素除各行>,结果存入A。
for j=1:nfor i=2:ja ij=a ij-a i1a1j-a i2a2j -…-a i,i-1a i-1,jendfor i=j+1:na ij=(a ij-a i1a1j-a i2a2j -…-a i,j-1a j-1,j>/a jjendfor i=1:j-1a ij= a ij/a iiendend④列主元素的LU分解对照顺序消元和LU分解,列主元素法也可得列主元素的LU分解:PA=LUP是行交换结果的排列阵,L和U同前.例4.列主元素法解方程组并写出系数矩阵的LU分解.括号内是乘数,k=2时2,3行交换.因而有直接作列主元素LU分解,因为在k步要先选主元素,所以作如下改变:for k=1:n-1for i=k:na ik=a ik-l i1u1k-l i2u2k -…-l i,k-1u k-1,kend找p:p行k行i k=pfor j=k+1:nu kj=a kj-l k1u1j-l k2u2j -…- l k,k-1u k-1,jendfor i=k+1:nl ik= a ik/u kkendend可将l ik存于a ik,u kj存于a kj.二实验部分本章实验内容：实验题目：Gauss消元法,追赶法,范数。
实验内容：①编制用Gauss消元法求解线性方程组Ax=f的程序。
②编制用追赶法求解线性方程组Ax=f的程序。
③编制向量和矩阵的范数程序。
实验目的：①了解Gauss消元法原理及实现条件,熟练掌握Gauss消元法解方程组的算法,并能计算行列式的值。
②掌握追赶法，能利用追赶法求解线性方程组。
③理解向量和矩阵范数定义,性质并掌握其计算方法. 编程要求：利用Gauss消元法,追赶法解线性方程组。
分析误差。
计算算法:①Gauss消元法:1.消元过程设,对计算⒉回代过程②追赶法：1.分解Ax=f:( >2.解Lg=f,求g:(>3．解Ux=g,求x:<）③范数：常用向量范数有：<令x=( x,x2,…,x n>）11-范数：║x║1=│x1│+│x2│+…+│x n│2-范数：║x║2=(│x1│2+│x2│2+…+│x n│2>1/2∞-范数：║x║=max(│x│,│x2│,…,│x n│>1常用的三种向量范数导出的矩阵范数是：1-范数:║A║1= max{║Ax║1/║x║1=1}=2-范数:║A║2=m ax{║Ax║2/║x║2=1}=,λ1是A T A的最大特征值.∞-范数:║A║=max{║Ax║/║x║=1}=实验例题⑴：用Gauss消元法解方程组实验例题⑵:用追赶法解三对角方程组Ax=b,其中实验例题⑶：设，计算A的行列范数.程序①:Gauss消元法function x=Gauss(A,b>%A是线性方程组的系数矩阵,b为自由项.n=length(A>for k=1:n-1m(k+1:n,k>=A(k+1:n,k>/A(k,k>。
A(k+1:n,k+1:n>=A(k+1:n,k+1:n>-m(k+1:n,k>*A(k,k+1:n>。
b(k+1:n>=b(k+1:n>-m(k+1:n,k>*b(k>。
endx=zeros(n,1>。
x(n>=b(n>/A(n,n>。
for k=n-1:-1:1x(k>=(b(k>-A(k,k+1:n>*x(k+1:n>>/A(k,k>。
endx=x'。
disp(sprintf('k x(k>'>>。
for i=0:ndisp(sprintf('%d %f ',i,x(i+1>>>。
end数值结果：x=Gauss(A,b>n =3k x(k>0 1.1111111 0.7777782 2.555556程序②:追赶法function [x,y,beta]=zhuiganfa(a,b,c,f>%a,b,c是三对角阵的对角线上的元素,f是自由项.n=length(b>。
beta(1>=c(1>/b(1>。
for i=2:nbeta(i>=c(i>/(b(i>-a(i>*beta(i-1>>。
endy(1>=f(1>/b(1>。
for i=2:ny(i>=(f(i>-a(i>*y(i-1>>/(b(i>-a(i>*beta(i-1>>。
endx(n>=y(n>。
for i=n-1:-1:1x(i>=y(i>-beta(i>*x(i+1>。
enddisp(sprintf('k x(k> y(k> beta(k>'>>。
for i=0:ndisp(sprintf('%d %f ',i,x(i+1>,y(i+1>,beta(i+1>>>。
end数值结果:a=[0 -1 -1 -1 -1]'。
b=[2 2 2 2 2]'。
c=[-1 -1 -1 -1 0]'。
f=[1 0 0 0 0]'。
[x,y,beta]=zhuiganfa(a,b,c,f>k x(k> y(k> beta(k>0 0.833333 5.000000e-001 -0.5000001 0.666667 3.333333e-001 -0.6666672 0.500000 2.500000e-001 -0.7500003 0.333333 2.000000e-001 -0.8000004 0.166667 1.666667e-001 0.000000程序③：1.列范数:function fan=lie(A>%A为已知矩阵n=length(A>for j=1:nx(j>=0for i=1:nx(j>=x(j>+abs(A(i,j>>。
endendfan=max(x>disp(sprintf('n x(n>'>>。
for i=0:ndisp(sprintf(' %d %f',i,x(i+1>>>。
end数值结果:fan=lie(A>fan =0.8000n x(n>00.70000010.8000002.行范数:function fan=hang(A>%A为已知矩阵n=length(A>for i=1:nx(i>=0for j=1:nx(i>=x(i>+abs(A(i,j>>。
endendfan=max(x>disp(sprintf('n x(n>'>>。
for i=0:ndisp(sprintf(' %d %f',i,x(i+1>>>。
end数值结果:fan=hang(A>fan =1.1000n x(n>0 1.10000010.400000总结:从代数上看,直接分解法和Gauss消去法本质上一样,但如果我们采用”双精度累加”计算,那么直接三角分解法的精度要比Gauss消去法为高.求线性方程组的直接法,其算式繁杂,给人以枯燥沉闷的感觉.为了改善教案效果,本章着重介绍了三对角方程组的追赶法.三对角方程组以及其拓广形式的带状方程组有着广泛的实际应用.追赶法是解三对角线方程组(对角元占优势>的有效方法,它具有计算量少,方法简单,算法稳定等优点,具有鲜明的对称美.复习思考题1. 用消去法解线性方程组为什么最好选主元？怎样的方程组可以不用选主元？2. 用高斯—约当消去法求矩阵的逆，其理论根据是什么？3. 求矩阵A的LU分解的紧凑格式有何规律？写出LU分解法解Ax = b的算法。