fork download
copy 
  1. #include <stdio.h>
  2. #include <stdlib.h>
  3. #include <math.h>
  4.  
  5. #define N 4
  6.  
  7. // 4x4行列を表示する関数
  8. void printMatrix(float matrix[N][N]) {
  9.     for (int i = 0; i < N; i++) {
  10.         for (int j = 0; j < N; j++) {
  11.             printf("%10.4f ", matrix[i][j]);
  12.         }
  13.         printf("\n");
  14.     }
  15.     printf("\n");
  16. }
  17.  
  18. int main() {
  19.     // 1. データの初期化
  20.     float A_GE[N][N + 1] = {
  21.         {3.0,  1.5, -6.0,  4.8,  1.2},
  22.         {1.0,  1.5, -2.0, -2.4,  0.6},
  23.         {0.0, -1.5, -2.0, -1.0, -2.4},
  24.         {2.0,  4.0, -1.8, -0.6,  0.0}
  25.     };
  26.  
  27.     float A_GJ[N][N + 1] = {
  28.         {3.0,  1.5, -6.0,  4.8,  1.2},
  29.         {1.0,  1.5, -2.0, -2.4,  0.6},
  30.         {0.0, -1.5, -2.0, -1.0, -2.4},
  31.         {2.0,  4.0, -1.8, -0.6,  0.0}
  32.     };
  33.  
  34.     float A_Original[N][N] = {
  35.         {3.0,  1.5, -6.0,  4.8},
  36.         {1.0,  1.5, -2.0, -2.4},
  37.         {0.0, -1.5, -2.0, -1.0},
  38.         {2.0,  4.0, -1.8, -0.6}
  39.     };
  40.  
  41.     float x[N];
  42.     float ratio, sum;
  43.  
  44.     // ==========================================
  45.     // 2. ガウスの消去法 (Gaussian Elimination)
  46.     // ==========================================
  47.     for (int i = 0; i < N - 1; i++) {
  48.         for (int j = i + 1; j < N; j++) {
  49.             ratio = A_GE[j][i] / A_GE[i][i];
  50.             for (int k = 0; k < N + 1; k++) {
  51.                 A_GE[j][k] -= ratio * A_GE[i][k];
  52.             }
  53.         }
  54.     }
  55.  
  56.     // 後退代入 (Back-Substitution)
  57.     x[N - 1] = A_GE[N - 1][N] / A_GE[N - 1][N - 1];
  58.     for (int i = N - 2; i >= 0; i--) {
  59.         sum = 0;
  60.         for (int j = i + 1; j < N; j++) {
  61.             sum += A_GE[i][j] * x[j];
  62.         }
  63.         x[i] = (A_GE[i][N] - sum) / A_GE[i][i];
  64.     }
  65.  
  66.     printf("--- 1. ガウスの消去法による計算結果 ---\n");
  67.     for (int i = 0; i < N; i++) {
  68.         printf("x%d = %10.4f\n", i + 1, x[i]);
  69.     }
  70.     printf("\n");
  71.  
  72.     // ==========================================
  73.     // 3. ガウス・ジョルダン法 (Gauss-Jordan Method)
  74.     // ==========================================
  75.     for (int i = 0; i < N; i++) {
  76.         float pivot = A_GJ[i][i];
  77.         for (int j = 0; j < N + 1; j++) {
  78.             A_GJ[i][j] /= pivot;
  79.         }
  80.  
  81.         for (int k = 0; k < N; k++) {
  82.             if (k != i) {
  83.                 float factor = A_GJ[k][i];
  84.                 for (int j = 0; j < N + 1; j++) {
  85.                     A_GJ[k][j] -= factor * A_GJ[i][j];
  86.                 }
  87.             }
  88.         }
  89.     }
  90.  
  91.     // ガウス・ジョルダン法の結果を出力
  92.     printf("--- 2. ガウス・ジョルダン法による計算結果 ---\n");
  93.     for (int i = 0; i < N; i++) {
  94.         printf("x%d = %10.4f\n", i + 1, A_GJ[i][N]);
  95.     }
  96.     printf("\n");
  97.  
  98.     // ==========================================
  99.     // 4. 逆行列の計算 (ガウス・ジョルダン法を使用)
  100.     // ==========================================
  101.     float A_Inv[N][N] = {
  102.         {1, 0, 0, 0},
  103.         {0, 1, 0, 0},
  104.         {0, 0, 1, 0},
  105.         {0, 0, 0, 1}
  106.     };
  107.  
  108.     // Aの値を操作用の一時行列にコピー
  109.     float A_Temp[N][N];
  110.     for (int i = 0; i < N; i++) {
  111.         for (int j = 0; j < N; j++) {
  112.             A_Temp[i][j] = A_Original[i][j];
  113.         }
  114.     }
  115.  
  116.     for (int i = 0; i < N; i++) {
  117.         float pivot = A_Temp[i][i];
  118.         for (int j = 0; j < N; j++) {
  119.             A_Temp[i][j] /= pivot;
  120.             A_Inv[i][j] /= pivot;
  121.         }
  122.         for (int k = 0; k < N; k++) {
  123.             if (k != i) {
  124.                 float factor = A_Temp[k][i];
  125.                 for (int j = 0; j < N; j++) {
  126.                     A_Temp[k][j] -= factor * A_Temp[i][j];
  127.                     A_Inv[k][j] -= factor * A_Inv[i][j];
  128.                 }
  129.             }
  130.         }
  131.     }
  132.  
  133.     printf("--- 3. 算出された逆行列 (A^-1) ---\n");
  134.     printMatrix(A_Inv);
  135.  
  136.     // ==========================================
  137.     // 5. 検算 (A * A^-1)
  138.     // ==========================================
  139.     float Identity_Check[N][N] = {0};
  140.     for (int i = 0; i < N; i++) {
  141.         for (int j = 0; j < N; j++) {
  142.             for (int k = 0; k < N; k++) {
  143.                 Identity_Check[i][j] += A_Original[i][k] * A_Inv[k][j];
  144.             }
  145.         }
  146.     }
  147.  
  148.     printf("--- 4. 検算結果 (A * A^-1) ---\n");
  149.     printMatrix(Identity_Check);
  150.  
  151.     return 0;
  152. }
Success #stdin #stdout 0s 5304KB
comments (?)
stdin
copy 
Standard input is empty
stdout
copy 
--- 1. ガウスの消去法による計算結果 ---
x1 =    81.0521
x2 =   -29.4800
x3 =    27.0200
x4 =    -7.4200

--- 2. ガウス・ジョルダン法による計算結果 ---
x1 =    81.0521
x2 =   -29.4800
x3 =    27.0200
x4 =    -7.4200

--- 3. 算出された逆行列 (A^-1) ---
    5.5933    26.6200   -24.3200   -21.2000 
   -2.0667    -9.8000     8.8000     8.0000 
    1.7667     8.7000    -8.2000    -7.0000 
   -0.4333    -2.7000     2.2000     2.0000 

--- 4. 検算結果 (A * A^-1) ---
    1.0000     0.0000    -0.0000    -0.0000 
    0.0000     1.0000    -0.0000    -0.0000 
   -0.0000    -0.0000     1.0000     0.0000 
    0.0000     0.0000    -0.0000     1.0000
https://ideone.com/H4srS5
language:
C (gcc 8.3)
created:
23 hours ago
visibility:
public

Share or Embed source code

Discover > Sphere Engine API

The brand new service which powers Ideone!

Discover > IDE Widget

Widget for compiling and running the source code in a web browser!