Doolittle分解
编程实现doolittle分解方法解方程组
Doolittle分解方法是一种用于解决线性方程组的数值方法，它可以将一个矩阵分解为一个下三角矩阵和一个上三角矩阵的乘积，从而可以方便地求解线性方程组。
在本文中，我们将介绍Doolittle分解方法的原理和实现过程，并用编程语言实现该方法来解方程组。
一、Doolittle分解方法原理1.1 Doolittle分解方法是一种LU分解的特例，它将一个矩阵A分解为一个下三角矩阵L和一个上三角矩阵U的乘积，即A=LU。
其中，L 的主对角线元素全为1，U的主对角线以上的元素全为0。
这样的分解可以方便地求解线性方程组Ax=b，其中b是一个已知的列向量。
1.2 Doolittle分解方法的具体实现过程是通过高斯消元法来实现的。
将矩阵A分解为一个下三角矩阵L和一个上三角矩阵U，然后通过回代法求解线性方程组Ax=b。
具体来说，我们首先将矩阵A分解为L 和U，然后用L和U的乘积代替原来的矩阵A，将原来的线性方程组Ax=b变为LUx=b，然后通过两次回代法求解线性方程组Ly=b和Ux=y，最终得到线性方程组的解x。
1.3 Doolittle分解方法的优点是可以方便地求解多个方程组，因为一旦矩阵A被分解为L和U，就可以通过多次回代法来求解不同的线性方程组，而不需要重新分解矩阵A。
1.4 Doolittle分解方法的缺点是需要对原始的矩阵A进行分解，这需要一定的计算量，特别是对于比较大的矩阵来说。
Doolittle分解方法在实际应用中往往需要结合其他数值方法来提高求解线性方程组的效率。
二、Doolittle分解方法的实现过程2.1 我们需要定义一个函数来实现Doolittle分解。
该函数的输入是一个矩阵A，输出是矩阵A的下三角矩阵L和上三角矩阵U。
2.2 接下来，我们需要通过高斯消元法来实现Doolittle分解。
具体来说，我们首先对矩阵A进行行变换和列变换，使得矩阵A的主对角线元素非零，然后逐步消去矩阵A的非主对角线元素，得到下三角矩阵L和上三角矩阵U。
gauss消元法和doolittle乘法运算次数
《深入探讨高斯消元法和Doolittle分解的乘法运算次数》在数学和计算机科学领域，高斯消元法和Doolittle分解是两种常见的线性代数运算方法。
它们被广泛用于解决线性方程组和矩阵求逆等问题。
本文将从深度和广度的角度对这两种方法进行全面评估，并进一步探讨它们的乘法运算次数的比较。
1. 高斯消元法简介高斯消元法是一种用于解决线性方程组的方法，通过矩阵变换将其转化为上三角矩阵，从而求解方程组。
其基本思想是通过一系列的行变换，将系数矩阵变换为上三角矩阵，再通过回代求解出未知数的值。
在实际应用中，高斯消元法通常需要进行大量的乘法和加法运算，其乘法运算次数随矩阵的大小而增加。
2. Doolittle分解简介Doolittle分解是将一个矩阵分解为一个下三角矩阵和一个上三角矩阵的乘积，这种分解方法可以简化矩阵的求逆和解线性方程组的计算。
与高斯消元法相比，Doolittle分解在某些情况下可以更加高效地解决线性方程组的问题，尤其是对于大型矩阵的计算。
其乘法运算次数与矩阵的大小和稀疏程度密切相关。
3. 乘法运算次数比较在实际应用中，我们常常需要比较高斯消元法和Doolittle分解的乘法运算次数，以确定哪种方法更适合特定的问题。
根据理论分析和实际测试，我们可以得出以下结论：- 对于小型矩阵，通常情况下高斯消元法的乘法运算次数略少于Doolittle分解。
- 对于大型矩阵，Doolittle分解的乘法运算次数通常比高斯消元法少很多，尤其是在矩阵稀疏的情况下。
- 对于需要多次求解的问题，Doolittle分解可以通过分解一次，多次使用的方式，进一步减少总体的乘法运算次数。
4. 个人观点和理解从个人观点来看，高斯消元法和Doolittle分解都是非常重要的线性代数运算方法，它们各有优劣。
在实际应用中，我们需要根据具体的问题特点来选择合适的方法。
对于小型矩阵或需一次性解决问题的情况，高斯消元法可能更加便捷；而对于大型矩阵或需要多次使用的情况，Doolittle分解可能更具优势。
杜立特尔三角分解法
杜立特尔三角分解法杜立特尔三角分解法- 全面评估、深度探讨与个人观点解读1. 引言杜立特尔三角分解法是线性代数中一种重要的数值方法，用于解决线性方程组。
它的广泛应用于科学计算、工程学和金融建模等领域中，具有高精度、快速稳定的特点。
本文将对杜立特尔三角分解法进行全面评估、深度探讨，并分享个人对该方法的理解与观点。
2. 杜立特尔三角分解法的基本原理杜立特尔三角分解法是一种将矩阵分解为下三角矩阵和上三角矩阵之乘积的方法。
它的基本原理是通过高斯消元法将线性方程组转化为上三角矩阵的形式，再通过前向、后向代替求解出线性方程组的解。
3. 杜立特尔三角分解法的步骤3.1 高斯消元法在杜立特尔三角分解法中，首先需要使用高斯消元法将线性方程组转化为上三角矩阵的形式。
高斯消元法通过消元和回代的操作，将方程组变换为上三角矩阵形式，即将主元以下的元素全部消为0。
3.2 前向代替在杜立特尔三角分解法的前向代替步骤中，从上到下，逐行将下三角矩阵的元素表示为主线元素与其他元素的线性组合，以求得下三角矩阵L。
3.3 后向代替在杜立特尔三角分解法的后向代替步骤中，从下到上，逐行将上三角矩阵的元素表示为主线元素与其他元素的线性组合，以求得上三角矩阵U。
3.4 解线性方程组通过前向代替和后向代替的步骤，我们可以求解出上三角矩阵和下三角矩阵。
我们可以通过反向代替的方式，将线性方程组表示为矩阵形式进行求解，即解出线性方程组的解向量。
4. 杜立特尔三角分解法的优势与应用杜立特尔三角分解法相较于其他解线性方程组的方法，具有以下优势： - 高精度：杜立特尔三角分解法在求解线性方程组时能够提供高精度的解。
- 快速稳定：该方法的时间复杂度较低，求解效率高，且对数值误差具有较高的容忍度。
4.1 科学计算应用杜立特尔三角分解法在科学计算领域中应用广泛，特别是在数值模拟和计算物理学等领域。
它可以在矩阵形式下求解微分方程、计算特征值和特征向量等。
4.2 工程学应用在工程学中，杜立特尔三角分解法常用于求解稀疏矩阵的线性方程组。
矩阵的Doolittle递归分解算法及符号程序设计
A22 = L22 U 22 .
即对 A22 进行 Doolit tle 分解
3 算法求精
由上述推导知 , 若矩阵可进行 Doolit tle 分解 , 使用一次分解就可求出 L 的前 d 行与 U 的前 d 行 , 接下来对 ( m - d) 行 ( n - d) 列矩阵 A′ 22 进行 Doolit 2
1) 由上述公式得到 A 的 Doolit tle 分解算法 for ( k = 0 ; k < m 且 k < n ; k + + ) {for ( j = k ;j < n ;j + + ) r kj ← a kj ;
l kk ← 1; fo r (i = k + 1 ;i < m ;i + + ) lik ← aik / a kk ; fo r (i = k + 1 ;i < m ;i + + ) for ( j = k + 1 ;j < n ;j + + ) aij ← aij - lik u kj ; 矩阵 L 和 U 的其他元素为 0 .
}
2 D ool i tt le 分解递归算法的推导
设 A 是 m 行 n 列的矩阵 , L 是 m 行 p 列矩阵 ,
U 是 p 行 n 列矩阵 , 可以进行 Doolittle 分解 , 则 Am ×n = L m ×p U p ×n , p = min ( m , n) . A11 d ×d A21 ( m- d) ×d L11 d ×d L21 ( md) × d
Doolittle分解
矩阵数值分析实验报告专业信息与计算科学班级学号姓名指导教师Doolittle 分解法一、实验目的在Gauss 消元法中,对于n 阶方程组,应用消去发经过n-1步消元之后,得到一个与Ax=b 等价的代数线性方程组)1()1(--=n n b x A ,而且)1(-n A 为一个上三角矩阵.所以我们想是否能把矩阵A 分解成一个下三角阵与一个上三角阵的乘积 A=LR,其中L 为下三角阵,R 为上三角阵.就变成了两个三角形方程组⎩⎨⎧==yRx b Ly , 的求解问题。
二、算法思想Setp1：利用for 循环求出r[k][j]=a[k][j]-∑-=1k 1p ]kp [r ]kp [l ，l[ik]=（a[ik]-∑-=1k 1p ]ip [r ]ip [l ）/r[k][k]。
Step2：y[i]=b[i]-∑-=1i 1k ]k [y ]ik [l ,得出x[i]=(y[i]-∑+=n1i k ]k [x ]ik [r )/r[i][i].三、程序代码#include <stdio.h>#include <stdlib.h>#define N 10 //矩阵大小范围float getmx(float a[N][N], float x[N], int i, int n) {float mx = 0;int r;for(r=i+1; r<n; r++){mx += a[i][r] * x[r];}return mx;}float getmy(float a[N][N], float y[N], int i, int n) {float my = 0;int r;for(r=0; r<n; r++){if(i != r) my += a[i][r] * y[r];}return my;}float getx(float a[N][N], float b[N], float x[N], int i, int n) {float result;if(i==n-1) //计算最后一个x的值result = (float)(b[i]/a[n-1][n-1]);else //计算其他x值（对于公式中的求和部分，需要调用getmx()函数） result = (float)((b[i]-getmx(a,x,i,n))/a[i][i]);return result;}float gety(float a[N][N], float b[N], float y[N], int i, int n) {float result;if(i==0) //计算第一个y的值result = float(b[i]/a[i][i]);else //计算其他y值（对于公式中的求和部分，需要调用getmy()函数） result = float((b[i]-getmy(a,y,i,n))/a[i][i]);return result;}void main(){float l[N][N]={0}; //定义L矩阵float u[N][N]={0}; //定义U矩阵float y[N]={0}; //定义数组Yfloat x[N]={0}; //定义数组Xfloat a[N][N]={{0},{0},{0}}; //定义系数矩阵float b[N]={0}; //定义右端项float sum=0;int i,j,k;int n=0;int flag=1;//用户手工输入矩阵while(flag){printf("请输入系数矩阵的大小：");scanf("%d", &n);if(n>N){printf("矩阵过大！\n");continue;}flag=0;}printf("请输入系数矩阵值：\n");for(i=0; i<n; i++){for(j=0; j<n; j++){printf("a[%d][%d]: ", i, j);scanf("%f", &a[i][j]);}}printf("请输入右端项数组：\n");for(i=0; i<n; i++){printf("b[%d]: ", i);scanf("%f", &b[i]);}/*显示原始矩阵*/printf("原始矩阵：\n");for(i=0; i<n; i++){for(j=0; j<n; j++)printf("%0.3f ",a[i][j]);printf("\n");}printf("\n\n");/*初始化矩阵l*/for(i=0; i<n; i++){for(j=0; j<n; j++){if(i==j) l[i][j] = 1;}}/*开始LU分解*//*第一步：对矩阵U的首行进行计算*/for(i=0; i<n; i++){u[0][i] = (float)(a[0][i]/l[0][0]); }/*第二步：逐步进行LU分解*/for(i=0; i<n; i++){/*对“L列”进行计算*/for(j=i+1; j<n; j++){for(k=0,sum=0; k<n; k++){if(k != i) sum += l[j][k]*u[k][i];}l[j][i] = (float)((a[j][i]-sum)/u[i][i]); }/*对“U行”进行计算*/for(j=i+1; j<n; j++){for(k=0,sum=0; k<n; k++){if(k != i+1) sum += l[i+1][k]*u[k][j]; }u[i+1][j] = (float)((a[i+1][j]-sum));}}/*输出矩阵l*/printf("矩阵L：\n");for(i=0; i<n; i++){for(j=0; j<n; j++){printf("%0.3f ", l[i][j]);}printf("\n");}/*输出矩阵u*/printf("\n矩阵U：\n");for(i=0; i<n; i++){for(j=0; j<n; j++){printf("%0.3f ", u[i][j]);}printf("\n");}/*回代方式计算数组Y*/for(i=0; i<n; i++){y[i] = gety(l,b,y,i,n);}/*显示数组Y*/printf("\n\n数组Y：\n");for(i=0; i<n; i++){printf("y%d = %0.3f\n", i+1,y[i]); }/*回代方式计算数组X*/for(i=n-1; i>=0; i--){x[i] = getx(u,y,x,i,n);}/*显示数组X*/printf("\n\n数组X：\n");for(i=0; i<n; i++){printf("x%d = %0.3f\n", i+1,x[i]); }}四、运行结果五、参考文献[1]刑志栋，矩阵数值分析，陕西：陕西科学技术出版社， 2005。
线性方程组直接解法实验
实验一 线性方程组直接解法实验一、实验目的1.运用matlab 软件完成线性方程组的直接实验；2.通过实验,了解Doolittle 分解方法和列主元消去法解方程组的过程，并比较两种方法的优点。
二、实验题目分别用Doolittle 分解方法和列主元消去法解方程组123410-7018-3 2.09999962 5.9000015-15-1521021⎛⎫⎛⎫⎛⎫⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪⎝⎭⎝⎭⎝⎭x x x x . 输出A ，b ；Doolittle 分解方法的L 和U ；解向量x,det A ;列主元方法的行交换次序，解向量x,det A ；比较两种方法所得的结果。
三、实验原理1) Doolittle 分解方法的原理算法原理：应用高斯消去法解n 阶线性方程Ax b =经过1n -步消去后得出一个等价的上三角形方程组()()n n A x b =，对上三角形方程组用逐步回代就可以求出解来。
这个过程也可通过矩阵分解来实现。
将非奇异阵分解成一个下三角阵L 和上三角阵U 的乘积称为对矩阵A 的三角分解，又称LU 分解。
根据LU 分解,将Ax b =分解为Ly bUx y =⎧⎨=⎩形式,简化了求解问题。
程序框图：变量说明：ij a 为系数矩阵元素，i b 为常数矩阵系数，,ij ij l u 分别为下、上三角矩阵元素。
2)列主元消去法解方程组的原理算法原理：列选主元是当变换到第k步时，从k列的kk a及以下的各元素中选取绝对值a的位置上，然后再进行消元过程。
交换系数矩阵中最大的元素，通过行交换将其交换到kk的两行（包括常数项），相当于两个方程的位置交换了。
程序框图：Array变量说明：k表示消元到a为消元第k步时第k步，kk主对角线元素3)四、实验过程及结果1)Doolittle分解方法的输出结果----------计算实习题----------Page64 第1题用Doolittle分解方法解方程组A =10.0000 -7.0000 0 1.0000-3.0000 2.1000 6.0000 2.00005.0000 -1.0000 5.0000 -1.00002.0000 1.0000 0 2.0000b =8.00005.90005.00001.0000L =1.0e+006 *0.0000 0 0 0-0.0000 0.0000 0 00.0000 -2.5000 0.0000 00.0000 -2.4000 0.0000 0.0000 U =1.0e+007 *0.0000 -0.0000 0 0.00000 -0.0000 0.0000 0.00000 0 1.5000 0.57500 0 0 0.0000 X =-0.0000-1.00001.00001.0000det(A)值为-762.00009000----------输出完毕----------2)列主元消去法输出结果----------计算实习题----------Page64 第1题列主元消去法解方程组A =10.0000 -7.0000 0 1.0000-3.0000 2.1000 6.0000 2.00005.0000 -1.0000 5.0000 -1.00002.0000 1.0000 0 2.0000b =8.00005.90005.00001.0000X =0.0000-1.00001.00001.0000detA值为-762.00009000----------输出完毕----------五、实验分析1.运用LU分解法可以成批地解方程组,且速度快.若c先求LU=A3,再解(LU)x=b,则要重新计算,计算量增加;如果按照上述方法计算,能够减少运算次数,加快运算速度.3. ⑴无论当n=10、n=100、n=1000时,x1与x2的值都相等,且随着n的增大,变化的只是解的中间部分数字,头、前后几位数都没有变化.⑵高斯消去法应用于三对角方程组得到的就是所谓的“追赶法”.追赶法不需要对零元素计算，只有6n-5次乘除法计算量，求解速度快.且当系数矩阵对角占优时数值稳定，是解三对角方程组的优秀解法.⑶用LU分解法解此方程组速度慢.顺序高斯消去法实际上就是将方程组的系数矩阵分解成单位下三角矩阵与上三角矩阵的乘积.顺序高斯消去法的消元过程相当于LU分解过程和Ly=b的求解,回代过程则相当于解线性方程组Ux=y,故其求解速度慢.六、附原程序1)Doolittle分解方法原程序fprintf('----------计算实习题----------\n')fprintf('Page64 第1题用Doolittle分解方法解方程组\n')A=[10 -7 0 1 ; -3 2.099999 6 2 ;5 -1 5 -1 ; 2 1 0 2];b=[8;5.900001;5;1];n=length(A);U=zeros(n,n);L=eye(n,n);U(1,:)=A(1,:);L(2:n,1)=A(2:n,1)/U(1,1);for i=2:n;U(i,i:n)=A(i,i:n)-L(i,1:i-1)*U(1:i-1,i:n);L(i+1:n,i)=(A(i+1:n,i)-L(i+1:n,1:i-1)*U(1:i-1,i))/U(i,i); endY=zeros(n);Y(1)=b(1);for i=2:nY(i)=b(i)-L(i,1:i-1)*Y(1:i-1,1);endX=zeros(n,1);if det(U)==0;X=0;elseX(n)=Y(n)/U(n,n);for i=n-1:-1:1X(i)=(Y(i)-U(i,i+1:n)*X(i+1:n,1))/U(i,i);endendAbLUXfprintf('det(A)值为%9.8f\n',det(A))fprintf('----------输出完毕 ----------\n')2)列主元消去法原程序fprintf('----------计算实习题----------\n')fprintf('Page64 第1题列主元消去法解方程组\n')A=[10 -7 0 1 ; -3 2.099999 6 2 ;5 -1 5 -1 ; 2 1 0 2];b=[8;5.900001;5;1];C=[A b];n=length(A);D=zeros(n,n+1);l=zeros(n,1);for i=1:nD=C;k=min(find(C(i:n,i)==max(C(i:n,i))));C(i,i:n+1)=D(k+i-1,i:n+1);C(k+i-1,i:n+1)=D(i,i:n+1);l(i+1:n,1)=C(i+1:n,i)/C(i,i);C(i+1:n,i:n+1)= C(i+1:n,i:n+1)- l(i+1:n,1)*C(i,i:n+1); endX=zeros(n,1);X(n)=C(n,n+1)/C(n,n);for i=n-1:-1:1X(i)=(C(i,n+1)-C(i,i+1:n)*X(i+1:n,1))/C(i,i); endAbXfprintf('detA值为%9.8f\n',det(A))fprintf('----------输出完毕----------\n')。
Doolittle分解法的构造过程
Doolittle 分解法的构造过程把矩阵A 写成下列两个矩阵相乘的形式：A=LU ，其中L 为下三角矩阵，U 为上三角矩阵。
这样我们可以把线性方程组 Ax=b 写成Ax=（LU)x=L( U x ) = b令 U x=y,则原线性方程组 Ax=b 化为两个简单的三角方程组：Ux=y 和Ly=b 。
于是可首先 求解Ly=b 得到向量y ，然后求解 Ux=y,从而求解线性方程组 Ax=b 的目的。
设矩阵A 的Doolittle 分解为： 为求出矩阵L 和U ，根据矩阵乘法规则，有用公式写法有：()ijjj j ii i i i j i a u u u l l l l U L j =⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎭⎫⎝⎛=⋅-000,,0,1,,,,213211 LU u u u u u u l l l a a a a a a a a a A nn n n n n nn n n n n =⎪⎪⎪⎪⎪⎭⎫ ⎝⎛⎪⎪⎪⎪⎪⎭⎫ ⎝⎛=⎪⎪⎪⎪⎪⎭⎫ ⎝⎛= 222112112121212222111211111{}∑∑==⋅=⋅=nk j i k kjik kj ik ij u l u l a 1,min 1ni u l u l u l a nj u u l u l u l a i nk k k ik k ik i jj nk k kj k kj k j ,,2,,2,11111111111111111111 =⋅=⋅=⋅===⋅=⋅=⋅=∑∑∑∑====iii k ki jkn k i k ki jkki jkji iji k kj ik nk ik kj ikkj kij u l u lu lu la u u l u lu la j i ji ⋅+⋅=⋅=⋅=+⋅=⋅=⋅=≤∑∑∑∑∑∑-===-===111111111,有时当于是可得Doolittle 分解公式：u 1j =a 1j j=1,2,…,n l i1=a i1 / u 11 i=2,3,…,n2）Doolittle 分解法算法1.输入变量个数n 、系数矩阵A 、常数项b2 如果a 11=0,则输出“LU 分解失败”提示并终止，否则a j1 ⇐ a j1/a 11 （j=2,…,n ）ij u u l a l ij u l a u iii k ki jk ji i k kjik ij ij ji >⋅-=≥⋅-=∑∑-=-= /)(1111注：因为 L 和U 中的三角零元素都不必存储，L 的对角元素也因为都是1没有必要再记录在程序中。
用matlab编程做doolittle 分解
用matlab编程做doolittle 分解Doolittle分解是一种用于识别系统频率响应中模态重叠的方法。
在模态分析中，模态重叠是一个常见问题，因为它可能导致错误的模态参数估计。
为了解决这个问题，研究人员开发了Doolittle分解方法。
Doolittle分解基于频率响应函数的局部平均，以识别独立的模态分量。
在本文中，我们将使用MATLAB编程语言演示如何使用Doolittle分解来识别系统中的独立模态分量。
首先，我们需要安装和加载所需的MA TLAB工具箱。
要执行Doolittle分解，我们需要使用Control System Toolbox中的信号处理工具。
要安装此工具箱，请按照MATLAB的说明进行操作。
安装完成后，加载Control System Toolbox和Signal Processing Toolbox。
接下来，我们需要创建一个模拟系统。
我们可以使用Transfer Function Designer工具箱创建一个简单的一阶系统，其传递函数为：G(s) = 1 / (s + 2)将此传递函数添加到MATLAB工作空间中，以便稍后使用。
现在，我们可以编写MATLAB脚本来执行Doolittle分解。
以下是一个简单的脚本示例：```matlab```% 创建模拟系统sys = ss(1, [1, 2], [], []);```% 添加已知的传递函数到系统中sys = addTerminalBlock(sys, [[1 0 0 1]], 1);```% 添加未知的传递函数到系统中sys = addTerminalBlock(sys, [[0 1 0 1]], 2);```% 计算频率响应函数h = tf(sys);% 将频率响应函数转换为数值表示形式num, den = numden(h);```% 设置频率向量freqs = logspace(-1, 2, 500);% 初始化Doolittle分解变量modes = zeros(size(den, 1), 1);weights = zeros(size(den, 1), 1);local_averages = zeros(size(den, 1), 1);num_modes = zeros(size(den, 1), 1);```% 执行Doolittle分解for i = 1:size(den, 1)local_averages(i) = mean(abs(zeros(1, length(freqs)) .* num(i, :) / den(i, :)));weights(i) = sum(abs(num(i, :)).^2) / sum(abs(den(i, :)).^2);modes(i) = local_averages(i) / weights(i);num_modes(i) = i;end```% 可视化结果figure;plot(freqs, modes);xlabel('Frequency (rad/s)');ylabel('Normalized mode shapes');title('Doolittle Decomposition of Simulated System');```在这个脚本中，我们首先创建了一个模拟系统，并添加了已知的和未知的传递函数。
矩阵数值分析doolittle分解
《矩阵数值分析》一、目的意义（1）.学会将方程组的系数矩阵分解成单位上三角和下三角；（2）.学会用c 语言解决矩阵分解。
二、内容要求Doolitte 分解问题：将给定的方程的系数矩阵和右端项输入即可得出方程的解。
三、问题解决的方法与算法1.方法第一步：LR 分解Doolittle 分解报告1若n 阶线性方程组系数矩阵非奇异，则有A=LR 其中L 为单位下三角，R 为上三角，按照矩阵相乘的运算法则，容易求出L 及R 的元素，A=LR 可以写成 ⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡ann an an n a a a n a a a ....21................2....22211....1211 = ⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡1....2ln 1ln ......1211l ⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡rnn n r r n r r r ..........2.....221.....1211 比较等号两边第I 行和第j 列的元素，可以得到L 和R 矩阵的元素满足下列数学公式ni j r r l a l ni i j r l a r i k ii ki jk ji ji i k kj ik ij ij ,...1,/)(,....1,,1111+=-=+=-=∑∑-=-=第二步： 求解由于将A 分解为LR ，故,b LRx Ax == 分两步求解 ⎩⎨⎧==yRx b Ly ，先求出y ，再求出x ，其计算公式：⎪⎪⎩⎪⎪⎨⎧-=-==-=∑∑+=-=n i k ii k ik i i i k k ik i i n n i r x r y x n i y l b y 1111,...1,,/)(,...2,1, 对于矩阵是用数组表示，将矩阵的分解过程转换成程序语言，最后用代 码实现出来，然后进行矩阵的各种基本运算。
2算法：设A 为方程组的系数矩阵，则A=LU ，其中L 为单位上三角矩阵，U 为 单位 下三角矩阵，其中L 中的在对角线上方的元素等于零，U 中的元素在对角线下方的等于零，L 中的不为零除对角线的元素为ik l =ik a -∑-=11k j jk ij r l ,(i=k ，k+1,……，n);U中的不为零的元素为kk k m mj km kj kj l r l a /)(r 11∑-=-=（j=k+1,k+2,……，n ）;计算公式：⎪⎪⎩⎪⎪⎨⎧-=-=∑∑-=+=ii i j j ij i i n j k kjk j l y l b x r y /)(y x 111j （i=1,2,......，n:j=n,n+1, (1)四、流程图五、计算程序#include<stdio.h>#include<math.h>void main(){double a[100][100],r[100][100],l[100][100],x[100],y[100],b[100];printf("请输入矩阵的阶数：");int n,j,i,k,p;double sum,sum2,sum3,sum4;scanf("%d",&n);printf("请输入系数矩阵：\n");for(i=1;i<=n;i++){for(j=1;j<=n;j++){scanf("%lf",&a[i][j]);}}printf("请输入矩阵的右端项:\n");for(i=1;i<=n;i++){scanf("%lf",&b[i]);}for(k=1;k<=n;k++){for(i=k;i<=n;i++){p=1;sum=0;while(p<=k-1){sum=sum+l[k][p]*r[p][i];p++;}r[k][i]=a[k][i]-sum;}for(i=k+1;i<=n;i++){p=1;sum2=0;while(p<=k-1){sum2=sum2+l[i][p]*r[p][k];p++;}l[i][k]=(a[i][k]-sum2)/r[k][k];}}for(i=1;i<=n;i++){k=1;sum3=0;while(k<=i-1){if(i==k)l[i][k]=1;else{sum3=sum3+l[i][k]*y[k];k++;}}y[i]=b[i]-sum3;}for(i=n;i>=1;i--){k=i+1;sum4=0;while(k<=n){sum4=sum4+r[i][k]*x[k];k++;}x[i]=(y[i]-sum4)/r[i][i];printf("%lf ",x[i]);}}六、计算结果与分析计算题目⎪⎩⎪⎨⎧=++=++=++52242342321321321x x x x x x x x x计算结果见下图：结果分析：可以得出，计算精度还是很准确的，计算起来比较方便，但是鉴于能力有限，还有很多不足，还应多学习改进。
矩阵分解及应用
引言数学是人类历史中发展最早，也是发展最为庞大的基础学科。
许多人说数学是万理之源，因为许多学科的研究都是以数学做为基础，有了数学的夯实基础，人类才铸就起了众多学科的高楼大厦，所以数学的研究和发展一直在不断的发展壮大。
在数学中有一支耀眼的分支，那就是矩阵。
在古今矩阵的研究发展长河中产生了许多闪耀星河的大家。
英国数学大家詹姆斯·约瑟夫·西尔维斯特，一个数学狂人，正是他的孜孜不倦的研究使得矩阵理论正式被确立并开启了矩阵发展的快速发展通道。
凯莱和西尔维斯特是非常要好的朋友，他也是一位非常伟大的数学大师，正是他们伟大的友谊，加上两人的齐心协力最后他们共同发展了行列式和矩阵的理论。
后来高斯在矩阵方面的研究取得重要的成就，尤其是高斯消去法的确立，加速了矩阵理论的完善和发展。
而在我国，矩阵的概念古已有之。
从最早的数学大家刘徽开始我们古代数学大家都已或多或少的研究了矩阵。
尤其在数学大家刘徽写的《九章算术》中，它最早提出了矩阵的类似定义。
而且是将矩阵的类似定义用在了解决遍乘直除问题里了。
这已经开始孕育出了最早的矩阵形式。
随着时间转移，矩阵的理论不断的完善，在对于那些大型矩阵的计算中如果用基本方法显得过于繁重，于是发展出了矩阵的分解，随着对矩阵分解的不断研究完善，矩阵分解方法和理论也日趋成熟矩阵经常被当做是数学工具，因为在数学问题中要经常用上矩阵的知识。
矩阵是一个表格，要掌握其运算法则，作为表格的运算与数的运算既有联系又有差别，在所有矩阵的运算方法中，矩阵的分解是他们中一种最重要并且也是应用最广泛。
矩阵分解主要是对高斯消去法的延续和拓展。
在一些大型的矩阵计算中，其计算量大，化简繁杂，使得计算非常复杂。
如果运用矩阵的分解，将那些大型矩阵分解成简单的矩阵的乘积形式，则可大大降低计算的难度以及计算量。
这就是矩阵分解的主要目的。
而且对于矩阵的秩的问题，特征值的问题，行列式的问题等等，通过矩阵的分解后都可以清楚明晰的反应出来。
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
Doolittle 分解
班级：计算062 姓名：王保翔 学号：3060811028
目的意义：在Gauss 消元法中,对于n 阶方程组,应用消去发经过n-1步消元之后,
得到一个与Ax=b 等价的代数线性方程组)1()1(--=n n b x A ,而且)1(-n A 为一个上三角矩阵.所以我们想是否能把矩阵A 分解成一个下三角阵与一个上三角阵的乘积 A=LR,其中L 为下三角阵,R 为上三角阵.就变成
了两个三角形方程组⎩⎨⎧==y
Rx b Ly , 的求解问题. 算法：
Setp1：利用for 循环求出r[k][j]=a[k][j]-∑-=1k 1p ]kp [r ]kp [l ，l[ik]=（a[ik]-∑-=1
k 1p ]ip [r ]ip [l ）
/r[k][k]。
Step2：y[i]=b[i]-∑-=1i 1k ]k [y ]ik [l ,得出x[i]=(y[i]-
∑+=n
1i k ]k [x ]ik [r )/r[i][i]. 源程序：
#include<stdio.h>
#define N 10
void main()
{
int i,j,k,n;
double a[N+1][N+1],l[N+1][N+1],r[N+1][N+1];
double x[N+1],y[N+1],b[N+1];
double p,q,t,s;
printf("Input n:");
scanf("%d",&n);
for(i=1;i<=n;i++)
{
printf("Input b[%d]:",i);
scanf("%lf",&b[i]);
}
for(i=1;i<=n;i++)
for(j=1;j<=n;j++)
{
printf("Input a[%d][%d]:",i,j);
scanf("%lf",&a[i][j]);
}
for(i=1;i<=n;i++)
{
for(j=i;j<=n;j++)
{
for(k=1,p=0;k<=i-1;k++)
p=p+l[i][k]*r[k][j];
r[i][j]=a[i][j]-p;
}
for(j=i+1;j<=n;j++)
{
for(k=1,q=0;k<=i-1;k++)
q=q+l[j][k]*r[k][i];
l[j][i]=(a[i][j]-q)/r[i][i];
}
l[i][i]=1;
}
for(i=1;i<=n;i++)
{
for(j=i;j<=n;j++)
printf("r[%d][%d]:%lf\n",i,j,r[i][j]);
for(j=i+1;j<=n;j++)
printf("l[%d][%d]:%lf\n",j,i,l[j][i]);
}
for(i=1;i<=n;i++)
{
for(k=1,t=0;k<=i-1;k++)
t=t+l[i][k]*y[k];
y[i]=b[i]-t;
}
for(i=1;i<=n;i++)
printf("Output y[%d],b[%d]:%lf,%lf\n",i,i,y[i],b[i]); for(i=n;i>=1;i--)
{
for(k=i+1,s=0;k<=n;k++)
s=s+r[i][k]*x[k];
x[i]=(y[i]-s)/r[i][i];
}
for(i=1;i<=n;i++)
printf("Output x[%d]:%lf\n",i,x[i]);
}
算例及结果分析：
分析：输入矩阵A对应输出单位下三角矩阵L和上三角矩阵R，将L乘以R得到的就是A，说明结果可靠。
输入b[i]解方程组Ax=b。
参考文献
[1]刑志栋，矩阵数值分析，陕西：陕西科学技术出版社，2005。
[2]谭浩强，C语言程序设计，北京：清华大学出版社，2005。
[3]翁惠玉，c语言程序设计思想与方法，北京：人民邮电出版社，2008。