緊急!今日が締切りの課題です。

フォーラム(掲示板)ルール
フォーラム(掲示板)ルールはこちら  ※コードを貼り付ける場合は [code][/code] で囲って下さい。詳しくはこちら
unittk501st
記事: 11
登録日時: 13年前

緊急!今日が締切りの課題です。

#1

投稿記事 by unittk501st » 12年前

キーボードから次元N を入力した後,N N 行列A とN 次元ベクトル⃗b のデータを要素
ごとにキーボードから入力し,その行列・ベクトルを使って,連立1 次方程式A⃗x = ⃗b を
部分ピボット選択付きガウスの消去法を用いて解け.
その計算課程も表示せよ.

表示は以下のようにせよ.
N=3 【Enter】
A[0][0] = 1 【Enter】
A[0][1] = 1 【Enter】
A[0][2] = 2 【Enter】
A[1][0] = 3 【Enter】
A[1][1] = 2 【Enter】
A[1][2] = -1 【Enter】
A[2][0] = 2 【Enter】
A[2][1] = -1 【Enter】
A[2][2] = 1 【Enter】
b[0] = 4 【Enter】
b[1] = -1 【Enter】
b[2] = 5 【Enter】
A, b =
1.0 1.0 2.0 | 4.0
3.0 2.0 -1.0 | -1.0
2.0 -1.0 1.0 | 5.0
column 0: row 1 is max
A, b =
3.0 2.0 -1.0 | -1.0
0.0 0.3 2.3 | 4.3
0.0 -2.3 1.7 | 5.7
column 1: row 2 is max
A, b =
3.0 2.0 -1.0 | -1.0
0.0 -2.3 1.7 | 5.7
0.0 -0.0 2.6 | 5.1
A, b =
3.0 2.0 -1.0 | -1.0
0.0 -2.3 1.7 | 5.7
0.0 -0.0 2.6 | 2.0
A, b =
3.0 2.0 -1.0 | -1.0
0.0 -2.3 1.7 | -1.0
0.0 -0.0 2.6 | 2.0
A, b =
3.0 2.0 -1.0 | 1.0
0.0 -2.3 1.7 | -1.0
0.0 -0.0 2.6 | 2.0

コード:

# include <stdio.h>
# include <stdlib.h>
# include <math.h>


void PrintAb(double *A, double *b, int n){
        int i, j;
	printf("\nA, b =\n ");
	for(i=0; i<n; i++){
	        for(j=0; j<n; j++){
			printf("%.2f\t", A[j]);
		}
		printf("|  %.2f\t\n", b[i]);
		A +=n;
		printf(" ");
	}
	printf("\n");

return;
} 


void swap(double *x, double *y){
	double z = *x;
	*x=*y;
	*y=z;
}

int main(void){
    int N, k, i, j, m;
	double *A, *ai, *ak, *b, r;

	/*N入力*/
	printf("N = ");
	scanf("%d", &N);
	printf("\n");
	
	/*A入力*/
	A=(double*)malloc(N*N*sizeof(double));
	if(A==NULL){
        printf("Can't allocate memory. \n");
        exit(1);
	}
	ai=A;
	for(i=0; i<N; i++){
		for(j=0; j<N; j++){
			printf("A[%d][%d] = ", i, j);
			scanf("%lf", ai+j);
			printf("\n");
		}
		ai +=N;
	}

	ai=A;

	/*b入力*/
	b=(double*)malloc(N*sizeof(double));
	if(b==NULL){
        printf("Can't allocate memory. \n");
        exit(1);
	}
	for(i=0; i<N; i++){
		printf("b[%d]    = ",i);
		scanf("%lf",b+i);
		printf("\n");
	}
		
	/*表示1*/
	printf("\n\n");
	
	PrintAb(A, b, N);

	
	/*ak定義*/
	ak=A;

	/*前進消去*/
	for (k=0; k<N-1; k++) {   /* k=0,...,N-2 */
		/* 絶対最大要素の検索 */
		ak+=N*k;
		m=k;  /*A[k][k]が最大と仮定*/
		for (i=k; i<N; i++){ /* i=k,...,N-1 */
			ai+=N*i;
			if( fabs(ak[k]) < fabs(ai[k]) ){
				m=i;
			}
			ai=A;
		}
		printf("column %d: row %d is max\n", k, m);

		/* m≠kのとき入れ替え */
		if(m!=k){
			ai+=N*m;
			swap(&ak[k], &ai[k]);
			swap(&b[k], &b[m]);
		}
		ai=A;
		ak=A;
		/*第k列消去*/
	  for (i=k+1; i<N; i++) {   /* i=k+1,...,N-1 */
	    ai+=N*i;
	    ak+=N*k;
	    r = *(ai+k)/ *(ak+k);
	    for (j=k; j<N; j++) {   /* j=k,...,N-1 */
	      *(ai+j) -= *(ak+j) * r;
	    }
	    b[i] -= b[k] * r;
	    ai=A;
	    ak=A;
	  }
	  PrintAb(A, b, N);
	}
	

	/*後退代入*/
	ai=A+(N * (N-1));
	for (i=N-1; i>=0; i--) { /* i=N-1,...,0 */
	  for (j=i+1; j<N; j++) { /* j=i+1,...,N-1 */
	    b[i]-=*(ai+j) * b[j];
	  }
	  b[i]=b[i]/ *(ai+i);
	  PrintAb(A, b, N);
	  ai-=N;
	}


	/*メモリ解放*/
	free(A);
	free(b);

return 0;
}
コンパイルして実行するとこうなってしまいます。

N = 3

A[0][0] = 1

A[0][1] = 1

A[0][2] = 2

A[1][0] = 3

A[1][1] = 2

A[1][2] = -1

A[2][0] = 2

A[2][1] = -1

A[2][2] = 1

b[0] = 4

b[1] = -1

b[2] = 5




A, b =
1.00 1.00  2.00 | 4.00
3.00 2.00  1.00 | -1.00
2.00 -1.00 1.00 | 5.00

column 0: row 2 is max

A, b =
2.00 1.00  2.00 | 5.00
0.00 0.50  -4.00 | -8.50
0.00 -1.50 0.00 | 1.50

column 1: row 2 is max

A, b =
2.00 1.00  2.00 | 5.00
0.00 -1.50 -4.00 | 1.50
0.00 0.00  -1.33 | -8.00


A, b =
2.00 1.00  2.00 | 5.00
0.00 -1.50 -4.00 | 1.50
0.00 0.00  -1.33 | 6.00


A, b =
2.00 1.00  2.00 | 5.00
0.00 -1.50 -4.00 | -17.00
0.00 0.00  -1.33 | 6.00


A, b =
2.00 1.00  2.00 | 5.00
0.00 -1.50 -4.00 | -17.00
0.00 0.00  -1.33 | 6.00

どこを正せばいいでしょうか。ご教示ください。宜しくお願いします。

beatle
記事: 1281
登録日時: 13年前
住所: 埼玉
連絡を取る:

Re: 緊急!今日が締切りの課題です。

#2

投稿記事 by beatle » 12年前

部分ピボット選択付きガウスの消去法を知らないのでプログラムを読んでもよくわからないのですが,
計算途中の行列の値などを表示させて,手計算した場合と異なっている場所を探せば何か分かるのではないでしょうか.

かずま

Re: 緊急!今日が締切りの課題です。

#3

投稿記事 by かずま » 12年前

unittk501st さんが書きました:

コード:

			swap(&ak[k], &ai[k]);
どこを正せばいいでしょうか。ご教示ください。宜しくお願いします。
ここかな?

コード:

			for (j = 0; j < N; j++) swap(&ak[j], &ai[j]);

unittk501st
記事: 11
登録日時: 13年前

Re: 緊急!今日が締切りの課題です。

#4

投稿記事 by unittk501st » 12年前

かずま さんが書きました:
unittk501st さんが書きました:

コード:

			swap(&ak[k], &ai[k]);
どこを正せばいいでしょうか。ご教示ください。宜しくお願いします。
ここかな?

コード:

			for (j = 0; j < N; j++) swap(&ak[j], &ai[j]);
swapはfor (k=0; k<N-1; k++)の中で回すので、再度繰り返すのはおかしいと思います。

かずま

Re: 緊急!今日が締切りの課題です。

#5

投稿記事 by かずま » 12年前

unittk501st さんが書きました:swapはfor (k=0; k<N-1; k++)の中で回すので、再度繰り返すのはおかしいと思います。
m行目と k行目全体を交換するんだから、(N+1) 個の列全部を swap しないといけません。
k は行番号。 j は列番号なので、繰り返しではありません。

コード:

#include <stdio.h>
#include <stdlib.h>
#include <math.h>
 
#define A(i,j) (a[(i)*n + (j)])

void PrintAb(double *a, double *b, int n)
{
    int i, j;
    printf("A, b =\n");
    for (i = 0; i < n; i++) {
        for (j = 0; j < n; j++)
            printf("%6.1f", A(i,j));
        printf("  | %6.1f\n", b[i]);
    }
}


void swap(double *x, double *y)
{
    double z = *x; *x = *y; *y = z;
}

int main(void)
{
    int n, k, i, j, m;
    double *a, *b, r;

    /*N入力 */
    printf("N = "); if (scanf("%d", &n) != 1) return 1;

    /*A入力 */
    a = malloc(n * n * sizeof(double));
    if (!a) printf("Can't allocate memory. \n"), exit(1);

    for (i = 0; i < n; i++)
        for (j = 0; j < n; j++) {
            printf("A[%d][%d] = ", i, j);
            if (scanf("%lf", &A(i,j)) != 1) return 1;
        }
    /*b入力 */
    b = malloc(n * sizeof(double));
    if (!b) printf("Can't allocate memory. \n"), exit(1);
    for (i = 0; i < n; i++) {
        printf("b[%d]    = ", i);
        if (scanf("%lf", &b[i]) != 1) return 1;
    }
    printf("\n");
    PrintAb(a, b, n); /*表示1 */
    for (k = 0; k < n - 1; k++) { /*前進消去 */
        m = k;   /* 絶対最大要素の検索, A[k][k]が最大と仮定 */
        for (i = k + 1; i < n; i++)
            if (fabs(A(m,k)) < fabs(A(i,k))) m = i;
        printf("column %d: row %d is max\n", k, m);

        if (m != k) { /* m != k のとき入れ替え */
            for (i = 0; i < n; i++) swap(&A(k,i), &A(m,i));
            swap(&b[k], &b[m]);
        }
        for (i = k + 1; i < n; i++) { /*第k列消去 */
            r = A(i,k) / A(k,k);
            for (j = k; j < n; j++)
                A(i,j) -= A(k,j) * r;
            b[i] -= b[k] * r;
        }
        PrintAb(a, b, n);
    }
    for (i = n - 1; i >= 0; i--) { /*後退代入 */
        for (j = i + 1; j < n; j++) /* j=i+1,...,N-1 */
            b[i] -= A(i,j) * b[j];
        b[i] = b[i] / A(i,i);
        PrintAb(a, b, n);
    }
    free(a); free(b); /*メモリ解放 */
    return 0;
}

閉鎖

“C言語何でも質問掲示板” へ戻る