#includestdio.h#includestdlib.h#defineN10//矩阵大小范围/**使用已经求出的x，向前计算x（供getx()调用）*floata[][]矩阵U*floatx[]方程组解*inti解的序号（数组X元素序号）*intn矩阵大小*turn公式中需要的和*/floatgetmx(floata[N][N],floatx[N],inti,intn){floatmx=0;intr;for(r=i+1;rn;r++){mx+=a[r]*x[r];}turnmx;}/**使用已经求出的y，向前计算y（供gety()调用）*floata[][]矩阵L*floaty[]数组Y*inti数组Y元素序号*intn矩阵大小*turn公式中需要的和*/floatgetmy(floata[N][N],floaty[N],inti,intn){floatmy=0;intr;for(r=0;rn;r++){if(i!=r)my+=a[r]*y[r];}turnmy;}/**解方程组，计算某x*floata[][]矩阵U*floatx[]方程组解*inti解的序号*intn矩阵大小*turn方程组的第i个解（数组X的第i个元素值）*/floatgetx(floata[N][N],floatb[N],floatx[N],inti,intn){floatsult;if(i==n-1)//计算最后一个x的值sult=(float)(b/a[n-1][n-1]);else//计算其他x值（对于公式中的求和部分，需要调用getmx()函数）sult=(float)((b-getmx(a,x,i,n))/a);turnsult;}/**解数组Y，计算其中一元素值*floata[][]矩阵L*floaty[]数组Y*inti数组Y元素序号*intn矩阵大小*turn数组Y的第i个元素值*/floatgety(floata[N][N],floatb[N],floaty[N],inti,intn){floatsult;if(i==0)//计算第一个y的值sult=float(b/a);else//计算其他y值（对于公式中的求和部分，需要调用getmy()函数）sult=float((b-getmy(a,y,i,n))/a);turnsult;}voidmain(){floatl[N][N]={0};//定义L矩阵floatu[N][N]={0};//定义U矩阵floaty[N]={0};//定义数组Yfloatx[N]={0};//定义数组Xfloata[N][N];//定义系数矩阵floatb[N];//定义右端项floatsum=0;inti,j,k;intn;intflag=1;while(flag){printf("请输入系数矩阵的大小：");scanf("%d",n);if(nN){printf("矩阵过大！\n");continue;}flag=0;}printf("请输入系数矩阵值：\n");for(i=0;in;i++){for(j=0;jn;j++){printf("a[%d][%d]:",i,j);scanf("%f",a[j]);}}printf("请输入右端项数组：\n");for(i=0;in;i++){printf("b[%d]:",i);scanf("%f",b);}/*显示原始矩阵*/printf("\n原始矩阵：\n");for(i=0;in;i++){for(j=0;jn;j++)printf("%0.3f",a[j]);printf("\n");}printf("\n\n");/*初始化矩阵l*/for(i=0;in;i++){for(j=0;jn;j++){if(i==j)l[j]=1;}}/*开始LU分解*//*第一步：对矩阵U的首行进行计算*/for(i=0;in;i++){u[0]=(float)(a[0]/l[0][0]);}/*第二步：逐步进行LU分解*/for(i=0;in-1;i++){/*对“L列”进行计算*/for(j=i+1;jn;j++){for(k=0,sum=0;kn;k++){if(k!=i)sum+=l[j][k]*u[k];}l[j]=(float)((a[j]-sum)/u);}/*对“U行”进行计算*/for(j=i+1;jn;j++){for(k=0,sum=0;kn;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;in;i++){for(j=0;jn;j++){printf("%0.3f",l[j]);}printf("\n");}/*输出矩阵u*/printf("\n矩阵U：\n");for(i=0;in;i++){for(j=0;jn;j++){printf("%0.3f",u[j]);}printf("\n");}/*回代方式计算数组Y*/for(i=0;in;i++){y=gety(l,b,y,i,n);}/*显示数组Y*/printf("\n\n数组Y：\n");for(i=0;in;i++){printf("y%d=%0.3f\n",i+1,y);}/*回代方式计算数组X*/for(i=n-1;i=0;i--){x=getx(u,y,x,i,n);}/*显示数组X*/printf("\n\n数组X：\n");for(i=0;in;i++){printf("x%d=%0.3f\n",i+1,x);}}
