三次样条插值的Matlab实现(自然边界和第一边界条件)(精)

(第一边界条件源代码:

function y=yt1(x0,y0,f_0,f_n,x _____________(1 %第一类边界条件下三次样条插值;

%xi 所求点;

%yi所求点函数值;

%x 已知插值点;

%y 已知插值点函数值;

%f_0左端点一次导数值;

%f_n右端点一次导数值;

n = length(x0;

z = length(y0;

h = zeros(n-1,1;

k=zeros(n-2,1;

l=zeros(n-2,1;

S=2*eye(n;

fori=1:n-1

h(i= x0(i+1-x0(i;

end

fori=1:n-2

k(i= h(i+1/(h(i+1+h(i;

l(i= 1-k(i;

end

%对于第一种边界条件:

k = [1;k]; _______________________(2 l = [l;1]; _______________________(3 %构建系数矩阵 S :

fori = 1:n-1

S(i,i+1 = k(i;

S(i+1,i = l(i;

end

%建立均差表:

F=zeros(n-1,2;

fori = 1:n-1

F(i,1 = (y0(i+1-y0(i/(x0(i+1-x0(i;

end

D = zeros(n-2,1;

fori = 1:n-2

F(i,2 = (F(i+1,1-F(i,1/(x0(i+2-x0(i;

D(i,1 = 6 * F(i,2;

end

%构建函数 D :

d0 = 6*(F(1,2-f_0/h(1; ___________(4

dn = 6*(f_n-F(n-1,2/h(n-1; ___________(5

D = [d0;D;dn]; ______________(6

m= S\D;

%寻找 x 所在位置,并求出对应插值:

fori = 1:length(x

for j = 1:n-1

if (x(i<=x0(j+1&(x(i>=x0(j

y(i =( m(j*(x0(j+1-x(i^3/(6*h(j+...

(m(j+1*(x(i-x0(j^3/(6*h(j+...

(y0(j-(m(j*h(j^2/6*(x0(j+1-x(i/h(j+... (y0(j+1-(m(j+1*h(j^2/6*(x(i-x0(j/h(j ; break; else continue;

end

end

end

(2 (自然边界条件源代码:

仅仅需要对上面部分标注的位置做如下修改 :

__(1:function y=yt2(x0,y0,x __(2:k=[0;k]

__(3:l=[l;0]

__(4+(5:删除

— (6:D=[0:D:0]