LU分解

解：首先求出Q矩阵，计算过程如下。 2.下一个分块矩阵的消元变换。 =L*U

1. L矩阵第k列变换

由此可以分离出L和U如下：

1.L矩阵的第k列变换

将L矩阵这一列的全部元素除以Q矩阵这一列的对角线元素。 2.下一个分块矩阵的消元变换

分块矩阵的每一行元素减去（同一行最近的非分块矩阵元素与同一列最近的非分块矩阵元素的乘积）

程序4.1 矩阵LU分解的C语言程序 Int lluu(n,a,l,u) Int n //矩阵A的阶数

double a[n][n] //存放n阶实矩阵A；返回时存放Q矩阵 double l[n][n] //返回下三角矩阵L double u[n][n] //返回上三角矩阵U

int lluu() //函数返回一个整型标志。若返回的标志值为0，则表

示程序失败，输出信息fail；若不为0，则表示正常

返

回。

#include #include Int lluu(n,a,l,u) Int n;

double a[],l[],u[]; {

int i,j,w,v,ll;//ll是变量不是位数

for(k=0;k<=n-2;k++) //计算Q矩阵。分析序型

{ll=k*n+k;

if(fab(a[ll]+1.0==1.0) //主对角元素为0，错误返回 ,printf(“0 state fail\\n”);return 0;-

for(i=k+1;i<=n-1;i++) //计算L矩阵第k列元素。分析序型

{w=i*n+k;a[w]=a[w]/a[ll];}

for(i=k+1;i<=n-1;i++) //下一个分块矩阵消元变换。分析序型 {w=i*n+k;

for(j=k+1;j<=n-1;j++) { v=i*n+j;

a[v]=a[v]-a[w]*a[k*n+j]; } } }

for(i=0;i<=n-1;i++) //分析序型 {for(j=0;j{w=i*n+j;l[w]=a[w];u[w]=0.0;} w=i*n+I;

l[w]=1.0;u[w]=a[w]; for(j=i+1;j<=n-1;j++)

{w=i*n+j;l[w]=0.0;u{w}=a[w];} } return 1; } main() { int I,j;

static double l[4][4],u[4][4];

static double a[4][4]={{2.0,4.0,4.0,2.0}, {3.0,3.0,12.0,6.0}, {2.0,4.0,-1.0,2.0}, {4.0,2.0,1.0,1.0}}; i=lluu(4,a,1,u); if(i!=0)

,printf(“MAT L IS:\\n”); for(i=0;i<=3;i++) {for(j=0;j<=3;j++) printf(“%13.7e”,l*i+*j+); printf(“\\n”); }

printf(“\\n”); printf(“MAT U IS:\\n”); for(i=0;i<=3;i++) {for(j=0;j<=3;j++) printf(“%13.7e”,u[i][j]); printf(“\\n”); } printf(“\\n”); }

}