fork download
copy 
  1. #include <stdio.h>
  2.  
  3. //_連立方程式に表示
  4. void printEqs(
  5.   double  a[3][3],  //_係数行列
  6.   double  b[3])     //_右辺ベクトル
  7. {
  8.   int i,j;
  9.   printf("(a)[3][3],{b}[3])\n");
  10.   for( i=0;i<3;i++ ) {
  11.     for( j=0;j<3;j++ ) {
  12.       printf("_%9f",a[i][j]);
  13.     } printf("_%9f\n",b[i]);
  14.   }
  15. }
  16.  
  17. //_３元１次連立方程式のガウス消去
  18. int gausselim(
  19.   double  a[3][3],  //_係数行列
  20.   double  b[3])     //_右辺ベクトル→解
  21. {
  22.   double  Akk;      //_ピボット要素の保存値
  23.   double  Aik;      //_ピボット列要素の保存値
  24.   int     k;        //_ピボットインデクス
  25.   int     i;        //_行インデクス
  26.   int     j;        //_列インデクス
  27.   int     n=3;      //_方程式の数
  28.  
  29.   //_前進消去(k)
  30.   for( k=0;k<n;k++ ) {
  31.   	Akk = a[k][k];
  32.     //_デバッグ表示開始
  33.     printf("ピボット(k = %d)\n",k);
  34.     printEqs(a,b);
  35.     printf("pivot=a[%d][%d]=%f\n",k,k,Akk);
  36.     //_デバッグ表示終了
  37.     if( Akk == 0 ) return -1;
  38.     //_ピボット行正規化
  39.     for( j=k;j<n;j++ ) a[k][j] /= Akk;
  40.     b[k] /= Akk;
  41.     //_非ピボット行消去
  42.     for( i=k+1;i<n;i++ ) {
  43.       Aik = a[i][k];
  44.       for( j=k;j<n;j++ ) {
  45.         a[i][j] -= Aik*a[k][j];
  46.       } b[i] -= Aik*b[k];
  47.     }
  48.   }
  49.   //_後退代入(i)
  50.   for( i=n-2;i>=0;i-- ) {
  51.     for( j=i+1;j<n;j++ ) {
  52.       b[i] -= a[i][j]*b[j];
  53.     }
  54.   }
  55.   return 0;
  56. }
  57. /*
  58. 3x+1y+1z=8
  59. 6x+2y+3z=4
  60. 9x+4y+5z=12
  61. を解く。
  62. */
  63. int main(){
  64.   double  a[3][3]={{3.,1.,1.},{6.,2.,3.},{9.,4.,5.}};
  65.   double  b[3]={8.,4.,12.};
  66.   int     ret;
  67.   ret = gausselim(a,b);
  68.   if(ret==0) {
  69.      printf("{x}=%f,%f,%f\n",b[0],b[1],b[2]);
  70.   } else {
  71.     printf("ピボットがゼロ！、中止\n");
  72.     printf("ret = %d\n",ret);
  73.   }
  74.   return 0;
  75. }
Success #stdin #stdout 0s 5472KB
comments (?)
stdin
copy 
Standard input is empty
stdout
copy 
ピボット(k = 0)
(a)[3][3],{b}[3])
_ 3.000000_ 1.000000_ 1.000000_ 8.000000
_ 6.000000_ 2.000000_ 3.000000_ 4.000000
_ 9.000000_ 4.000000_ 5.000000_12.000000
pivot=a[0][0]=3.000000
ピボット(k = 1)
(a)[3][3],{b}[3])
_ 1.000000_ 0.333333_ 0.333333_ 2.666667
_ 0.000000_ 0.000000_ 1.000000_-12.000000
_ 0.000000_ 1.000000_ 2.000000_-12.000000
pivot=a[1][1]=0.000000
ピボットがゼロ！、中止
ret = -1
https://ideone.com/2sU3T3
language:
C (gcc 8.3)
created:
1 year 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!