| 查看: 845 | 回复: 5 | |||
[交流]
【求助】急求帮忙解决MATLAB 样条差值求解微分代码出现的问题，使代码正常运行
|
|
function y0=hillsplines3(x0,x,y) if size(x,1)~=1 x=x'; end if size(y,1)~=1 y=y'; end [a,t]=hillsplinex(x,y);a=a';t=t'; for i=1:length(t)-4 phi(i)=powerplus3(x0,t(i:i+4)); end y0=sum(a.*phi); function [ax,tx]=hillsplinex(t,y) h=diff(t);delta=diff(y)./h;n=length(h)+1 d1=pchipendpoint(h(1),h(2),delta(1),delta(2)); dn=pchipendpoint(h(n-1),h(n-2),delta(n-1),delta(n-2)); r=[d1;y';dn]; tx=[t(1)-[1:3]*h(1),t,t(length(t)+[1:3])*h(length(h))]; a=zeros(length(tx)-4,length(tx)-4); for i=1:3 a(1,i)=depowerplus3(t(1),tx(i:i+4)); end for i=length(tx)-6:length(tx)-4 a(length(tx)-4,i)=depowerplus3(t(length(t)),tx(i:i+4)); end for j=2:length(tx)-5 for i=j-1:j+1 a(j,i)=powerplus3(t(j-1),tx(i:i+4)); end end ax=a\r;tx=tx'; function y=powerplus3(x,t) c=t; for i=1:5 c(i)=[]; beta(i)=24/prod(t(i)-c); c=t; end powerplus=abs(x-t).^3; y=0.5*sum(beta.*powerplus); function y=depowerplus3(x,t) c=t; for i=1:5 c(i)=[]; beta(i)=24/prod(t(i)-c); c=t; end powerplus=3*sign(x-t).*(x-t).^2; y=0.5*sum(beta.*powerplus); function d=pchipendpoint(h1,h2,del1,del2) d=((2*h1+h2)*del1-h1*del2)/(h1+h2); if sign(d)~=sign(del1) d=0; elseif(sign(del1)~=sign(del2))&(abs(d)>abs(3*del1)) d=3*del1; end [ Last edited by huameitang05 on 2011-1-4 at 17:21 ] |
回复此楼» 猜你喜欢
售SCI-T0P文章，我:8O.5.5.1.O.54,科目齐全,可+急
已经有14人回复
-
售SCI一区T0P文章，我:8.O.55.1.O.54,科目齐全,可+急
已经有13人回复
-
越来越担心硕士毕业后没有竞争力了
已经有15人回复
-
售SCI一区文章，我:8O5.5.1.O5.4,科目全,可伽急
已经有6人回复
-
售SCI一区T0P文章，我:8.O.55.1.O54,科目全,可伽急
已经有3人回复
-
售SCI一区T0P文章，我:8.O.55.1.O.5.4,科目全,可+急
已经有5人回复
-
售SCI一区文章，我:8.O.55.1.O.54,科目齐全,可伽急
已经有7人回复
-
售SCI文章，我:8O5.5.1.O.54,科目齐全,可+急
已经有7人回复
-
售SCI-T0P文章，我:8O.5.5.1.O.54,科目齐全,可+急
已经有6人回复
-
售SCI一区T0P文章，我:8.O.55.1.O.54,科目齐全,可+急
已经有4人回复
2楼2011-01-04 19:06:32
3楼2011-01-04 19:31:42
4楼2011-01-04 19:47:16
5楼2011-01-04 19:49:24
6楼2011-02-28 21:36:36
50