crout分解法1Word文档格式.docx
《crout分解法1Word文档格式.docx》由会员分享，可在线阅读，更多相关《crout分解法1Word文档格式.docx（11页珍藏版）》请在冰点文库上搜索。
=li1（i=1,2,·
n）
a1j=（l11,0,0,0·
=l11·
u1j（j=1,2,·
n）
得：
li1=ai1;
u1j=a1j/l11（i,j=1,2,·
2当i<
j时
aij=（l11,li2,·
lii,0·
=
当i
得：
lij=aij-
（j=2,3,·
i，i=2,3,·
uij=（aij-
）/lii（i=2,3,·
n，j=i+1,i+2,·
（2）具体算法如下:
①对A作LU分解：
由A=LU及矩阵的乘法原理可得:
lij=aij-
（j=1,2,…,i,i=1,2,…n）
uij=（aij-
）/lii（j=i+1,i+2,…,n，i=1,2,…n）
②解两个三角型方程组
由A=LU及Ax=b可得L（Ux）=b，令Y=（y1,y2,…,yn）T=Ux，则LY=b,于是求解Ax=b就被化为求解下三角型方程组LY=b及单位上三角型方程组Ux=Y。
a）先解下三角型方程LY=b。
由l11y1=b1,l21y1+l22y2=b2l11,……ln1y1+ln2y2+…+lnnyn=bn
所以y1=b1/l11,yi=（bi-
）/lii,i=2,3,…,n
b）再解单位上三角型方程组Ux=Y。
由UX=Y得x1+u12x2+…+u1nxn=y1，x2+…+u2nxn=y2
………xn-1+un-1nxn=yn-1，xn=yn，利用回代解法可得方程组AX=b的解为
xn=yn,xi=yi-
i=n-1,…,2,1
二、Crout方法解线性方程组的程序
（1）程序代码：
#include"
stdio.h"
math.h"
//头文件
#defineN20//自定义N=20
intmain（）//主函数
{
inti,j,k;
intsize;
floata[N][N],l[N][N],u[N][N];
floatb[N],x[N],y[N];
//定义变量
printf（"
\t\t\tCrout分解法解方程组\n"
）;
请输入方阵A的n：
"
scanf（"
%d"
&
size）;
\n"
请输入方程组的系数:
for（i=0;
i<
size;
i++）
{
for（j=0;
j<
j++）
{
scanf（"
%f"
a[i][j]）;
//输入方程组系数矩阵a[][]
}
}
\n请输入方程组的y:
scanf（"
b[i]）;
//输入结果矩阵b[]
\n方阵A[][]为：
printf（"
%f"
a[i][j]）;
//输出a[][]
printf（"
\n方程组y为:
b[i]）;
//输出b[]
u[i][i]=1;
//定初始值令u[i][i]=1
for（j=i+1;
l[i][j]=0;
//定初始值令l[i][j]=0
for（j=0;
for（i=j+1;
u[i][j]=0;
//定初始值令u[i][j]=0
l[0][0]=a[0][0];
for（i=1;
l[i][0]=a[i][0];
//计算第一行的l[][]
u[0][i]=a[0][i]/l[0][0];
//计算第一列的u[][]
size-1;
for（j=1;
=i;
j++）//计算第2行到第size-1行的l[][]
l[i][j]=a[i][j];
for（k=0;
k<
j;
k++）
{
l[i][j]=l[i][j]-l[i][k]*u[k][j];
}
j++）//计算第2行到第size行的u[][]
u[i][j]=a[i][j];
=i-1;
u[i][j]=u[i][j]-l[i][k]*u[k][j];
u[i][j]=u[i][j]/l[i][i];
for（j=1;
j++）//计算第size行的l[][]
l[size-1][j]=a[size-1][j];
for（k=0;
=j-1;
l[size-1][j]=l[size-1][j]-l[size-1][k]*u[k][j];
输出矩阵L[i][j]\n"
i++）{
j++）
l[i][j]）;
"
//输出下三角矩阵l[][]
输出矩阵U[i][j]\n"
u[i][j]）;
//输出单位上三角矩阵u[][]
y[0]=b[0]/l[0][0];
//给y[0]初始值
i++）//计算y[i]的值
y[i]=b[i];
y[i]=y[i]-l[i][k]*y[k];
//计算公式
y[i]=y[i]/l[i][i];
y值:
y[%d]=%f"
i+1,y[i]）;
//输出y[i]的结果
\n\n"
x的值:
x[size-1]=y[size-1];
//给x[size-1]赋值
for（i=size-2;
i>
=0;
i--）
x[i]=y[i];
for（k=i+1;
x[i]=x[i]-u[i][k]*x[k];
//计算x[]
x[%d]=%f\n"
i+1,x[i]）;
//输出x[i]的结果
}
（2）程序应用
升级会员 