资源描述:
《线性方程矩阵三角分解法课件.ppt》由会员上传分享,免费在线阅读,更多相关内容在教育资源-天天文库。
1、计算物理矩阵的直接三角分解法列主元三角分解法`1````````````````````````````````````````1计算物理矩阵三角分解法Gauss消元法:增广矩阵(A
2、b)经过若干次行初等变换,使之化为(U
3、g)形式,其中矩阵U为上三角矩阵由线性代数理论可知,对一个矩阵进行一次初等行变换,相当于该矩阵左乘一个相应的初等矩阵,因此Gauss消元过程可表示为其中为单位下三角矩阵,即若记(仍为下三角矩阵),则Gauss过程隐含了系数矩阵的一个矩阵三角分解计算物理三角分解定义设A是n阶方阵(n>=2),则称A
4、=LU为矩阵A的一个三角分解,其中L是下三角矩阵,U是上三角矩阵Doolittle分解L是单位下三角矩阵,U是上三角矩阵Crout分解L是下三角矩阵,U是单位上三角矩阵定理若n阶矩阵A的k阶顺序主子式不为零,则A有唯一的Doolittle分解与Crout分解计算物理证明可逆得证计算物理直接三角分解法设矩阵A存在Doolittle分解A=LU,则方程组计算物理计算物理为得到方程组的解,需先求出矩阵A的Doolittle分解由矩阵乘法得U的第一行L的第一列计算物理U的第k行、L的第k列计算物理为了节约存储空间,方便编程计
5、算,将每一步中的结算结果u、l值存放在对应的a所占的空间位置计算物理直接三角分解法解题步骤先求U矩阵的行元素,再求L矩阵的列元素一行一列交叉进行计算,将计算值存放到矩阵A的对应位置得到完整的U、L矩阵系数矩阵A计算列向量y得到方程组解x计算物理列主元三角分解法目标:与Gauss列主元消元法的原理相同,为了消除因除数过小而产生的误差扩散第K步计算物理为了消除误差扩散,应当满足但分解公式中不一定满足该条件,此时选行号,使若,则交换矩阵A中第k行和第行的元素此时一定满足(是中的主元)计算物理直接三角分解法解题步骤求U矩阵的
6、第一行元素与L矩阵的第一列元素一行一列交叉进行计算,将计算值存放到矩阵A的对应位置得到完整的U、L矩阵增广矩阵A
7、b计算列向量y得到方程组解x列主元三角分解法解题步骤第一步:求出第一列的主元,交换行第k步:求出主元,交换行求U矩阵的第k行元素与L矩阵的第k列元素计算物理例:用列主元三角分解法解下列方程组Matlab程序计算物理a=[2,1,2;4,5,4;6,-3,5];b=[6;18;5];n=3;A=[a,b];u=[];l=[];s=[];y=[];x=[];p=A(1,1);q=1;fori=1:nifp>=
8、A(i,1)p=p;q=q;elsep=A(i,1);q=i;end;end;forj=1:n+1t=A(1,j);A(1,j)=A(q,j);A(q,j)=t;end;forj=1:nu(1,j)=A(1,j);A(1,j)=u(1,j);end;fori=2:nl(i,1)=A(i,1)/u(1,1);A(i,1)=l(i,1);end;fork=2:nfori=k:nd=0;form=1:k-1d=d+A(i,m)*A(m,k);end;s(i,k)=A(i,k)-d;end;p=s(k,k);q=k;fori
9、=k:nifp>=s(i,k)p=p;q=q;elsep=s(i,k);q=i;end;end;forj=1:n+1t=A(k,j);A(k,j)=A(q,j);A(q,j)=t;end;计算物理forj=k:nd=0;form=1:k-1d=d+A(k,m)*A(m,j);end;u(k,j)=A(k,j)-d;A(k,j)=u(k,j);end;ifk<=n-1fori=k+1:nd=0;form=1:k-1d=d+A(i,m)*A(m,k);end;l(i,k)=1/u(k,k)*(A(i,k)-d);A(i,
10、k)=l(i,k);end;end;end;y(1,1)=A(1,n+1);fork=2:nd=0;form=1:k-1d=d+A(k,m)*y(m,1);end;y(k,1)=A(k,n+1)-d;end;x(n,1)=y(n,1)/A(n,n);fork=n-1:-1:1d=0;form=k+1:nd=d+A(k,m)*x(m,1);end;x(k,1)=1/A(k,k)*(y(k,1)-d);end;计算物理Thankyou!