三次样条插值—三弯矩法（C语言实现）

目录

1、概述

2、案例分析 

（1）案例

（2）分析 

 3、代码及结果

---

1、概述

2、案例分析 

（1）案例

（2）分析 

 

 

 3、代码及结果

#include<stdio.h> 
#include<math.h> 
#define N 9 
void main() 
{ 
	double x[N+1]={1,2,3,4,5,6,7,8,9,10}, 
		y[N+1]={0,0.69314718,1.0986123,1.3862944,1.6094378,1.7917595,1. 9459101,2.079445,2.1972246,2.3025851}, 
		h[N+1],d[N+1],a[N+1],c[N+1],b[N+1]={2,2,2,2,2,2,2,2,2,2},s[N+1] ,t[N+1],l[N+1],M[N+1],f,f1; 
	int i; 
	for(i=1;i<=N;i++) /*使用对任意分化的三弯矩插值法*/ 
		h[i-1]=x[i]-x[i-1]; 
	d[0]=6/h[0]*((y[1]-y[0])/h[0]-1); 
	d[N]=6/h[N-1]*(0.1-(y[N]-y[N-1])/h[N-1]); 
	for(i=1;i<=N-1;i++) 
	{
     d[i]=6/(h[i-1]+h[i])*((y[i+1]-y[i])/h[i]-(y[i]-y[i-1])/h[i-1]) ; 
	 a[i]=h[i-1]/(h[i-1]+h[i]); 
	 c[i]=1-a[i]; 
	}
	c[0]=1;
	a[N]=1;
	s[0]=b[0];
	t[0]=d[0]; /*用追赶法求解三对角方程组*/ 
	for(i=0;i<=N-1;i++) 
	{ 
		l[i+1]=a[i+1]/s[i]; 
		s[i+1]=b[i+1]-l[i+1]*c[i]; 
		t[i+1]=d[i+1]-l[i+1]*t[i]; 
	}
	M[N]=t[N]/s[N]; 
	for(i=N-1;i>=0;i--) 
		M[i]=(t[i]-c[i]*M[i+1])/s[i];
	/*求算 4.563 这点的函数值*/
	f=M[3]*pow((x[4]-4.563),3)/6/h[3] 
		+M[4]*pow((4.563-x[3]),3)/6/h[3] 
		+(y[3]-M[3]*h[3]*h[3]/6)*(x[4]-4.563)/h[3] +(y[4]-M[4]*h[3]*h[3]/6)*(4.563-x[3])/h[3];
	/*求算 4.563 这点的一阶导 数值*/
	f1=-3*M[3]*pow((x[4]-4.563),2)/6/h[3] 
		+3*M[4]*pow((4.563-x[3]),2)/6/h[3] -(y[3]-M[3]*h[3]*h[3]/6)/h[3] 
		+(y[4]-M[4]*h[3]*h[3]/6)/h[3]; 
	printf("f(4.563)=%lf f'(4.563)=%lf\n",f,f1);