数值分析-课程设计论文—样条求导.docVIP

课程设计报告 课程设计题目：样条求导 姓 名 ： 专 业 ： 班 级 ： 指导教师 ： 时 间 ： 题目： 找一个已知函数，分别用样条求导方法和简化样条求导方法求其数值导数，和精确值比较，画出导数的图象，并指出它们的异同。原函数： 一、摘要 在matlab环境下熟悉用熟悉的计算机语言编程上机完成在某一点u处的导数，必须将u作为节点，其缺点显而易见。若，根据三次样条插值，有 ， --------(1)[1] ，其实不必也有，，即先做的三次样条插值多项式，再求在任意一点的一阶甚至二阶导数，，，作为在相应点的一阶或二阶数值导数在区间作的三次样条函数，,有 ，其中，假设节点是等距的，则 ，--------(2) [1] 由于 因此基本方程组（三对角方程组）为 ， --------(3) [1] 加上第一类边界条件 将代入（2）式，然后对求导，得 ，--------(4) [1] 上式再对求导，得 ，--------(5) [1] 从而 ， --------(6) [1] ， --------(7) [1] （4）、（5）、（6）、（7）统称为样条求导公式 2、简化样条求导法 简化（j=0,1, n）的求解，避免解三对角方程组，首先用三点或五点公式求出节点处的导数以替代（j=0,1, n）直接代入（4）、（5）式中便可求得一阶和二阶导数。 ①分段三点求导公式 ，，对于任取的相邻三点由三点公式，有 [1] ②分段五点求导公式【1】 由五点公式， 三、程序代码及运算结果 3次样条m固支边界法 h=0.2; x=-1:h:1; y=x.^3+x.^2; [m,n]=size(x); m0=1;mn=5; H=[]; landa=[]; mu=[]; for i=1:n-1 H(i)=x(i+1)-x(i); end; for i=1:n-2 landa(i)=H(i+1)/(H(i)+H(i+1)); mu(i)=1-landa(i); end; c=[]; for i=1:n-2 c(i+1)=3*(landa(i)*(y(i+1)-y(i))/H(i)+mu(i)*(y(i+2)-y(i+1))/H(i+1)); end; c(1)=m0; c(n)=mn; a=zeros(n,n); a(1,1)=1; a(n,n)=1; for i=2:n-1 a(i,i)=2; a(i,i-1)=landa(i-1); a(i,i+1)=landa(i-1); end; m=inv(a)*c; for j=1:n-1 t=x(j):h/(n-1):x(j+1); X(j,:)=t; Y1(j,:)=2*((t-x(j+1)).*(3*t-2*x(j)-x(j+1)+h)*y(j)-(t-x(j)).*(3*t-x(j)-2*x(j+1)-h).*y(j+1))/(x(j+1)-x(j))^3+((t-x(j+1)).*(3*t-2*x(j)-x(j+1))*m(j)+(t-x(j)).*(3*t-x(j)-2*x(j+1))*m(j+1))/(h)^2; Y2(j,:)=((6*t-2*x(j)-4*x(j+1))*m(j)+(6*t-4*x(j)-2*x(j+1))*m(j+1))/h^2+6*(x(j)+x(j+1)-2*t)*(y(j+1)-y(j))/h^3; end syms t; for j=1:n-1 S1(j)=vpa(simple(2*((t-x(j+1))*(3*t-2*x(j)-x(j+1)+h)*y(j)-(t-x(j))*(3*t-x(j)-2*x(j+1)-h)*y(j+1))/(x(j+1)-x(j))^3+((t-x(j+1))*(3*t-2*x(j)-x(j+1))*m(j)+(t-x(j))*(3*t-x(j)-2*x(j+1))*m(j+1))/(h)^2),5); S2(j)=vpa(simple(((6*t-2*x(j)-4*x(j+1))*m(j)+(6*t-4*x(j)-2*x(j+1))*m(j+1))/h^2+6*(x(j)+x(j+1)-2*t)*(y(j+1)-y(j))/h^3),5); endr=-1:0.001:1; hold on plot(r,3*r.^2+2