LU和QR分解法解线性方程组(3)

void solution(double A[][N],double b[],int n); //求解上三角方程的解 void QR(double A[][N],double b[],int n);

void main() { double A[4][4]={{4,2,1,5},{8,7,2,10},{4,8,3,6},{12,6,11,20}}; double b[4]={-2,-7,-7,-3}; int i,j; int n=4; printf(\待求解的方程组系数矩阵A为：\\n\ for(i=0;i

void matrix_time(double A[][N],double B[][N],double C[][N],int n) //两个矩阵相乘，结果存在矩阵C[][N]

{ int i,j,k; for(i=0;i

void matrix_vec(double A[][N],double B[N],double C[N],int n) //矩阵和向量相乘,结果存在向量C[N]

{ int i,j; for(i=0;i

} }

double vec_value(double A[],int n) //求向量的模 { double Sum=0; int i; for(i=0;i

void vec_time(double a[],double H[][N],int n) //两个向量相乘得一个矩阵； { int i,j; for(i=0;i

void householder(double *a,double H[][N],int n, int m) //计算Householder矩阵 { int i,j; double first; //存放首元素 double vec_v; //存放向量的模 double a1[N]={0}; for(i=0;i

for(j=m;j

void matrix_turn(double A[][N],int n) //求矩阵的转置 { double a[N][N]={0}; int i,j; for(i=0;i

void solution(double A[][N],double b[],int n) //求解上三角方程的解 { int i,j; double x[N]={0}; double sum=0; for(i=n-1;i>=0;i--) { for(j=n-1;j>i;j--) sum+=A[i][j]*x[j]; x[i]=(b[i]-sum)/A[i][i]; sum=0; } printf(\矩阵的QR分解求解结果为：\\n\ for(i=0;i

void QR(double A[][N],double b[],int n) { int i,j,k,t; double H1[N][N]={0},H2[N][N]={0},H3[N][N]={0}; //H1,H2存放相邻两次的Householder矩阵，H3为中间变量矩阵

double temp[N][N]={0}; double Q[N][N]={0},R[N][N]={0},Q1[N][N]={0};

double b1[N]={0}; for(i=0;i

for(t=0;t

}

}

for(j=0;j

for(k=0;k

matrix_vec(Q,b,b1,n);//Q矩阵为正交阵，故Q的逆就等于Q装置，求解Q的逆

matrix_turn(Q,n); //转置，求解Q

printf(\printf(\系数矩阵A进行QR分解后的Q矩阵为：\\n\for(i=0;i

solution(R,b1,n); //求解上三角方程