三角分解解线性方程组报告

解线性方程组的Doolittle分解

一、 实验目的

1、 2、

学习和掌握线性代数方程组的Doolittle分解法。 运用Doolittle分解法进行计算。

二、方法原理

在Gauss消元法中,对于n阶方程组,应用消去法经过n-1步消元之后,得到一个与Ax=b等价的代数线性方程组A(n1)xb(n1),而且A(n1)为一个上三角矩阵.所以

我们想是否能把矩阵A分解成一个下三角阵与一个上三角阵的乘积 A=LR,其中

Lyb,L为下三角阵,R为上三角阵.就变成了两个三角形方程组Rxy 的求解问题。

三、算法描述：

Setp1：利用for循环求出rkii1akilkprpi，lik(aikp1k1lp1k1ippkr)/rkk。

Step2：

yibilikyk,得出xi(yik1ki1rx)/r

ikkiin四、流程图：

开始 输入方程右边的数和方程组的系数矩阵 r[k][j]=a[k][j]-l[kp]r[kp] p1k1l[ik]=（a[ik]-l[ip]r[ip]）/r[k][k] p1k1y[i]=b[i]-l[ik]y[k] k1i1x[i]=(y[i]-ki1r[ik]x[k])/r[i][i] n输出最终运算结果 结束 五、程序原代码如下：

#include #define N 10 void main() { int i,j,k,p,n; double l[N+1][N+1]={0},r[N+1][N+1]={0},x[N+1],y[N+1],b[N+1],a[N+1][N+1]; double s1,s2,s3,s4; printf(\"请输入矩阵的阶数n:\");

scanf(\"%d\ printf(\"请输入矩阵的右端项:\\n\"); for(i=1;i<=n;i++) {

scanf(\"%lf\ } printf(\"请输入矩阵的系数:\\n\"); for(i=1; i<=n; i++) {

for(j=1; j<=n; j++) {

scanf(\"%lf\ } } for(k=1;k<=n;k++) { for(i=k;i<=n;i++) { s1=0; for(p=1;p<=k-1;p++) s1=s1+l[k][p]*r[p][i]; r[k][i]=a[k][i]-s1; printf(\"r[%d][%d]=%lf\\n\ } for(i=k+1;i<=n;i++) { s2=0; for(p=1;p<=k-1;p++) s2=s2+l[i][p]*r[p][k]; l[i][k]=(a[i][k]-s2)/r[k][k]; printf(\"l[%d][%d]=%lf\\n\ } l[k][k]=1; }

printf(\"输出上三角矩阵r：\\n\"); for(i=1;i<=n;i++) { for(j=1;j<=n;j++) printf(\"%f \ printf(\"\\n\"); } printf(\"输出下三角单位矩阵l：\\n\"); for(i=1;i<=n;i++) {

}

}

for(j=1;j<=n;j++) printf(\"%f \ printf(\"\\n\");

for(i=1;i<=n;i++) { s3=0; for(k=1;k<=i-1;k++) s3=s3+l[i][k]*y[k]; y[i]=b[i]-s3; }

for(i=1;i<=n;i++) printf(\"y[%d]=%lf\\n\for(i=n;i>=1;i--) {

s4=0; for(k=i+1;k<=n;k++) s4=s4+r[i][k]*x[k]; x[i]=(y[i]-s4)/r[i][i]; }

for(i=1;i<=n;i++) printf(\"x[%d]=%lf\\n\

六、数值计算:

2x1x2x34x13x22x36 x2x2x5231

七、参考文献：

[1]刑志栋，矩阵数值分析， 陕西： 陕西科学技术出版社， 2005。 [2]谭浩强，C语言程序设计，北京：清华大学出版社，2005。

[3]翁惠玉，c语言程序设计思想与方法，北京：人民邮电出版社，2008