computational-geometry - 三次样条插值：如何计算第二个样条 S2

Question

哪个公式对 S2 有用？

到目前为止，这是我的代码：

program spline_1_2;
uses crt;
const
k=1000;
var
input_file:text;
out1_file:text;
out2_file:text;
n,i:integer;
a_int,b_int,h:real;
delta_1, delta_2:real;
x, alpha, beta, y, a_coef, b_coef, c_coef, d_coef:array[0..k] of real;
S, S1, S2:array[0..k] of real;
A, B, C, F, z:array[0..k] of real;

Function f_given(x:real):real;
 begin
  f_given:=Cos((x*x)/3)-Exp(Sin(2*x));
 end;
Function f_prime(x:real):real;
 begin
  f_prime:=-Cos((x*x)/3)*(1/3)*2*x-Exp(Sin(2*x))*Cos(2*x)*2;
 end;

begin
clrscr;

Assign(input_file,'in.txt');
Reset(input_file);

Assign(out1_file,'out1.txt');
Rewrite(out1_file);

Assign(out2_file,'out2.txt');
Rewrite(out2_file);

while not EOF(input_file) do 
 begin
  readln(input_file,a_int,b_int);
  readln(input_file,h);
  readln(input_file,n);
end;  
 for i:=0 to n do
  begin
    x[i]:=a_int+(i-0.5)*h;
    f_given(x[i]); 
    y[i]:= f_given(x[i]);
    f_prime(x[i]);

   A[i]:= h;
  C[i]:= 4*h;
  B[i]:= h;

 end;

for i:=1 to n-1 do
begin 
 a_coef[i]:=y[i];
 F[i]:=((y[i+1]-2*y[i]+y[i-1])*3)/(h*h);
 A[1]:=0;
  C[n-1]:=0;
end;
alpha[1]:=-C[1]/B[1];
beta[1]:=F[1]/B[1];

  for i:=1 to n-1 do
 begin
  z[i]:=(A[i]*alpha[i-1] + C[i]);
   alpha[i]:= -B[i]/z[i];
  beta[i]:= ( F[i]-A[i]*beta[i-1])/z[i];
end;

c_coef[n-1]:=(F[i] - A[i] * beta[n-1])/(C[i] + A[i] * alpha[n-1]);

for i:=n-2 downto 1 do
begin
c_coef[i]:=alpha[i]*c_coef[i+1] + beta[i];
c_coef[0]:=0;
c_coef[n]:=0;

 end;



for i:=0 to n-1 do
 begin 
  d_coef[i]:=(c_coef[i+1]-c_coef[i])/(3*h);
 end;
  b_coef[0]:=f_prime(x[0]);
  b_coef[n-1]:=f_prime(x[n-1]);

for i:=1 to n-2 do
 begin 
   b_coef[i]:=(h/3)*(c_coef[i+1]+2*c_coef[i]) + (y[i+1] - y[i])/h;
 end;

 //our spline
for i:=1 to n do
 begin
   S[i]:= a_coef[i] + b_coef[i]*(x[i-1]-x[i]) + c_coef[i]*sqr(x[i-1]-x[i]) +    d_coef[i]*sqr(x[i-1]-x[i])*(x[i-1]-x[i]);
   S1[i]:=b_coef[i]+2*c_coef[i]*(x[i-1]-x[i])+3*d_coef[i]*(x[i-1]-x[i])*(x[i-1]-x[i]);
 end;

for i:=1 to n-1 do
 begin 
   Writeln(out1_file, x[i],'  ',y[i]:3:6,'  ',S[i]:3:6);
   delta_1:=abs(y[i]-S[i]);
   Writeln(out2_file,x[i],'  ',f_prime(x[i]):3:6,'  ', S1[i]:3:6);
   delta_2:=abs(f_prime(x[i])-S1[i]);
  end;
   Writeln(out1_file,'"delta1" = ',delta_1);
   Writeln(out2_file,'"delta2" = ',delta_2);
   Readln;  
 Close(input_file);
 Close(out1_file);
 Close(out2_file);
end.

score 0 · Accepted Answer

很可能已经过时了，但无论如何先做一些笔记

1.据我了解，您有 f(x) 并想从中创建 SPLINE

一条三次曲线或多条加入的样条曲线？
S1,S2 的区别或含义是什么？
那么 f(x) 是什么？我看到 f_given 和 f_prime 是哪一个？

2. 在 SPLINE 计算之后，您必须计算偏差

在 SPLINE 和原始 f(x) 之间
以及它们的推导之间

3.你的输入是什么（a，h，n）？

我假设 a 是 x 间隔的开始
h 是点之间的 x 步长
n 是要采样的点数

如何解决：

1.得到样条

只需从 f(x) 获得 4 个点的控制点
并通过求解 SPLINE 方程计算 SPLINE 系数（您的 a、b、c、d）

在这之后：

t = < 0.0 , 1.0 >
x = x0 + (x1-x0)*t
y = a + b*t + c*t*t + d*t*t*t
y'=     b   + c*t   + d*t*t

其中 x0,x1 是 x 的范围
并且 t 是 SPLINE 参数
如果您有更多连接的 SPLINE 段，则分别计算每个段
不要忘记正确加入他们，所以如果你有点 p0,p1,p2,...

然后 SPLINE 调用序列将蜂

p0,p0,p0,p1 // 3x p0 to ensure that SPLINE starts at point p0
p0,p0,p1,p2
p0,p1,p2,p3 // here continue normally 
p1,p2,p3,p4
p2,p3,p4,p5
...
p5,p6,p7,p8
p6,p7,p8,p9 // at the end the last point is also 3x
p7,p8,p9,p9
p8,p9,p9,p9

2.计算差异（偏差）

只为 x = x0 .. x1
并计算差异
存储最大结果

computational-geometry - 三次样条插值：如何计算第二个样条 S2

1 回答 1

Related

Reference