#include"stdio.h"
#include"malloc.h"
#include"math.h"
#include"string.h"
main(){
int i,j,n;
printf("請輸入需計算中心差分階數:");
scanf("%d",&j);
while (j==0||j%2==1){
printf("中心差分階數不可為0或奇數\n請重新輸入:");
scanf("%d",&j);
}
i=j/2;
n=j+1;
float a;//為系數矩陣分配空間
a=(float)malloc(n*sizeof(float*));
a[0]=(float*)malloc(n*n*sizeof(float));
int c,d,e,h;
for(c=1;c
a[c]=a[0]+n*c;
}
float *b;//為常數分配空間
b=(float*)malloc(n*sizeof(float));
float *g;//為權系數分配空間
g=(float*)malloc(n*sizeof(float));
memset(g,0,n*sizeof(float));
e=1;
float f;
for(c=0;c
for(d=0;d
if(c==0){
a[d][c]=1.0;
}
else{
h=d-i;
a[d][c]=pow(h,c);
f=(float)e;
a[d][c]=a[d][c]/e;
}
}
e=e*(c+1);
}
for(c=0;c
if(c==1)
b[c]=1.0;
else
b[c]=0.0;
}
for(c=0;c
for(d=c+1;d
if(a[c][c]==0){
for(e=c;e
f=a[e][c];
a[e][c]=a[e][d];
a[e][d]=f;
}
f=b[c];
b[c]=b[d];
b[d]=f;
}
if(a[c][c]!=0)
break;
}
for(d=c+1;d
a[d][c]/=a[c][c];
b[c]/=a[c][c];
a[c][c]=1;
}
for(d=c+1;d
for(e=c;e
a[e][d]=a[e][d]-a[c][d]*a[e][c];
}
b[d]=b[d]-a[c][d]*b[c];
}
}//
for(c=n-1;c>=0;c--){//計算權系數的值
g[c]=b[c];
if(c>0)
b[c-1]=b[c-1]-a[c][c-1]*g[c];
}
for(c=0;c<n;c++){//打印權系數的值
d=c-i;
printf("w%d=%f(1/dx)\n",d,g[c]);
}
free(a[0]);
free(a);
free(b);
free(g);
}
-nan就是負的無窮大,檢查下公式用得對不對,除數是否接近0。還有整數不能作為除數、乘數,需要先轉換成浮點。