Hermite cubic spline absolute error depended on amount of data points

Question

I'm trying to understand how does Hermite spline absolute error change depending on number of data points(n).

#include <stdio.h>
#include <math.h>
#include <stdlib.h>

double f(double x)
{
    return 5*pow(x,9)-19*pow(x,5)+5 ;
}
double fd(double x)
{
    return 45*pow(x,8)-95*pow(x,4);
}

double cubespline(double *x, double *y, double *d, double t, int n)
{
    double P,c0,c1,c2,c3,h;
    int i ;
    for(i=0;i<n;i++)
    {
        if(t>=x[i] && t<=x[i+1])
          {
             h=x[i+1]-x[i];
             c0=y[i];
             c1=d[i];
             c2=(3*y[i+1]-3*y[i]-d[i+1]*h-2*d[i])/(h*h);
             c3=(2*y[i]-2*y[i+1]+d[i+1]*h+d[i]*h)/(h*h*h);
             P=c0+c1*(t-x[i])+c2*(t-x[i])*(t-x[i])+c3*(t-x[i])*(t-x[i])*(t-x[i]);
          }
}
    return P;
}
int main(void)
{
    double *x, *y,*d, t=0.7,S,a=0,b=1;
    int n=2,i,k,iter=0,M=15;
    //printf("t=") ; scanf("%le",&t);
    //printf("n=") ; scanf("%d",&n);
   
    x=(double*)malloc(n*2*pow(2,M)*sizeof(double));
    y=(double*)malloc(n*2*pow(2,M)*sizeof(double));
    d=(double*)malloc(n*2*pow(2,M)*sizeof(double));
   
    while(iter<M)
{
    iter++ ;
    printf("iter=%d\n",iter) ;
   
    x[0]=a ; x[n]=b ;
    for(i=1;i<n;i++)
    {
        x[i]=a+((b-a)*i)/n;  
        y[i]=f(x[i]);
        d[i]=fd(x[i]);
    }
   
    S=cubespline(x,y,d,t,n) ;
     
    printf("n=%d\n",n) ;
    printf("Error=%le\n",fabs( f(t)-S) )  ;
    printf("\n\n") ;
    n*=2 ;
}
    return 0;
}

It seems like error increases as we increase n.For example:

    n=2 Error=1.360171
    n=4 Error=5.530201
    n=8 Error=8.471831
    n=16 Error=1.42389
    n=32 Error=5.883009
    ...
    n=32678 Error=1.455327e+01

The same results i get when i compute maximum of difference between splines with n+1,2n+1 and 4n+1 points.I do this to calculate error when we don't know the function so we can't calculate difference between f(x) and spline(x).If the difference between splines decreases as we increase n,it means that error decreases as we increase n.But again,i get that the difference between splines also increases as we increase n

#include <stdio.h>
#include <math.h>
#include <stdlib.h>

double f(double x)
{
    return sin(x) ;
}

double fd(double x)
{
    return cos(x);
}


double cubespline(double *x, double *fx, double *dx, double t,int n )
{
    double P,c0,c1,c2,c3,h;
    int i ;
   
     for(i=1 ; i<n-1 ; i++) if(t>=x[i]) break ;
     //printf("%d    %le   %le    %le\n",i,x[i],fx[i],dx[i]) ;
     h=x[i+1]-x[i];
     c0=fx[i];
     c1=dx[i];
     c2=(3*fx[i+1]-3*fx[i]-dx[i+1]*h-2*dx[i])/(h*h);
     c3=(2*fx[i]-2*fx[i+1]+dx[i+1]*h+dx[i]*h)/(h*h*h);
     P=c0+c1*(t-x[i])+c2*(t-x[i])*(t-x[i])+c3*(t-x[i])*(t-x[i])*(t-x[i]);
     
     return P;
}

int main(void)
{
    double *x, *fx,*dx,*y,*fy,*dy,*z,*fz,*dz,t,r=0,maxr=0,a=0,b=2*M_PI ;
    int n=128,i;
    //printf("t=") ; scanf("%le",&t);
    //printf("n=") ; scanf("%d",&n);
   
    x=(double*)malloc((4*n+1)*sizeof(double));
    fx=(double*)malloc((4*n+1)*sizeof(double));
    dx=(double*)malloc((4*n+1)*sizeof(double));
   
    y=(double*)malloc((n*2+1)*sizeof(double));
    fy=(double*)malloc((n*2+1)*sizeof(double));
    dy=(double*)malloc((n*2+1)*sizeof(double));
   
    z=(double*)malloc((n+1)*sizeof(double));
    fz=(double*)malloc((n+1)*sizeof(double));
    dz=(double*)malloc((n+1)*sizeof(double));
   
    for(i=0;i<=4*n;i++)
    {
x[i]=a+((b-a)*i)/(4*n);  
        fx[i]=f(x[i]);
        dx[i]=fd(x[i]);
        //printf("%d    %le   %le    %le\n",i,x[i],fx[i],dx[i]) ;
    }
   
    for(i=0;i<=2*n;i++)
    {
y[i]=a+((b-a)*i)/(2*n);  
        fy[i]=f(y[i]);
        dy[i]=fd(y[i]);
        //printf("%d    %le   %le    %le\n",i,y[i],fy[i],dy[i]) ;  

    }
    for(i=0;i<=n;i++)
    {
z[i]=a+((b-a)*i)/n;  
        fz[i]=f(z[i]);
        dz[i]=fd(z[i]);
        //printf("%d    %le   %le    %le\n",i,z[i],fz[i],dz[i]) ;
    }
   
    maxr=fabs(cubespline(x,fx,dx,x[0],4*n+1)-cubespline(y,fy,dy,x[0],2*n+1) ) ;
    for(i=1 ; i<=4*n ; i++)
   {
       r=fabs(cubespline(x,fx,dx,x[i],4*n+1)-cubespline(y,fy,dy,x[i],2*n+1) ) ;
       if(r>maxr) maxr=r ;
   }
   
     printf("Max1=%le\n",maxr) ;
     
    maxr=fabs(cubespline(y,fy,dy,y[0],2*n+1)-cubespline(z,fz,dz,y[0],n+1) ) ;
    for(i=1 ; i<=2*n ; i++)  
   {
       r=fabs(cubespline(y,fy,dy,y[i],2*n+1)-cubespline(z,fz,dz,y[i],n+1) ) ;
       if(r>maxr) maxr=r ;
   }
   
     printf("Max2=%le\n",maxr) ;
    return 0;
}

Shouldn't error decrease as we increase n?

Answer 1

After having a look at the code I don't seem to find any implementation error, just beware that an estimate of error for such a spline is not directly dependent on the number of initial points given to the algorithm:

where h is the maximum of step size so your test should be developing accordingly.