计算二元样条上的线积分

Question

我可以用 scipy.interpolate.RectBivariateSpline 构造一些二维数据的样条并很容易地查询 area 积分。但是没有 easily/quickly 计算此类样条曲线上的线积分的函数。

from scipy.interpolate import RectBivariateSpline
import numpy as np

# example data
X, Y = np.meshgrid(np.arange(0, 100), np.arange(0, 100))
Z = np.sin(X / 10) + np.cos(Y / 10)

# construct the spline
rbs = RectBivariateSpline(X[0, :], Y[:, 0], Z)

# get integral
area_integral = rbs.integral(0, 10, 30, 50)

print(area_integral)

>>>40.9287528271697

我正在寻找的是 rbs.line_integral(0,10, 30,50) 之类的东西，它将给我从 [0,30] -> [10,50] 的行的积分;不是这两点组成的矩形的面积。

fitpack有这个功能吗？我可以看到，在幕后，scipy 调用了 dfitpack.dblint()，但我不确定从那里去哪里。

我还应该补充：我不想在一条线上制作一系列点，使用例如查询这些点rbs.ev(x, y)，然后对结果求和。这非常慢，可能会引入数值积分错误。

编辑

RectBivariateSpline.integral()给我的是：

\int_A f(x, y) dxdy

其中f(x,y)为样条逼近的函数，A为点p1=(x0,y0)p2=(x1,y1)之间矩形的面积。所以，dxdy是微分面积，A是面积，求出的积分是体积

我要找的是：

\int_C f(x, y) ds

其中f(x,y)是样条逼近的函数，C是p1=(x0,y0)和p2=(x1,y1)之间的直线。所以ds是微分长度，S是长度，求积分是面积。

Answer 1

二元样条曲面是分段多项式。例如，在双三次表达式的情况下，通过代入参数方程 X = X0 + t.δX、Y = Y0 + t.δY，您可以获得 t 中的六次多项式。通过跟踪线段并与网格相交，在每个网格单元中找到参数X0、Y0、δX和δY以及积分范围。然后如果 scipy 好心给你交叉的二元多项式的系数，你可以得到 t.

中的多项式的系数

如果ds的表达式是√dx²+dy²，则等于√δX²+δY².dt并且反导数是初等的。因此，“足以”枚举交叉单元格中的多项式并计算所有需要的系数。

如果表达式是√dx²+dy²+dz²，你卡住了，因为根据链式法则dz = (δX.P(t) + δY.Q(t)) dt，其中P、Q是多项式，积分没有一个封闭形式的表达式。

注意。由于多项式展开的复杂性，对于合理的步长，直接的 Simpson 积分计算可以更快、更准确并不是不可想象的。

计算二元样条上的线积分

Calculate Line Integral over Bivariate Spline

python

interpolation

scipy