用三次样条插值求离散点斜率 matlab程序

shouce

用三次样条插值求斜率 三次样条插值的matlabt程序 function x=followup(a,b,c,d)n=length(d); a(1)=0; 4 F. e; B8 D% U, V/ N: d% O%“追”的过程 # L( w9 U/ D( |: C& `L(1)=b(1); * w! G) ^( K# x/ f7 Q/ I' ey(1)=d(1)/L(1); u(1)=c(1)/L(1); 0 P4 Q) V; O, T& k, ^% Q# ~8 @* Gfor i=2n-1) + D0 d! s, o3 O9 X* j# Y5 H    L(i)=b(i)-a(i)*u(i-1); 0 U3 B) ?3 a& g6 {2 ]$ i( G" y8 L    y(i)=(d(i)-y(i-1)*a(i))/L(i);     u(i)=c(i)/L(i); end L(n)=b(n)-a(n)*u(n-1); y(n)=(d(n)-y(n-1)*a(n))/L(n); ' Y9 K6 g& j! \) ^%“赶”的过程 x(n)=y(n); 0 ?( @5 p7 m5 h2 y! Lfor i=(n-1):-1:1     x(i)=y(i)-u(i)*x(i+1); & k  p% p5 Q- x+ iend function[s,y0]=spline3 (x,y,x0) %x,y为数表x0为插值点s表示插值函数y0为x0对应的插值函数值 : @, [) \; o% b9 P: {) {. T1 jsyms t n=length(x); %得出n 2 x4 L; L. Z  V# f" W$ |for i=1:n-1; ( m' m% K3 s8 N    h(i)=x(i+1)-x(i); end for i=2:n-1; % c5 H8 p# U, U' F4 }; G    lamda(i)=h(i)/(h(i-1)+h(i));     miu(i)=1-lamda(i);     g(i)=3*(lamda(i)*((y(i)-y(i-1))/h(i-1))+miu(i)*((y(i+1)-y(i))/h(i))); end ) I6 s, z- y( K4 B$ S. H" ~& Ag(1)=3*((y(2)-y(1))/h(1)); ) W" ^; J) i% z0 }, x$ Y; Ng(n)=3*((y(n)-y(n-1))/h(n-1)); %前边求出lamda miu和g从而可以确定系数矩阵 4 W, J* A7 K. u* }* T' _1 B3 p! i! smiu(1)=1; miu(4)=0; lamda(n)=1; lamda(1)=0; %根据第二边界条件补充两个lamda和miu的值 $ l8 w$ \5 w# O* F* t" C$ ^2 ~. Wfor i=1:n 1 f2 o+ \4 g, V) m    beta(i)=2; + G% Q7 R( w( A* h5 O5 I6 G3 Hend ( D( M! L4 c, l9 }' Q6 Im=followup(lamda,beta,miu,g) 5 M: ]" g7 T7 K3 F. z8 ?%解出m的值从而可确定st st为各段的插值多项式 for i=1:n-1 " Z! F" a& s7 y/ O/ D1 q    st(i)=(t-x(i+1))^2*(h(i)+2*(t-x(i)))*y(i)/(h(i)^3)...     +(t-x(i))^2*(h(i)+2*(x(i+1)-t))*y(i+1)/(h(i)^3)...     +(t-x(i))^2*(t-x(i+1))*m(i+1)/(h(i)^2)... % Q( i& [4 m2 h/ {& M, w! d& W    +(t-x(i+1))^2*(t-x(i))*m(i)/(h(i)^2); end + d# |. c! S9 I% {% W%得到插值的结果各段的t的表达式 %接下来要将插值点x0代入首先确定x0所在的插值区间 # A: s, }. T  d3 jfor i=1:n-1     if (x(i)>x0) ( z" Q2 D9 G* b# f8 i, l9 q5 R        in=i;     end end s=st(in); + o$ k0 [6 [# d* n* ^& ds=expand(s); , R& M9 `2 G( Y9 ?: b' `! T7 rs=collect(s,'t'); / V# u( L0 ~- h- [8 _% Ny0=subs(s,'t',x0) ( o* @5 M* J6 m1 n2 e%s是插值多项式y0是插值点的函数值 , \1 P% g3 d2 n ( s" E- g; v" T6 E 在matlab中输入 ' h% c9 d0 A  S! k1 l/ |# j x=[1 2 4 5]； y=[1 3 4 2]； spline3(x,y,2) 会得到各点的斜率 1 H% D0 R, C5 Z! x% h/ n' a0 I 4 H% ?2 {* C: M ) t$ y+ Q' W7 M