1次元ポアソン方程式を差分法で解くプログラム

フォーラム(掲示板)ルール
フォーラム(掲示板)ルールはこちら  ※コードを貼り付ける場合は [code][/code] で囲って下さい。詳しくはこちら
兎汰李

1次元ポアソン方程式を差分法で解くプログラム

#1

投稿記事 by 兎汰李 » 14年前

1次元ポアソン方程式を差分法で解くプログラムをつくりたいです。
以下のようなプログラムを作成しました。これで、「1次元ポアソン方程式を差分法で解くプログラム」は完成でよいのでしょうか?
足りない部分や間違いの指摘お願いします。

コード:

#include<stdio.h>
double a[10][10],b[10]; double w[10];
int N=0; int i=0; int j=0; int k=0;  double s,p;
void decomp(void);
void solver(void);
void sakusei(int);
void toku(double[][10],int);

void sakusei(int N){
	for(i=1;i<=N;i++){
		for(j=1;j<=N;j++){
			a[i][j]=0;
		}
	}
	a[1][1]=2;
	a[1][2]=-2;
	for(i=2;i<=N-1;i++){
		a[i][i-1]=-1;
		a[i][i]=2;
		a[i][i+1]=-1;
	}
	a[N][N-1]=-1;
	a[N][N]=2;
	printf("---左辺ベクトル---\n");
	for(i=1;i<=N;i++){
		for(j=1;j<=N;j++){
			printf("tk[%d][%d]=%4.2f  ",i,j,a[i][j]);
		}
		printf("\n");
	}
	for(i=1;i<=N;i++){
		b[i]=(1.0/N)*(1.0/N);
	}
	printf("---右辺ベクトル---\n");
	for(i=1;i<=N;i++){
		printf("tf[%d]=%4.2f \n", i,b[i]);
	}
}
void choku(double[][10],int){
//tf[i],tu[i]の初期化
	int i=0;
	for(i=0;i<N;i++){
		b[i]=0;
		w[i]=0;
	}
	for(i=0;i<=N;i++){
		b[i]=(1.0/N)*(1.0/N);
	}
	
	//連立方程式の呼び出し
	decomp();
	solver();
}
//前進消去
void decomp(void)
{
	int i=0,j=0,k=0;
	double s=0,p=0;
	for(i=1;i<=N;i++){
		s=a[i][i];
		b[i]=b[i]/s;
		for(j=1;j<=N;j++){
			a[i][j]=a[i][j]/s;
		}
		for(k=i+1;k<=N;k++){
			p=a[k][i];
			for(j=1;j<=N;j++){
				a[k][j]=a[k][j]-p*a[i][j];
			}
			b[k]=b[k]-p*b[i];
		}
	}
}

//後退代入
void solver(void)
{
	
	double s[10];
	double sum=0;
	for(i=0;i<10;i++){
		s[i]=0;
	}
	for(i=N-1;i>0;i--){
		for(j=N;j>i;j--){
			s[j]=a[i][j]*b[j];
			sum+=s[j];
			a[i][j]=0;
		}
	b[i]=b[i]-sum;
	sum=0;
	}
	for(i=1;i<=N;i++){
		w[i]=b[i];
	}
}



	
int main(void){
	printf("---入力---\n");
	printf("刻み幅nn=");
	scanf("%d",&N);
	sakusei(N);
	choku(a,N);
	printf("---解---\n");
	for(i=1;i<=N;i++){
		printf("w[%d]=%4.2f \n", i,w[i]);
	}
	return 0;
}
 

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