C言語関係掲示板

過去ログ

No.910 部分選択法によるピボット選択

[戻る] [ホームページ]
No.11531

ピボット操作を加えたプログラムへの改良
投稿者---C/C++が発展途上の人(2004/01/06 17:36:54)


下のようなプログラムをピボット操作を加えて改良したいのですが、
どうすればいいのでしょうか?

/* Gaussの単純消去法 */
#include <stdio.h>

void main(void) {

  int n,i,j,k;
  double b, a[32][33];
  
   /* 行列の大きさをセット */
  printf("1次連立方程式の解を求めます。\n");
  printf("何元の連立方程式にしますか?\n n=");
  scanf("%d",&n);
  
   if(n > 32 || n < 1) {
    puts(" n の値は1〜32までです。\n");
    exit(1);
  }

  /* 方程式の入力 */
  printf("方程式の係数を入力してください。\n");
  for(k=0;k<n;k++)
    for(j=0;j<n+1;j++) {
      printf("a[%d][%d]=",k,j);
          scanf("%lf",&a[k][j]);
    }

  for(k=0;k<n;k++) {
    
    /* ピポットを1にする */
    b = a[k][k];
    for(j=0;j<n+1;j++) a[k][j] /= b;
    
    /* ピポット行よりも下の行について */
    for(i=k+1;i<n;i++){

      /* ピポットと同じ列を0に */
      b = -a[i][k];
      for(j=0;j<n+1;j++) a[i][j] += (b*a[k][j]);
    }
  }

  /* 次に後置代入 */
  for(k=n-2;k>=0;k--) {
    for(j=k+1;j<n;j++)
      a[k][n] -= (a[k][j] * a[j][n]);
  }
  
  /* 答を表示 */
  printf("Gauss単純消去法による解\n");
  for(i=0;i<n;i++)
    printf("x[%d]=%lf\n",i,a[i][n]);
}


No.11538

Re:ピボット操作を加えたプログラムへの改良
投稿者---たか(2004/01/07 09:33:11)


>下のようなプログラムをピボット操作を加えて改良したいのですが、
>どうすればいいのでしょうか?
>

ピボット選択には部分選択と完全選択があります。部分選択はガウスの
消去法第k段でk列の絶対値だけを比較して行を交換する方法、完全選択は
k列より右の列も比較して絶対値最大の列を交換する方法です。

部分選択は単純なのでちょっと元のプログラムに手を加えるだけで済み
ますが、完全選択は行でなく列を交換するので後から逆変換が必要です。
しかし精度的には完全選択の方が有利なので、科学技術計算には完全選択
がよく使われます。

ここでは改造が簡単な(それゆえ見通しが良い)部分選択法をご紹介しま
す。

/* Gaussの単純消去法(部分ピボット選択) */
#include <stdio.h>
#include <stdlib.h>
#include <math.h>

#define MAXN 32

int main(void)
{
  int n, i, j, k, l, imax;
  double b, t, a[MAXN][MAXN + 1];
  
   /* 行列の大きさをセット */
  printf("1次連立方程式の解を求めます。\n");
  printf("何元の連立方程式にしますか?\nn = ");
  scanf("%d", &n);
  
   if (n > MAXN || n < 1) {
    printf(" n の値は1〜%dまでです。\n", MAXN);
    exit(1);
  }

  /* 方程式の入力 */
  printf("方程式の係数を入力してください。\n");
  for (k = 0; k < n; k++)
    for (j = 0; j < n + 1; j++) {
      printf("a[%d][%d]=", k, j);
          scanf("%lf", &a[k][j]);
    }

  for (k = 0; k < n; k++) {
    /* ピポットを走査する */
    t = 0; imax = k;
    for (l = k; l < n; l++) {
      if (fabs(a[l][k]) > t) {
        imax = l;
        t = a[l][k];
      }
    }
    /* 見つかったピボット行をk行と入れ替える */
    for (l = 0; l < n + 1; l++) {
      t = a[k][l]; a[k][l] = a[imax][l]; a[imax][l] = t;
    }
    /* 後は同じ */
    b = a[k][k];
    for (j = 0; j < n + 1; j++) a[k][j] /= b;
    
    /* ピポット行よりも下の行について */
    for (i = k + 1; i < n; i++) {
      /* ピポットと同じ列を0に */
      b = -a[i][k];
      for (j = 0; j < n + 1; j++) a[i][j] += b * a[k][j];
    }
  }

  /* 次に後置代入 */
  for (k = n - 2; k >= 0; k--)
    for (j = k + 1; j < n; j++)
      a[k][n] -= a[k][j] * a[j][n];
  
  /* 答を表示 */
  printf("Gauss単純消去法(部分ピボット選択)による解\n");
  for (i = 0; i < n; i++)
    printf("x[%d]=%lf\n", i, a[i][n]);

  return 0;
}



No.11540

Re:ピボット操作を加えたプログラムへの改良
投稿者---C/C が発展途上の人(2004/01/07 13:30:07)


ありがとうございます。参考になりました。