x=[-1,0,1,2],y=[0,1,0,0]のデータで
区間x=0~1 をスプライン補間で計算させています。

MuPAD でcubicSplineを用いた場合と
C言語によるアルゴリズム辞典から作ったソフトでは計算結果が
微妙に異なります。

どちらが3次スプライン補間として正しいのかお教え願えないでしょうか?
あるいはどちらも正しいとして、スプラインの種別が違うのでしょうか?

非常に漠然としていますが、よろしくお願いします。

「自分のツールだとこういう結果だった」というようなアドバイスでも大歓迎です。


X      MuPAD        C言語によるアルゴ

0      1            1
0.125  0.922851563  0.9488281
0.25   0.8203125    0.853125
0.375  0.698242188  0.7246094
0.5    0.5625       0.575
0.625  0.418945313  0.4160156
0.75   0.2734375    0.259375
0.875  0.131835938  0.1167969

A 回答 (2件)

x=[-1,0,1,2],y=[0,1,0,0]


この3次スプライン補間を、次のように定義し計算しました。
 (1)3次式で、節点での関数値はすべて一致
 (2)2回微分可能で、導関数、2次導関数が連続
 (3)両端 x=-1, x=2 で2次導関数が0

「'」は導関数、「"」は二次導関数です。

f1(-1)=0, f1(0)=1; f2(0)=1, f2(1)=0; f3(1)=0, f3(2)=0
f1'(0)=f2'(0); f2'(1)=f3'(1)
f1"(-1)=0; f1"(0)=f2"(0); f2"(1)=f3"(1); f3"(2)=0

これを解くと、条件を満たす3次式f1,f2,f3は一組だけ存在し、

f1(x) = -0.6 x^3 - 1.8 x^2 - 0.2 x + 1
f2(x) = x^3 - 1.8 x^2 - 0.2 x + 1
f3(x) = -0.4 x^3 + 2.4 x^2 - 4.4 x + 2.4

区間[-1,0],[0,1],[1,2] にそれぞれf1,f2,f3を当てはめてつなぐと、2回連続微分可能な3次スプライン補間が得られます。

f2(x)の[0,1]での値は次のようになります。ご質問の「C言語によるアルゴ」に一致しているようです。ご所望のスプライン補間が、上記の定義と異なる場合は補足してください。

f2(0/8)=1
f2(1/8)=0.948828125
f2(2/8)=0.853125
f2(3/8)=0.724609375
f2(4/8)=0.575
f2(5/8)=0.416015625
f2(6/8)=0.259375
f2(7/8)=0.116796875
f2(8/8)=0

この回答への補足

補足します。

MuPADのcubicSplineの計算結果は、この4点を通る3次多項式と一致しているようです。

そのようなこともあるのかな?と思いましたが、3次スプライン補間の定義を考えると、なんとなく釈然としません。

補足日時:2005/04/11 09:51
    • good
    • 0
この回答へのお礼

大変明確な回答をありがとうございます。

私もただ今、他のアルゴリズム参考書で3次スプラインを検算させようとしています。
ただ、相当時間が掛かりそうですので、まずはお礼をのべておきます。
他にも情報があればお願いいたします。

お礼日時:2005/04/11 09:42

スプラインというのは,とびとびの点に対して,


滑らかにつながる曲線を補完する方法ですよね.

私の記憶では,
「指定された点をすべて通るもの」

「指定された点を必ずしも通らない」
いろいろなスプラインがあります.

どちらのものか忘れましたが,「Bスプライン」というものがあります.


3次スプラインといえども,この2つはあると思います.

その点はよろしいでしょうか?

いろいろなスプラインがあるので,
市販のライブラリとの違いを論じる場合に,

あなたが,アルゴリズムのコード化に失敗しているのか,もともとねらっている式が違うのかをまずは,
判断する必要があります.

どのような補完をお望みですか?

それは,3次式でよいのですか?

二次式ではだめなのですか?

お望みの補完方法を明らかにした方が良い回答が得られると思います.
    • good
    • 0
この回答へのお礼

早速の回答ありがとうございます。

今回は各点を通る事が大前提ですのでB-スプラインではないです。

次数は3次、始点・終点の境界条件は「自然条件」と呼ばれるもので2次微分係数は0です。

Cのコードは上記にあるように、間違ってはなさそうです。

お礼日時:2005/04/11 09:46

お探しのQ&Aが見つからない時は、教えて!gooで質問しましょう!

このQ&Aを見た人が検索しているワード

このQ&Aと関連する良く見られている質問

Q3次スプライン補間?

3次スプライン補間ですが、3次関数で補間するので既知の3点より、
補間したい値を求めると思っていたのですが、
下記の資料を見ると、既知の2点から値を求めています。
接線を使って、計算しているみたいですが、イマイチ分かりません。
分かりやすく教えていただけないでしょうか?
http://www.caero.mech.tohoku.ac.jp/publicData/Daiguji/Chapter4.pdf

Aベストアンサー

> 3次関数で補間するので

 ここを誤解なさってるようです。そうはなくて、「2階微分まで連続な区分的3次関数」で補間するんです。

「n-1階微分まで連続な区分的n次関数」f(x)というのは、あらかじめ点列 k[j] (j=0,1,…, N) (これらを"knot" (「つなぎ目」という意味)と言います)が決めてあるものとし、また、f[j](x)をn次関数の列 (j=1,2,…,N)として、 k[0] ≦ x ≦ k[N]で定義される関数 f(x) が
(1) k[j-1] ≦ x ≦ k[j] のとき f(x) = f[j](x)
であって、しかも
(2) fのxによる0階微分f(x), 1階微分f'(x), 2階微分f''(x), …, n-1階微分(f^(n-1))(x)がどれもk[0] ≦ x ≦ k[N]の範囲で至る所連続
である、という意味。

 もちろん(1)から、k[j-1] < x < k[j] ならf(x), f'(x), f''(x), …, (f^(n-1))(x)が連続なのは自明。だから、xが丁度knot k[1], …, k[N-1]に等しいときにf(x)が(2)を満たすことが、f(x)が(2)を満たすための必要十分条件。つまり、
  (f[j]^(r))(k[j+1]) = (f[j+1]^(r))(k[j+1]) (r=0,1,2,…, n-1; j=1,…,N-1) …(★)
ということです。
 n次関数ならn+1個の係数がある。(ご質問には「3次関数で補間するので既知の3点」とありますが、勝手な3次関数なら既知の4点が必要。これはラグランジュ補間ですね。)
 しかし(2)の制約条件が付くため、f[j](x)は勝手なn次関数では駄目です。関数がN個あるんだから全部で(n+1)N個の係数がある。けれども、k[1], …, k[N-1]のN-1個のknotそれぞれについてn個の拘束条件(★)があるために、自由度は(n+1)N - n(N-1) = N+n だけしかない。この自由度を使って補間するんだから、結局「既知のN+n点 <x[i], y[i]> (i=1,2,…,N+n)を通るようにする」という補間をする。すなわち、
  f(x[i]) = y[i] (i=1,2,…,N+n)
と上記の式★とを連立した連立方程式を解くことによってfが決まるんです。

 なお問題によっては、一番端のknot k[0], k[N]におけるfのr階微分が予め指定されていたり(たとえば (f^(r))(k[0])=(f^(r))(k[N])=0 (r=2,3,…,n)とか)、あるいはfが周期的( (f^(r))(k[0])=(f^(r))(k[N]) (r=0,1,…,n) )であったり、といった条件が付く場合もあり、このときにはさらに自由度が減ることになります。

> 3次関数で補間するので

 ここを誤解なさってるようです。そうはなくて、「2階微分まで連続な区分的3次関数」で補間するんです。

「n-1階微分まで連続な区分的n次関数」f(x)というのは、あらかじめ点列 k[j] (j=0,1,…, N) (これらを"knot" (「つなぎ目」という意味)と言います)が決めてあるものとし、また、f[j](x)をn次関数の列 (j=1,2,…,N)として、 k[0] ≦ x ≦ k[N]で定義される関数 f(x) が
(1) k[j-1] ≦ x ≦ k[j] のとき f(x) = f[j](x)
であって、しかも
(2) fのxによる0階微分f(x), 1階微分f'(x), 2...続きを読む

Qエクセルにて、3次スプラインのマクロを組みました。

エクセルにて、3次スプラインのマクロを組みました。
結果がところどころしか表示されません。
何度見直しても自分では間違いがわかりません。
どなたか、間違い箇所がわかりませんでしょうか。


A列には5行目から、0-720までの721個の少数値
0
1.005586592
2.011173184

B列には5行目から、721個の少数値
5.83
5.69
5.66

D列には5行目から、0-720まで、1ずつ増加の整数値
0
1
2

E列には補間値をマクロから代入
以下マクロです。
お手数をおかけしますが、目を通して指摘していただけませんでしょうか。

Sub 補間_3次スプライン()

Dim i As Long
Dim h(1000) As Double
Dim dif1(1000) As Double
Dim dif2(1000) As Double
Dim data_count As Long
Dim lngDataEnd As Long
Dim dataX(1000) As Double
Dim dataY(1000) As Double
Dim x, y, yy0, yy1, yy2, yy3 As Double
lngDataEnd = ThisWorkbook.Sheets("3次スプライン").Range("A65536").End(xlUp).Row
data_count = lngDataEnd - 5

For i = 1 To data_count
dataX(i) = ThisWorkbook.Sheets("3次スプライン").Cells(i + 4, 1)
dataY(i) = ThisWorkbook.Sheets("3次スプライン").Cells(i + 4, 2)
Next i

h(0) = 0
dif2(0) = 0

On Error Resume Next
For i = 1 To data_count
h(i) = dataX(i) - dataX(i - 1) '//間隔を計算
dif1(i) = (dataY(i) - dataY(i - 1)) / h(i) '//一次微分を計算
Next i

For i = 1 To data_count
'二次微分を計算
dif2(i) = (dif1(i + 1) - dif1(i)) / (dataX(i + 1) - dataX(i - 1))
Next i

i = 1
For x = 0 To 720
If x < dataX(i) Then
yy0 = dif2(i - 1) / (6 * h(i)) * (dataX(i) - x) * (dataX(i) - x) * (dataX(i) - x) '第1項
yy1 = dif2(i) / (6 * h(i)) * (x - dataX(i - 1)) * (x - dataX(i - 1)) * (x - dataX(i - 1)) '第2項
yy2 = (dataY(i - 1) / h(i) - h(i) * dif2(i - 1) / 6) * (dataX(i) - x) '第3項
yy3 = (dataY(i) / h(i) - h(i) * dif2(i) / 6) * (x - dataX(i - 1)) '第4項
y = yy0 + yy1 + yy2 + yy3
ThisWorkbook.Sheets("3次スプライン").Cells(x + 4, 5) = y
Else: i = i + 1
End If
Next x
End Sub

エクセルにて、3次スプラインのマクロを組みました。
結果がところどころしか表示されません。
何度見直しても自分では間違いがわかりません。
どなたか、間違い箇所がわかりませんでしょうか。


A列には5行目から、0-720までの721個の少数値
0
1.005586592
2.011173184

B列には5行目から、721個の少数値
5.83
5.69
5.66

D列には5行目から、0-720まで、1ずつ増加の整数値
0
1
2

E列には補間値をマクロから代入
以下マクロです。
お手数をおかけしますが、目を通して指摘していただけま...続きを読む

Aベストアンサー

3次スプラインの計算式は確認していないので、それは置いておくとして、
結果がところどころしか表示されない原因は、最後のほうの
If x < dataX(i) Then
の部分ではないでしょうか。
この方法では、「If x < dataX(i) Then」を満たさなければ、そのxの行は表示されません。

次のようにしてはどうですか。

i = 1
For x = 0 To 720
Do Until x < dataX(i)
i = i + 1
If i > data_count Then Exit Do
Loop
If i > data_count Then Exit For
yy1 = ・・・・・
・・・・・・・
ThisWorkbook.Sheets("3次スプライン").Cells(x + 4, 5) = y
Next x

Qスプライン補間関数が実装されている数学ライブラリについて

スプライン補間関数が実装されているC++のライブラリがあれば、教えていただけないでしょうか?
Boostにありそうだったので見てみたのですが、探し方が悪かったせいか見つけることができなかったので、質問させていただきました。
ご存じの方がいらっしゃったら、ご教授いただければと存じます。
よろしくお願いいたします。

Aベストアンサー

以下はC++によるスプライン補間の実装例です。

http://people.sc.fsu.edu/~burkardt/cpp_src/spline/spline.html


人気Q&Aランキング

おすすめ情報