求c語言寫的雙三次插值函數

void SPL(int n, double *x, double *y, int ni, double *xi, double *yi)； 是你所要。

10年積累的做網站、成都做網站經驗，可以快速應對客戶對網站的新想法和需求。提供各種問題對應的解決方案。讓選擇我們的客戶得到更好、更有力的網絡服務。我雖然不認識你，你也不認識我。但先網站設計后付款的網站建設流程，更有葉城免費網站建設讓你可以放心的選擇與我們合作。

已知 n 個點 x,y; x 必須已按順序排好。要插值 ni 點，橫坐標 xi[], 輸出 yi[]。

程序里用double 型，保證計算精度。

SPL調用現成的程序。

現成的程序很多。端點處理方法不同，結果會有不同。想同matlab比較，你需 嘗試 調用 spline（）函數 時，令 end1 為 1, 設 slope1 的值，令 end2 為 1 設 slope2 的值。

#include stdio.h

#include math.h

int spline (int n, int end1, int end2,

double slope1, double slope2,

double x[], double y[],

double b[], double c[], double d[],

int *iflag)

{

int nm1, ib, i, ascend;

double t;

nm1 = n - 1;

*iflag = 0;

if (n 2)

{ /* no possible interpolation */

*iflag = 1;

goto LeaveSpline;

}

ascend = 1;

for (i = 1; i n; ++i) if (x[i] = x[i-1]) ascend = 0;

if (!ascend)

{

*iflag = 2;

goto LeaveSpline;

}

if (n = 3)

{

d[0] = x[1] - x[0];

c[1] = (y[1] - y[0]) / d[0];

for (i = 1; i nm1; ++i)

{

d[i] = x[i+1] - x[i];

b[i] = 2.0 * (d[i-1] + d[i]);

c[i+1] = (y[i+1] - y[i]) / d[i];

c[i] = c[i+1] - c[i];

}

/* ---- Default End conditions */

b[0] = -d[0];

b[nm1] = -d[n-2];

c[0] = 0.0;

c[nm1] = 0.0;

if (n != 3)

{

c[0] = c[2] / (x[3] - x[1]) - c[1] / (x[2] - x[0]);

c[nm1] = c[n-2] / (x[nm1] - x[n-3]) - c[n-3] / (x[n-2] - x[n-4]);

c[0] = c[0] * d[0] * d[0] / (x[3] - x[0]);

c[nm1] = -c[nm1] * d[n-2] * d[n-2] / (x[nm1] - x[n-4]);

}

/* Alternative end conditions -- known slopes */

if (end1 == 1)

{

b[0] = 2.0 * (x[1] - x[0]);

c[0] = (y[1] - y[0]) / (x[1] - x[0]) - slope1;

}

if (end2 == 1)

{

b[nm1] = 2.0 * (x[nm1] - x[n-2]);

c[nm1] = slope2 - (y[nm1] - y[n-2]) / (x[nm1] - x[n-2]);

}

/* Forward elimination */

for (i = 1; i n; ++i)

{

t = d[i-1] / b[i-1];

b[i] = b[i] - t * d[i-1];

c[i] = c[i] - t * c[i-1];

}

/* Back substitution */

c[nm1] = c[nm1] / b[nm1];

for (ib = 0; ib nm1; ++ib)

{

i = n - ib - 2;

c[i] = (c[i] - d[i] * c[i+1]) / b[i];

}

b[nm1] = (y[nm1] - y[n-2]) / d[n-2] + d[n-2] * (c[n-2] + 2.0 * c[nm1]);

for (i = 0; i nm1; ++i)

{

b[i] = (y[i+1] - y[i]) / d[i] - d[i] * (c[i+1] + 2.0 * c[i]);

d[i] = (c[i+1] - c[i]) / d[i];

c[i] = 3.0 * c[i];

}

c[nm1] = 3.0 * c[nm1];

d[nm1] = d[n-2];

}

else

{

b[0] = (y[1] - y[0]) / (x[1] - x[0]);

c[0] = 0.0;

d[0] = 0.0;

b[1] = b[0];

c[1] = 0.0;

d[1] = 0.0;

}

LeaveSpline:

return 0;

}

double seval (int n, double u,

double x[], double y[],

double b[], double c[], double d[],

int *last)

{

int i, j, k;

double w;

i = *last;

if (i = n-1) i = 0;

if (i 0) i = 0;

if ((x[i] u) || (x[i+1] u))

{

i = 0;

j = n;

do

{

k = (i + j) / 2;

if (u x[k]) j = k;

if (u = x[k]) i = k;

}

while (j i+1);

}

*last = i;

w = u - x[i];

w = y[i] + w * (b[i] + w * (c[i] + w * d[i]));

return (w);

}

void SPL(int n, double *x, double *y, int ni, double *xi, double *yi)

{

double *b, *c, *d;

int iflag,last,i;

b = (double *) malloc(sizeof(double) * n);

c = (double *)malloc(sizeof(double) * n);

d = (double *)malloc(sizeof(double) * n);

if (!d) { printf("no enough memory for b,c,d\n");}

else {

spline (n,0,0,0,0,x,y,b,c,d,iflag);

if (iflag==0) printf("I got coef b,c,d now\n"); else printf("x not in order or other error\n");

for (i=0;ini;i++) yi[i] = seval(ni,xi[i],x,y,b,c,d,last);

free(b);free(c);free(d);

};

}

main(){

double x[6]={0.,1.,2.,3.,4.,5};

double y[6]={0.,0.5,2.0,1.6,0.5,0.0};

double u[8]={0.5,1,1.5,2,2.5,3,3.5,4};

double s[8];

int i;

SPL(6, x,y, 8, u, s);

for (i=0;i8;i++) printf("%lf %lf \n",u[i],s[i]);

return 0;

}

求用c語言編寫牛頓插值法

牛頓插值法：

#includestdio.h

#includealloc.h

float Language(float *x,float *y,float xx,int n)

{

int i,j;

float *a,yy=0.0;

a=(float *)malloc(n*sizeof(float));

for(i=0;i=n-1;i++)

{

a[i]=y[i];

for(j=0;j=n-1;j++)

if(j!=i)a[i]*=(xx-x[j])/(x[i]-x[j]);

yy+=a[i];

}

free(a);

return yy;

}

void main()

{

float x[4]={0.56160,0.5628,0.56401,0.56521};

float y[4]={0.82741,0.82659,0.82577,0.82495};

float xx=0.5635,yy;

float Language(float *,float *,float,int);

yy=Language(x,y,xx,4);

printf("x=%f,y=%f\n",xx,yy);

getchar();

}

2.牛頓插值法#includestdio.h

#includemath.h

#define N 4

void Difference(float *x,float *y,int n)

{

float *f;

int k,i;

f=(float *)malloc(n*sizeof(float));

for(k=1;k=n;k++)

{

f[0]=y[k];

for(i=0;ik;i++)

f[i+1]=(f[i]-y[i])/(x[k]-x[i]);

y[k]=f[k];

}

return;

}

main()

{

int i;

float varx=0.895,b;

float x[N+1]={0.4,0.55,0.65,0.8,0.9};

float y[N+1]={0.41075,0.57815,0.69675,0.88811,1.02652};

Difference(x,(float *)y,N);

b=y[N];

for(i=N-1;i=0;i--)b=b*(varx-x[i])+y[i];

printf("Nn(%f)=%f",varx,b);

getchar();

}

留下個郵箱，我發給你：牛頓插值法的程序設計與應用

高分懸賞 求三維數據點C語言插值計算程序

問題補充，因字數限制，挪到這

1.拉格朗日插值簡介：

對給定的n個插值節點x1,x2,…,xn,及其對應的函數值y1=f(x1), y2=f(x2),…, yn=f(xn);使用拉格朗日插值公式，計算在x點處的對應的函數值f(x);

2.一維拉格朗日插值c語言程序：

Int lagrange(x0, y0, n, x, y)

Float xo[], yo[], x;

Int n;

Float *y

{

Int i, j;

Float p;

*y=0;

If (n1)

{

For(i=0;in;i++)

{

P=1;

For(j=1;jn;j++)

{

If(i!=J)

P=p*(x-x0[j]/x0[i]-x0[j]);

}

*y=*y+p*y0[i];

Return(0);

}

Else

Return(-1);

}

3.例題。已知函數如下表所示，求x=0.472處的函數值：

X 0.46 0.47 0.48 0.49

Y 0.484655 0.4903745 0.502750 0.511668

計算這個問題的c語言程序如下：

#minclude stdio

#includeMnath.h

Main()

{

Float x0[4]={ 0.46, 0.47,0.48,0.49};

Float y0[4]={ 0.484655 ,0.4903745 ,0.502750 ,0.511668};

Float x, y;

Int n, rtn;

N=4;

X=0.472;

Rth=lagrange(x0,y0,n,x,y);

If(rtn=0)

{

Prinf(“Y(0.472)=:%f\n”,y);

}

Else

{

Prinf(“n must be larger than 1.\n”);

}

}

計算結果：Y(0.472)=:0.495553

4.問題補充

我的問題與上面的例子類似，計算三維空間一點(x,y,z)對應的函數值（Vx,Vy,Vz）.不同的是自變量(point_coordinate.txt)為三維空間散亂點（不是正方體的頂點），因變量(point_data.txt)為矢量（向量 ）。插值算法比較多，常數法，拉格朗日插值，埃特金插值，三階樣條插值等。最簡單的就是常數法，查找離目標點（x,y,z）距離最近的已知自變量(Xi,Yi,Zi)，把該點的函數值賦給目標點做函數值，求高手幫忙寫寫。

用C語言實現拉格朗日插值、牛頓插值、等距結點插值算法

#includestdio.h

#includestdlib.h

#includeiostream.h

typedef struct data

{

float x;

float y;

}Data;//變量x和函數值y的結構

Data d[20];//最多二十組數據

float f(int s,int t)//牛頓插值法，用以返回插商

{

if(t==s+1)

return (d[t].y-d[s].y)/(d[t].x-d[s].x);

else

return (f(s+1,t)-f(s,t-1))/(d[t].x-d[s].x);

}

float Newton(float x,int count)

{

int n;

while(1)

{

cout"請輸入n值(即n次插值):";//獲得插值次數

cinn;

if(n=count-1)// 插值次數不得大于count－1次

break;

else

system("cls");

}

//初始化t，y，yt。

float t=1.0;

float y=d[0].y;

float yt=0.0;

//計算y值

for(int j=1;j=n;j++)

{

t=(x-d[j-1].x)*t;

yt=f(0,j)*t;

//coutf(0,j)endl;

y=y+yt;

}

return y;

}

float lagrange(float x,int count)

{

float y=0.0;

for(int k=0;kcount;k++)//這兒默認為count－1次插值

{

float p=1.0;//初始化p

for(int j=0;jcount;j++)

{//計算p的值

if(k==j)continue;//判斷是否為同一個數

p=p*(x-d[j].x)/(d[k].x-d[j].x);

}

y=y+p*d[k].y;//求和

}

return y;//返回y的值

}

void main()

{

float x,y;

int count;

while(1)

{

cout"請輸入x[i],y[i]的組數，不得超過20組:";//要求用戶輸入數據組數

cincount;

if(count=20)

break;//檢查輸入的是否合法

system("cls");

}

//獲得各組數據

for(int i=0;icount;i++)

{

cout"請輸入第"i+1"組x的值:";

cind[i].x;

cout"請輸入第"i+1"組y的值:";

cind[i].y;

system("cls");

}

cout"請輸入x的值:";//獲得變量x的值

cinx;

while(1)

{

int choice=3;

cout"請您選擇使用哪種插值法計算："endl;

cout" (0):退出"endl;

cout" (1):Lagrange"endl;

cout" (2):Newton"endl;

cout"輸入你的選擇：";

cinchoice;//取得用戶的選擇項

if(choice==2)

{

cout"你選擇了牛頓插值計算方法，其結果為：";

y=Newton(x,count);break;//調用相應的處理函數

}

if(choice==1)

{

cout"你選擇了拉格朗日插值計算方法，其結果為：";

y=lagrange(x,count);break;//調用相應的處理函數

}

if(choice==0)

break;

system("cls");

cout"輸入錯誤!!!!"endl;

}

coutx" , "yendl;//輸出最終結果

}