Matlab之三次样条画图和表达式

haidaowang

这一题是得到数据点(0,3),(1,5),(2,4),(3,1)并得到它的三次样条表达式和画出三次样条后的图图形。

    以及对数据点(-1,3),(0,5),(3,1),(4,1),(5,1)并得到它的三次样条表达式和画出三次样条后的图图形。

    用函数spline可以直接得到，都是如果是要求自然三次样条呢？那就可以在数组y的左右两侧添0。如：

csa = spline(ax,[0 ay 0]);再用xxa = linspace(0,3,1000);plot(xxa,ppval(csa,xxa));进行绘图。linespace是将0~3分成1000份，然后ppval是求三次样条在不同的xxa上的值。MATLAB中的插值函数为interp1，其调用格式为： yi= interp1(x,y,xi,'method')，其中x，y为插值点，yi为在被插值点xi处的插值结果；x,y为向量， 'method'表示采用的插值方法，MATLAB提供的插值方法有几种： 'method'是最邻近插值， 'linear'线性插值； 'spline'三次样条插值； 'cubic'立方插值．缺省时表示线性插值。注意：所有的插值方法都要求x是单调的，并且xi不能够超过x的范围。

    最后要输出表达式的话这个还是有点复杂的：使用以下函数。
       pp=interp1(ax,ay,'spline','pp');
# ^- K3 G- Y% r5 d  [9 Y2 ?       breaks=pp.breaks; %breaks的第i和i+1个元素为m和n
! D$ _* h' D- x* I: ~7 o       coefs=pp.coefs; %假设coefs的第i行为a b c d，
- B: m0 }9 E* ?    接着再用一个循环得出每个表达式输出各个表达式即可。

    一下是我的代码，写的有一点粗糙，希望别见怪啊！

   % use natural cubic spline to interpolate data point
- G! Y- P' v9 Q% q; G   % a、(0,3),(1,5),(2,4),(3,1)
+ M& Z" P) a% s8 U: I9 Y7 i   % b、(-1,3),(0,5),(3,1),(4,1),(5,1)
function page_178_1_natural_spline_script
6 r  m) L$ z  z4 \0 e4 q) Cax = [0 1 2 3];
ay = [3 5 4 1];  %对数据点(0,3),(1,5),(2,4),(3,1)进行三次样条建模，并输出表达式和图像
bx = [-1 0 3 4 5];
6 }  D9 S$ h5 O0 M- \) Jby = [3 5 1 1 1];%对数据点(-1,3),(0,5),(3,1),(4,1),(5,1)进行三次样条建模，并输出表达式和图像
csa = spline(ax,[0 ay 0]);
5 d- d1 j& O& M+ f  f# Nxxa = linspace(0,3,1000);
0 h4 C  r# f; ?* Fsubplot(1,2,1);
plot(ax,ay,'o',xxa,ppval(csa,xxa),'-');
8 ^0 O! Z$ r3 n  ?1 \& f3 j4 {2 O% d- Txlabel('a x 0~3');
1 x' V# w+ j, u% r- ^! Nylabel('a y');
) `8 v# t/ `( qtitle('equation a');
' q) G, v, b& C8 Y' m. gcsb = spline(bx,[0 by 0]);
% b; Z# u0 U. E: S3 w# [xxb = linspace(-1,5,10000);
subplot(1,2,2);
: r6 C. a/ f4 Z, qplot(bx,by,'o',xxb,ppval(csb,xxb),'-');
! z7 S* y4 A( q* |1 u% c6 `$ l# }4 Q2 jxlabel('b x -1~5');
. W% h( z3 w' t+ a% N% ~" Vylabel('b y');
title('equation b');

pp=interp1(ax,ay,'spline','pp');
: _; Z4 A& N- V: Zbreaks=pp.breaks; %breaks的第i和i+1个元素为m和n
coefs=pp.coefs; %假设coefs的第i行为a b c d，
$ r8 a+ H( H% p* ?     %breaks的第i和i+1个元素为m和n那么在区间[m,n]的函数表达式就是
2 l% G8 R8 J9 W' Z. P     %a(x-m)^3+b(x-m)^2+c(x-m)+d
; m/ a5 [, b$ g& |4 s# hsyms x
disp('对于数据点(0,3),(1,5),(2,4),(3,1)的三次样条表达式为：');
, R" ?# `# B  M% w1 M. J, O# sfor i = 1:3
    y = coefs(i,1)*((x - breaks(i))^3) + coefs(i,2)*((x - breaks(i))^2) + coefs(i,3)*((x - breaks

    (i))) + coefs(i,4)
3 Z/ y% R7 w1 _. x9 H, \2 n2 kend
ppb=interp1(bx,by,'spline','pp');
breaksb=ppb.breaks; %breaks的第i和i+1个元素为m和n
coefsb=ppb.coefs; %假设coefs的第i行为a b c d，
%breaks的第i和i+1个元素为m和n那么在区间[m,n]的函数表达式就是
3 A# i, k: O! d6 X& M, M%a(x-m)^3+b(x-m)^2+c(x-m)+d
syms x
disp('对于数据点(-1,3),(0,5),(3,1),(4,1),(5,1)的三次样条表达式为：');
% r$ m& j/ ~' b3 zfor i = 1:4
   y = coefsb(i,1)*((x - breaksb(i))^3) + coefsb(i,2)*((x - breaksb(i))^2) + coefsb(i,3)*((x -

   breaksb(i))) + coefsb(i,4)
7 Z9 ^. y: C$ Y7 M6 Q' b3 Cend

# U& t9 B4 n/ `1 O, ^

BeardRen

学习