10 #define I_LENGTH (I_END - I_START)
13 int main(int argc,char **argv) {
17 double x,x_i,x_j,x_k,w;
23 printf("usage: %s <N>\n",argv[0]);
28 step=(double)I_LENGTH/(N+1);
31 strcpy(file,"polynom_interpolation.plot");
34 /* calculate fixpoints */
35 buf=(double *)malloc((N+2)*sizeof(double));
38 buf[i]=(sin(x)-x*cos(x))/(x*x*x*x);
39 /* manually, gp_api sux! */
40 dprintf(fd,"%f %f\n",x,buf[i]);
43 /* do interpolation */
45 buf_i=(double *)malloc((N+1)*N_I*sizeof(double));
46 for(i=0;i<(N+1);i++) {
47 for(l=1;l<N_I+1;l++) {
49 x_i=I_START+i*step+l*step_i;
50 for(j=0;j<(N+2);j++) {
53 for(k=0;k<(N+2);k++) {
55 if(j!=k) w=w*((x_i-x_k)/(x_j-x_k));
57 buf_i[i*(N_I+1)+l]+=w*buf[j];
59 dprintf(fd,"%f %f\n",x_i,buf_i[i*(N_I+1)+l]);