杜立特分解法.doc

2013-2014（1）专业课程实践论文 题目：杜立特分解法 一、算法理论 阶线性方程组的系数矩阵非奇异且有分解式，其中为单位下三角矩阵，为上三角矩阵，即，当时，；，当时，，矩阵的这种分解方法为Doolittle的分解。 比较等号两边的第行和第列的元素，得。 因为，所以，， 从而当，，时，，从而， 于是就得到了计算LR分解的一般计算公式。 二、算法框图 三、算法程序 #includeiostream #includestdlib.h #define N 3 using namespace std; int main() { double A[N+1][N+1]={{0,0,0,0},{0,2,1,1},{0,1,3,2},{0,1,2,2}}; double L[N+1][N+1]={0}; double U[N+1][N+1]={0}; double b[N+1]={0,4,6,5}; double y[N+1]; double x[N+1]; int i,j,k,p; for(j=1;j=N;j++) U[1][j]=A[1][j]; L[1][1]=1; for(i=2;i=N;i++) {L[i][1]=A[i][1]/U[1][1]; L[i][i]=1; } for(k=2;k=N;k++) { for(j=k;j=N;j++) {double l=0; for(p=1;p=k-1;p++) l+=L[k][p]*U[p][j]; U[k][j]=A[k][j]-l; } for(i=k+1;i=N;i++) {double l=0; for(p=1;p=k-1;p++) l+=L[i][p]*U[p][k]; L[i][k]=(A[i][k]-l)/U[k][k]; } } y[1]=b[1]; for(k=2;k=N;k++) { double l=0; for(j=1;j=k-1;j++) l+=L[k][j]*y[j]; y[k]=b[k]-l; } x[N]=y[N]/U[N][N]; for(k=N-1;k=1;k--) { double l=0; for(j=k+1;j=N;j++) l+=U[k][j]*x[j]; x[k]=(y[k]-l)/U[k][k]; } printf(向量y为:); for(i=1;i=N;i++)printf(%.1lf\t,y[i]);printf(\n); coutL：endl; for(i=1;i=N;i++) {for(j=1;j=N;j++) coutL[i][j]\t; coutendl;} coutendl; coutU：endl; for(i=1;i=N;i++) {for(j=1;j=N;j++) coutU[i][j]\t; coutendl;} coutendl; cout方程组的解：endl; for(i=1;i=N;i++) coutxi=x[i]\t; return 1; } 四、算法实现 例1．用杜立特分解法，求解三元方程组。 解： 例2．用杜立特分解法，求解三元方程组 解： 开始 输入矩阵行数 输入矩阵A第各行元素 分解矩阵A得L 输入矩阵B第各行元素 分解矩阵A得R 求得Y 求得X