Thursday, September 30, 2010

【Fortran90】 3次スプライン補間法

関数を何らかの方法で補間する必要があるので、その方法として3次スプライン補間法を考える。よく使われているらしいから。

でも、これがよく使われる理由とは?データがたくさんあるときでも、ラグランジュ補間のように当てはまりの精度が問題にならないこと、なのかな。なぜ3次かはわからない。

B-スプラインというのも耳にしたことあるけど、これはデータの点を通過しないという点で、個人的には、問題がある。

スプラインのサブルーチンはFortran Source CodesSPLINEがある。

3次スプライン補間法(Piecewise Cubic Spline)を使うには、次のサブルーチンを使う。

  • spline_cubic_val(n,t,y,ypp,tval,yval,ypval,yppval)

  • spline_cubic_set(n,t,y,ibcbeg,ybcbeg,ibcend,ybcend,ypp)
    • 2階微分の情報を与えるサブルーチン。
    • スプライン補間の境界条件はここで設定する。自然境界(natural boundary。2つの境界点における近似関数の2階微分が0となる)を使う場合、ibcbeg=2, ybcbeg=0.0D0, ibcend=2, ybcend=0.0D0とする必要がある。

Codingの手順(自然境界条件の場合)

1. n, t, yにそれぞれデータ・ポイントの数(整数)、argumentの配列(倍精度、n個の要素)、関数値の配列(倍精度、n個の要素)を設定する。

2. ``spline_cubic_set''サブルーチンをコールして2階微分の配列yppを得る。

3. ``spline_cubic_val''サブルーチンをコールして、argumentがtval(倍精度、input)の時の関数の近似値yval(倍精度、output)、1階微分ypval(倍精度、output)、2階微分yppval(倍精度、output)を得る。


Note 1: 微分可能でない関数の近似精度は低い。特に、微分可能でない点の近傍においてそれは顕著。

No comments:

Post a Comment