LU分解を用いて、3次正方行列の逆行列の計算をする問題なのですが、分かりませんでした。
まず、上手くLとUに分解できたのですか、その後のLの逆行列とUの逆行列を出すことができませんでした。
そして、問題の答えである、Lの逆行列とUの逆行列の掛算がうまくできませんでした。
LとUの逆行列ができないので、それはそうなのですが、すみません。
LとUの分解を以下に載せておきます。
#include <stdio.h> #include <math.h> #define N 3 float A[N][N]={{2.0,1.0,1.0},{1.0,5.0,-1.0},{-1.0,1.0,-3.0}}; int main(void){ int i,j,k,p,q; float L[N][N]={0.0},U[N][N]={0.0}; float lu; for(j=0;j<N;j++) L[j][j]=1.0; for(j=0;j<N;j++) U[0][j]=A[0][j]; for(j=1;j<N;j++) L[j][0]=A[j][0]/U[0][0]; for(j=1;j<N-1;j++){ for(q=j;q<N;q++){ lu=0.0; for(k=0;k<j;k++) lu+=L[j][k]*U[k][q]; U[j][q]=A[j][q]-lu; } for(p=j+1;p<N;p++){ lu=0.0; for(k=0;k<j;k++) lu+=L[p][k]*U[k][j]; L[p][j]=(A[p][j]-lu)/U[j][j]; } } lu=0.0; for(k=0;k<N-1;k++)lu+=L[N-1][k]*U[k][N-1]; U[N-1][N-1]=A[N-1][N-1]-lu; for(i=0;i<N;i++){ printf("[ "); for(j=0;j<N;j++)printf("%9f ",A[j]); printf("]\n"); } printf("のLU分解は\n"); for(i=0;i<N;i++){ printf("[ "); for(j=0;j<N;j++)printf("%9f ",L[j]); printf("]\n"); } printf("\n"); for(i=0;i<N;i++){ printf("[ "); for(j=0;j<N;j++)printf("%9f ",U[j]); printf("]\n"); } return(0); }
同じような問題である、3次正方行列の方程式を解く問題は解けましたので、以下に載せておきます。
#include <stdio.h> #include <math.h> #define N 3 float A[N][N]={{2.0,1.0,1.0},{1.0,5.0,-1.0},{-1.0,1.0,-3.0}}; float B[N]={3.0,-6.0,-8.0}; int main(void){ int j,k,p,q; float L[N][N]={0.0},U[N][N]={0.0}; float lu; float x[N],y[N]; for(j=0;j<N;j++) L[j][j]=1.0; for(j=0;j<N;j++) U[0][j]=A[0][j]; for(j=1;j<N;j++) L[j][0]=A[j][0]/U[0][0]; for(j=1;j<N-1;j++){ for(q=j;q<N;q++){ lu=0.0; for(k=0;k<j;k++) lu+=L[j][k]*U[k][q]; U[j][q]=A[j][q]-lu; } for(p=j+1;p<N;p++){ lu=0.0; for(k=0;k<j;k++) lu+=L[p][k]*U[k][j]; L[p][j]=(A[p][j]-lu)/U[j][j]; } } lu=0.0; for(k=0;k<N-1;k++)lu+=L[N-1][k]*U[k][N-1]; U[N-1][N-1]=A[N-1][N-1]-lu; y[0]=B[0]; for(k=1;k<N;k++){ lu=0.0; for(j=0;j<k;j++)lu+=L[k][j]*y[j]; y[k]=B[k]-lu; } x[N-1]=y[N-1]/U[N-1][N-1]; for(k=N-2;k>=0;k--){ lu=0.0; for(j=N-1;j>k;j--)lu+=U[k][j]*x[j]; x[k]=(y[k]-lu)/U[k][k]; } printf("解:x=%f y=%f z=%f\n",x[0],x[1],x[2]); return(0); }
非常に長くなりましたが、暇のある人は挑戦してみてください。お願いします。