数学建模社区-数学中国

标题: 程序有点问题？？？？？？？帮帮忙 [打印本页]

---

作者: 潇湘夜雨123456    时间: 2012-12-15 11:22
标题: 程序有点问题？？？？？？？帮帮忙
function [f,f0] = ThrSample1(x,y,y_1,y_n,x0)
%求已知数据点的第一类三次样条差值多项式及其插值点处的值
%已知数据点的x坐标向量：x
%已知数据点的y坐标向量：y
) A8 [; h9 \# v7 u7 y%左端点的一阶导数：y_1
%右端点的一阶导数：y_n
%插值点的x坐标：x0
%求得的三次样条差值多项式：f
/ U9 |% o* M' k, Q%求得的x0处的插值：f0
# f4 t6 B6 ?; U, `8 n+ V' Esyms t
f = 0.0000;
8 d1 L7 r( u* X1 Z/ ?, w! o. Wf0 = 0.0000;
8 r. Z' r1 }. `# }5 ~if(length(x)==length(y))
    n = length(x);
else
1 x2 ~$ O, k* i* K    disp('x和y的维数不相等！');
1 Y$ t! k7 a& E5 p% G9 M) `    return;
1 J% @5 Z" o7 V4 ?end              %维数检查
1 u9 M, Z: W+ z3 C9 O3 Y7 Z# e+ J9 Cfor i=1:n
    if(x(i)<=x0)&&(x(i+1)>=x0)
        ihdex = i;
        break;
    end
end              %找到x0所在区间
' D2 J% ?: A+ X% c% f4 fA = diag(2*ones(1,n));   %求解m的系数矩阵
u = zeros(n-2,1);
8 l! L6 m- Z* l8 v: `lamda = zeros(n-1,1);
* U1 e1 T+ Z1 D' G  Y! Jc = zeros(n,1);
+ ?" D& |! |* p/ |2 t4 Rfor i = 2:n-1
    u(i-1) = (x(i)-x(x-1))/(x(i+1)-x(i-1));
    lamda(i) = (x(i+1)-x(i))/(x(i+1)-x(i-1));
# G( m. C. c' v    c(i) = 3*lamda(i)*(y(i)-y(i-1))/(x(i)-x(i-1))+3*u(i-1)*(y(i+1)-y(i))/(x(i+1)-x(i));
/ w: m( g9 s$ z5 y1 `# ~ A(i,i+1) = u(i-1);
A(i,i-1) = lamda(i);   %构造系数矩阵及向量c
  S7 R7 Z$ N0 ]end
  m& }$ i& M4 f& g; e' i+ \c(1) = 2*y_1;
c(n) = 2*y_n;
9 `9 n- z0 Z3 `0 H% J+ V- I+ Dm = followup(A,c);   %用追赶法求解方程组
' z0 G) R: A% d+ ph = x(index+1) - x(index);  %x0所在区间长度
9 k1 O; J: E) w2 Qf = y(index)*(2*(t-x(index))+h)*(t-x(index+1))^2/h/h/h+...
    y(index+1)*(2*(x(index+1)-t)+h)*(t-x(index))^2/h/h/h+...
    m(index)*(t-x(index))*(x(index+1)-t)^2/h/h-...
    m(index+1)*(x(index+1)-t)*(t-x(index))^2/h/h; %x0所在区间的插值函数
1 V' m1 e! y" t: p1 Wf0 = abs(f,'t',x0)  %x0处的插值
& m/ X5 ]6 J5 [2 C# ?- ]   
   

---

欢迎光临 数学建模社区-数学中国 (http://www.madio.net/)

Powered by Discuz! X2.5