大学の課題で、コンパイルは通ったのですが、うまく読み込んでいないようです。
ガウス消去法を用いて Ax = B を x について解くというもので、行列Aと行列Bはそれぞれ.csvファイルに保存されています。
行列Aはコマンドライン関数として、行列Bはプログラム内で読み込むことがそれぞれ指定されています。
プログラムを実行してみると、行列A,B共に1列目の数字を読み込んだ際に値が変わってしまっているようです。その原因がどこにあるのかがわからずに困っています。
行列Aは4行4列の以下のようなもので、
{3.1 -5.3 3.8 -3.8}
{4.1 5.9 4.6 3.2}
{-5.9 7.9 -2.6 7.9}
{2.6 -3.2 -4.3 -5}
行列Bは1行4列の以下のようなものです。
{3}
{2}
{5}
{-1}
原因を教えていただけないでしょうか。よろしくお願いします!!
以下にプログラムと実行結果を載せます。
C
1#include <stdio.h> 2#include <stdlib.h> 3#include<string.h> 4#include <math.h> 5#define N 4 6#define EPS 1e-6 7typedef double Matrix[N][N+1]; 8 9 10//行列を読み込む関数 11void scan_data1(FILE* ifp, Matrix a) 12{ 13 int num[N] = { 0 }; //列の数を数える 14 for (int i = 0; i<N; i++) { 15 for (int j = 0; j<N + 1; j++) { 16 char data = 0; //文字を読み込むための 17 fscanf(ifp, "%c", &data); 18 if (data == '\n') { //改行していたら、ループを抜ける。 19 break; 20 } 21 else if (data == '\0') { //読み込みが終了したら、ループを抜ける。 22 break; 23 } 24 fscanf(ifp, "%lf", &a[i][j]); //数値を読み込む 25 } 26 } 27} 28 29void scan_data2(FILE* ifp, Matrix a) 30{ 31 for (int i = 0; i < N; i++) { 32 char data = 0; //文字を読み込むための 33 fscanf(ifp, "%c", &data); 34 if (data == '\n') { //改行していたら、ループを抜ける。 35 break; 36 } 37 else if (data == '\0') { //読み込みが終了したら、ループを抜ける。 38 break; 39 } 40 fscanf(ifp, "%lf", &a[i][N]); //数値を読み込む 41 } 42} 43 44//行列を表示する 45void print_matrix(Matrix a) 46{ 47 for (int i = 0; i<N; i++) { 48 for (int j = 0; j<N; j++) { 49 printf("%8.3f ", a[i][j]); 50 } 51 printf("| %8.3f\n", a[i][N]); 52 } 53 printf("\n"); 54} 55 56void change_low(Matrix a, int i) 57{ 58 int k = i; //最大となる行 59 double max_value = fabs(a[i][i]); 60 61 for (int j = i + 1; j<N; j++) { 62 if (max_value<fabs(a[j][i])) { 63 k = j; 64 max_value = fabs(a[j][i]); 65 } 66 } 67 if (max_value == fabs(a[i][i])) { //入れ替えがなかったら 68 printf("==== No change!! ====\n\n"); 69 } 70 else { 71 for (int j = 0; j<N + 1; j++) { 72 double tmp = a[i][j]; 73 a[i][j] = a[k][j]; 74 a[k][j] = tmp; 75 } 76 printf("==== change!! low:%d low:%d ====\n\n", i + 1, k + 1); 77 print_matrix(a); 78 } 79} 80 81void forward_elimination(Matrix a) 82{ 83 //前進消去 84 printf("====Foward elimination====\n"); 85 for (int i = 0; i<N; i++) { 86 change_low(a, i); 87 double pivot = a[i][i]; //ピボット 88 89 for (int j = i; j<N + 1; j++) 90 a[i][j] /= pivot; 91 92 //掃き出し計算 93 for (int k = i + 1; k<N; k++) { 94 if (k != i) { 95 double d = a[k][i]; 96 for (int j = i; j<N + 1; j++) //k行目をすべて揃えている 97 a[k][j] -= d*a[i][j] / a[i][i]; 98 } 99 } 100 print_matrix(a); 101 } 102} 103 104void backward_substitution(Matrix a) 105{ 106 //後退代入 107 printf("====Back substitution====\n"); 108 for (int i = N - 1; i >= 0; i--) { 109 for (int j = i + 1; j < N; j++) { 110 a[i][N] -= a[i][j] * a[j][N]; 111 a[i][j] = 0; 112 } 113 a[i][N] = a[i][N] / a[i][i]; 114 print_matrix(a); 115 } 116} 117 118// ガウス消去法を使って連立方程式を解く 119void gauss_elimination(Matrix a) 120{ 121 forward_elimination(a); 122 backward_substitution(a); 123} 124 125//計算結果を表示する 126void print_result(Matrix a) 127{ 128 printf("====results====\n"); 129 for (int i = 0; i<N; i++) { 130 printf(" x%d = %9.3lf\n", i + 1, a[i][N]); 131 } 132} 133 134int main(int argc, char* argv[]) 135{ 136 FILE* ifp1 = NULL; //入力用ファイルポインタ 137 char* ifile1 = NULL; //入力ファイル名 138 ifile1 = argv[1]; //パラメータの1番目を入力ファイル名 139 140 /*コマンドライン関数で読み込んだファイルを開く*/ 141 if (argc != 2) { //実行時引数が1個なかったらエラーとする 142 fprintf(stderr, "usage: %s inputfile outputfile\n", argv[0]); 143 exit(1); 144 } 145 if ((ifp1 = fopen(ifile1, "rt")) == NULL) { //入力ファイルを開く 146 fprintf(stderr, "Can't open file %s\a\n", ifile1); //開けなかったら終了 147 exit(2); 148 } 149 150 /*解の書き込まれているファイルを開く。*/ 151 FILE* ifp2 = NULL; 152 char* ifile2 = "b_vec.csv"; 153 if ((ifp2 = fopen(ifile2, "rt")) == NULL) { //入力ファイルを開く 154 fprintf(stderr, "Can't open file %s\a\n", ifile2); //開けなかったら終了 155 exit(3); 156 } 157 158 Matrix A = { 0 }; 159 160 scan_data1(ifp1, A); 161 print_matrix(A); //ファイルA確認 162 scan_data2(ifp2, A); 163 164 print_matrix(A); //ファイルAB確認 165 gauss_elimination(A); 166 print_result(A); 167 168 fclose(ifp1); //ファイルを閉じる 169 fclose(ifp2); 170 171 return (0); 172} 173
0.100 -5.300 3.800 -3.800 | 0.000 0.100 5.900 4.600 3.200 | 0.000 5.900 7.900 -2.600 7.900 | 0.000 0.600 -3.200 -4.300 -5.000 | 0.000 0.100 -5.300 3.800 -3.800 | 2.000 0.100 5.900 4.600 3.200 | 0.000 5.900 7.900 -2.600 7.900 | 0.000 0.600 -3.200 -4.300 -5.000 | 0.000 ====Foward elimination==== ==== change!! low:1 low:3 ==== 5.900 7.900 -2.600 7.900 | 0.000 0.100 5.900 4.600 3.200 | 0.000 0.100 -5.300 3.800 -3.800 | 2.000 0.600 -3.200 -4.300 -5.000 | 0.000 1.000 1.339 -0.441 1.339 | 0.000 0.000 5.766 4.644 3.066 | 0.000 0.000 -5.434 3.844 -3.934 | 2.000 0.000 -4.003 -4.036 -5.803 | 0.000 ==== No change!! ==== 1.000 1.339 -0.441 1.339 | 0.000 0.000 1.000 0.805 0.532 | 0.000 0.000 0.000 8.221 -1.044 | 2.000 0.000 0.000 -0.811 -3.675 | 0.000 ==== No change!! ==== 1.000 1.339 -0.441 1.339 | 0.000 0.000 1.000 0.805 0.532 | 0.000 0.000 0.000 1.000 -0.127 | 0.243 0.000 0.000 0.000 -3.778 | 0.197 ==== No change!! ==== 1.000 1.339 -0.441 1.339 | 0.000 0.000 1.000 0.805 0.532 | 0.000 0.000 0.000 1.000 -0.127 | 0.243 0.000 0.000 0.000 1.000 | -0.052 ====Back substitution==== 1.000 1.339 -0.441 1.339 | 0.000 0.000 1.000 0.805 0.532 | 0.000 0.000 0.000 1.000 -0.127 | 0.243 0.000 0.000 0.000 1.000 | -0.052 1.000 1.339 -0.441 1.339 | 0.000 0.000 1.000 0.805 0.532 | 0.000 0.000 0.000 1.000 0.000 | 0.237 0.000 0.000 0.000 1.000 | -0.052 1.000 1.339 -0.441 1.339 | 0.000 0.000 1.000 0.000 0.000 | -0.163 0.000 0.000 1.000 0.000 | 0.237 0.000 0.000 0.000 1.000 | -0.052 1.000 0.000 0.000 0.000 | 0.392 0.000 1.000 0.000 0.000 | -0.163 0.000 0.000 1.000 0.000 | 0.237 0.000 0.000 0.000 1.000 | -0.052 ====results==== x1 = 0.392 x2 = -0.163 x3 = 0.237 x4 = -0.052
回答3件
あなたの回答
tips
プレビュー
バッドをするには、ログインかつ
こちらの条件を満たす必要があります。
2017/06/09 05:39