三次样条插值算法

#include
#include
#define N 19
#define X 12
int Find（double *adouble x）;
double calcudet（double a[N][N]）;
int main（void）
{
	int ijpqk;
	int x[]={24612163060110180280400515};
	double a[N]b[N]h[N]u[N]l[N]g[N];
	double y[X];
	double b_[]={1.86548-0.046115};
	double A[N][N]={0}B[N][N]={0}C[N][N]={0};
	double M[N];
	double denonume;
	FILE *fp;
	
	/*从文件读取数据到数组*/
	fp=fopen（“Points.txt““r“）;
	for（i=0;i		fscanf（fp“%lf%lf“&a[i]&b[i]）;
	fclose（fp）;
	
	/*计算所需要的数据*/
	for（i=1;i		h[i]=a[i]-a[i-1];
	for（i=1;i		u[i]=h[i]/（h[i]+h[i+1]）;
	for（i=1;i	    l[i]=1-u[i];
	for（i=1;i		g[i]=6*（（b[i+1]-b[i]）/h[i+1]-（b[i]-b[i-1]）/h[i]）/（h[i]+h[i+1]）;
	g[0]=6*（（b[1]-b[0]）/h[1]-b_[0]）/h[1];
	g[N-1]=6*（b_[1]-（b[N-1]-b[N-2]）/h[N-1]）/h[N-1];
	
	/*系数行列式的初始化*/
	for（i=0;i	    for（j=0;j	        A[i][j]=0;
	for（i=0;i	    A[i][i]=2;
	A[0][1]=A[N-1][N-2]=1;
	for（i=1;i		A[i][i+1]=l[i];
		A[i][i-1]=u[i];
	}
	for（i=0;i	    for（j=0;j	        C[i][j]=A[i][j];
	
	/*Cramer法则计算线性方程组*/
	deno=calcudet（A）;
	for（j=0;j		for（p=0;p	        for（q=0;q	            B[p][q]=C[p][q];
		for（i=0;i		    B[i][j]=g[i];
		nume=calcudet（B）;
		M[j]=nume/deno;
	}
	
	/*计算S（x）*/
	for（i=0;i		k=Find（ax[i]）;
		y[i]=M[k-1]*pow（（a