数值代数上机报告

一、 目的意义：

把矩阵A 分解成一个下三角阵与一个上三角阵的乘积，即 A=LR ，其中L 为下

三角阵，R 为上三角阵，这样原线性方程组就可以化为⎩

⎨⎧==y Rx b Ly 的求解问题，方

便求解。

二、 算法：

1) 输入系数矩阵 A ； 2) 利用公式

ij

i k ij ij ij r l a r ∑-=-=1

1

和

ki i k jk ji ji r

r l a l /)(1

1

∑-=-=交错进

行，计算得出矩阵L 和R ；

3) 回带到⎩

⎨⎧==y Rx b

Ly 中得出原线性方程组的解。

三、 源程序：

#include

#include #include #include #define N 100 main() {

int i,j,k,s,n;

printf("请输入系数矩阵A 的阶数：n= "); scanf("%d",&n); float

a[N][N]={0},L[N][N]={0},R[N][N]={0},sigma1,sigma2,b[N],y[N],x[N];

/*为L 主对角线元素赋1*/ for(i=0;i

printf("请输入系数矩阵A ：\n"); /*输入系数矩阵A*/

报告

1

Doolittle 分解

for(i=0;i

{

for(j=0;j

scanf("%f",&a[i][j]);

}

for(k=0;k

{

for(j=k;j

{

sigma1=0;

for(s=0;s<=k-1;s++)

sigma1+=L[k][s]*R[s][j];

R[k][j]=a[k][j]-sigma1;

}

for(i=k;i

{

sigma2=0;

for(s=0;s<=k-1;s++)

sigma2+=L[i][s]*R[s][k];

L[i][k]=(a[i][k]-sigma2)/R[k][k];

}

}

printf("\n A矩阵为：\n");/*输出矩阵L、R*/ for(i=0;i

{

for(j=0;j

printf("%5.1f ",a[i][j]);

printf("\n");

}

printf("\n L矩阵为：\n");

for(i=0;i

{

for(j=0;j

printf("%5.1f ",L[i][j]);

printf("\n");

}

printf("\n R矩阵为：\n");

for(i=0;i

{

for(j=0;j

printf("%5.1f ",R[i][j]);

printf("\n");

}

printf("请输入b矩阵\n",i+1);

for(i=0;i

scanf("%f",&b[i]);

for(i=0;i

{

sigma1=0;

for(k=0;k<=i-1;k++)

sigma1+=L[i][k]*y[k];

y[i]=b[i]-sigma1;

}

for(i=n-1;i>=0;i--)

{

sigma2=0;

for(k=i+1;k

sigma2+=R[i][k]*x[k];

x[i]=(y[i]-sigma2)/R[i][i];

}

printf("解得x为：\n");

for(i=0;i

printf("%5.1f ",x[i]);

printf("\n");

}

四、计算结果与分析：

分析：运行结果与预想的结果相近，误差对结果的影响不是很大，比较理想

五、 参考文献：

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

一、 目的意义：

对称正定矩阵是实践中经常遇到的一种特殊矩阵类型矩阵，由于矩阵本身的流量好兴致，使得cholesky 分解在存储和运算量上较一般消去法节省一半左右，且解的精度高，cholesky 分解方法是目前计算机求解该类问题最有效的方法之一。

二、 算法：

1) 输入系数矩阵 A ；

2) 利用公式 ik

i k ik ii ii jji

ij ij jk

j k ik ij ij l s a d d s l l s a s ∑∑-=-=-==-=1

1

1

1

/交错进行，计算得出矩阵L 和D ；

3) 回带到⎩⎨⎧==y

D x L b

Ly T

T 中得出原线性方程组的解 三、 源程序：

#include

#include #include #define EPS 1.0e-8 #define N 20

double a[N][N], b[N], x[N];

报告

2

cholesky 分解