GSL中的插值

最新推荐文章于 2024-12-12 09:41:40 发布

没有啥的昵称

最新推荐文章于 2024-12-12 09:41:40 发布

阅读量2.3k

点赞数

CC 4.0 BY-SA版权

分类专栏： GSL 文章标签：算法

本文链接：https://blog.youkuaiyun.com/yxmlhc/article/details/120583410

GSL 专栏收录该内容

48 篇文章

订阅专栏

插值

本章描述的是进行插值的函数。本库提供了各种插值方法，包括Cubic、Akima和Steffen样条。插补类型是可互换的，不需要重新编译就可以使用不同的方法。插值可以定义为正规边界条件和周期边界条件。附加函数可用于计算插值函数的导数和积分。例程提供了插值一维和二维数据集。

这些插值方法产生通过每个数据点的曲线。用平滑曲线插值噪声数据参见基样条。

本节中描述的函数声明在头文件gsl_interp.h和gsl_spline.h中。

30.1 一维插值介绍

给定一组数据点(x1, y1)…(xn, yn)，本节描述的例程计算一个连续插值函数y(x)，使y(xi) = yi。插补是分段光滑的，它在端点的行为取决于所使用的插补的类型。

30.2 一维插值函数

给定数据集的插值函数存储在gsl_interp对象中。是由以下函数创建的。

gsl_interp

一维差值的工作空间。

gsl_interp * gsl_interp_alloc(const gsl_interp_type * T, size_t size)

本函数返回一个指针，指向一个新分配的类型为T的，size大小的数据点的插值对象。

int gsl_interp_init(gsl_interp * interp, const double xa[], const double ya[],

size_t size)

本函数初始化数据(xa, ya)的插值对象interp，其中xa和ya是大小为size的数组。插值对象(gsl_interp)不保存数据数组xa和ya，只存储从数据计算出的静态状态。xa数据数组总是假定是严格有序的，且x值不断增加；没有定义其他约定的行为。

void gsl_interp_free(gsl_interp * interp)

本函数释放了差值对象interp。

30.3 一维插值类型

插值库提供以下插值类型：

gsl_interp_type

gsl_interp_linear

线性插值，这种插值方法不需要任何额外的内存。

gsl_interp_polynomial

多项式插值，这种方法只用于插值少量的点，因为多项式插值对较好的数据集也引入了大振荡。插值多项式的项数等于点数。

gsl_interp_cspline

具有自然边界条件的三次样条，得到的曲线在每个区间上是分段三次的，在提供的数据点上具有匹配的一阶导数和二阶导数。二阶导数在第一个点和最后一个点都为零。

gsl_interp_cspline_periodic

具有周期边界条件的三次样条，得到的曲线在每个区间上是分段三次的，在提供的数据点上具有匹配的一阶导数和二阶导数。第一个点和最后一个点的导数也是匹配的。注意，数据中的最后一个点必须与第一个点具有相同的y值，否则生成的周期性插值将在边界处具有不连续。

gsl_interp_akima

具有自然边界条件的非圆形Akima样条，本方法采用了Wodicka的非圆角算法。

gsl_interp_akima_periodic

具有周期边界条件的非圆形Akima样条，本方法采用了Wodicka的非圆角算法。

gsl_interp_steffen

Steffen方法保证了给定数据点间插值函数的单调性。因此，极小值和极大值只能准确地出现在数据点上，并且数据点之间不可能存在虚假的振荡。插值函数在每个区间上是分段三次的。得到的曲线及其一阶导数保证是连续的，但二阶导数可能是不连续的。

30.4 一维索引查找和加速

搜索的状态存储在gsl_interp_accel对象中，它是一种用于插值查找的迭代器。

gsl_interp_accel

本工作空间存储用于插值查找的状态变量。它缓存索引查找之前的值。当后续插值点落在同一区间时，可以立即返回其索引值。

size_t gsl_interp_bsearch(const double x_array[], double x, size_t index_lo,

size_t index_hi)

本函数返回数组x_array的索引i，使x_array[i] <= x < x_array[i+1]。搜索索引的范围在[index_lo, index_hi]。定义HAVE_INLINE时，使用本函数的内联版本。

gsl_interp_accel * gsl_interp_accel_alloc(void)

本函数返回一个指向加速器对象的指针，它是一种用于插值查找的迭代器。它跟踪查找的状态，从而允许应用各种加速策略。

size_t gsl_interp_accel_find(gsl_interp_accel * a, const double x_array[],

size_t size, double x)

本函数使用给定的加速器a，对大小为size的数据数组x_array执行查找操作。这就是在计算插值期间执行查找的方式。函数返回一个索引i，使x_array[i] <= x < x_array[i+1]。定义HAVE_INLINE时，使用本函数的内联版本。

int gsl_interp_accel_reset(gsl_interp_accel * acc);

本函数重新初始化加速器对象acc。当缓存的信息不再适用时——例如，当切换到一个新的数据集时，应该使用它。

void gsl_interp_accel_free(gsl_interp_accel* acc)

本函数释放加速器对象acc。

30.5 一维插值函数的求值

double gsl_interp_eval(const gsl_interp * interp, const double xa[],

const double ya[], double x, gsl_interp_accel * acc)

int gsl_interp_eval_e(const gsl_interp * interp, const double xa[],

const double ya[], double x, gsl_interp_accel * acc, double * y)

这些函数使用插值对象interp、数据数组xa、ya和加速器acc返回给定点x的插值y值。当x超出xa的范围时，返回错误码GSL_EDOM, y的值为GSL_NAN。

double gsl_interp_eval_deriv(const gsl_interp * interp, const double xa[],

const double ya[], double x, gsl_interp_accel * acc)

int gsl_interp_eval_deriv_e(const gsl_interp * interp, const double xa[],

const double ya[], double x, gsl_interp_accel * acc,

double * d)

这些函数使用插值对象interp、数据数组xa、ya和加速器acc返回插值函数对给定点x的导数d。

double gsl_interp_eval_deriv2(const gsl_interp * interp, const double xa[],

const double ya[], double x, gsl_interp_accel * acc)

int gsl_interp_eval_deriv2_e(const gsl_interp * interp, const double xa[],

const double ya[], double x, gsl_interp_accel * acc,

double * d2)

这些函数使用插值对象interp、数据数组xa、ya和加速器acc返回对给定点x的插值函数的二阶导数d2。

double gsl_interp_eval_integ(const gsl_interp * interp, const double xa[],

const double ya[], double a, double b,

gsl_interp_accel * acc)

int gsl_interp_eval_integ_e(const gsl_interp * interp, const double xa[],

const double ya[], double a, double b,

gsl_interp_accel * acc, double * result)

这些函数使用插值对象interp、数据数组xa、ya和加速器acc返回被插值函数在[a, b]范围内的数值积分result。

----------------------------------------

图30.4:二维插值示例