Scipy 对多个列表进行数值积分

Question

我有一张图表，其中包含跨度步长和相应的旋转值。我需要对每个步骤执行数值积分以获得斜率值。我想知道因为 scipy 集成中已经有内置函数，例如梯形规则或辛普森规则。如何在没有任何附加功能的情况下在两个数组或数据列表上实现？

import scipy
fraction_of_span = [0, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1]
rotation = [0.33, 1.34, 2.62, 3.41, 3.87, 4.02, 3.87, 3.41, 2.62, 1.34, 0] 
result = scipy.trapz(fraction_of_span, rotation, 10)

预期结果：

result = [x0, x1, .........xn]

Answer 1

如上提议

import scipy
fraction_of_span = [0, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1]
rotation = [0.33, 1.34, 2.62, 3.41, 3.87, 4.02, 3.87, 3.41, 2.62, 1.34, 0]
result = scipy.trapz(fraction_of_span, rotation, 10)
print(result)
-2.6665

上面已经构建的解决方案示例，使用 simp

from scipy.integrate import simps
y = rotation
x = fraction_of_span


result_simps = simps(y, x)
print(result_simps)
2.6790000000000003

请注意，结果非常相似，因为方法略有不同。请注意符号应为正，因为积分仅在正值之间（旋转元素均为正）

有很好的资料符合你的要求，或许你想看看这里：docs.scipy.org/doc/scipy/reference/tutorial/integrate.html

让我们尝试获得类似的矢量结果。

为此，您可以转到上述函数并修改它们以获得结果。所以我去 https://github.com/numpy/numpy/blob/master/numpy/lib/function_base.py#L4081-L4169 和 modify/create 一个新的函数如下：

def trapz_modified(y, x=None, dx=1.0, axis=-1):
    """
    Integrate along the given axis using the composite trapezoidal rule.
    Integrate `y` (`x`) along given axis.
    Parameters
    ----------
    y : array_like
        Input array to integrate.
    x : array_like, optional
        The sample points corresponding to the `y` values. If `x` is None,
        the sample points are assumed to be evenly spaced `dx` apart. The
        default is None.
    dx : scalar, optional
        The spacing between sample points when `x` is None. The default is 1.
    axis : int, optional
        The axis along which to integrate.
    Returns
    -------
    trapz : float
        Definite integral as approximated by trapezoidal rule.
    See Also
    --------
    sum, cumsum
    Notes
    -----
    Image [2]_ illustrates trapezoidal rule -- y-axis locations of points
    will be taken from `y` array, by default x-axis distances between
    points will be 1.0, alternatively they can be provided with `x` array
    or with `dx` scalar.  Return value will be equal to combined area under
    the red lines.
    References
    ----------
    .. [1] Wikipedia page: https://en.wikipedia.org/wiki/Trapezoidal_rule
    .. [2] Illustration image:
           https://en.wikipedia.org/wiki/File:Composite_trapezoidal_rule_illustration.png
    Examples
    --------
    >>> np.trapz([1,2,3])
    4.0
    >>> np.trapz([1,2,3], x=[4,6,8])
    8.0
    >>> np.trapz([1,2,3], dx=2)
    8.0
    >>> a = np.arange(6).reshape(2, 3)
    >>> a
    array([[0, 1, 2],
           [3, 4, 5]])
    >>> np.trapz(a, axis=0)
    array([1.5, 2.5, 3.5])
    >>> np.trapz(a, axis=1)
    array([2.,  8.])
    """
    y = asanyarray(y)
    if x is None:
        d = dx
    else:
        x = asanyarray(x)
        if x.ndim == 1:
            d = diff(x)
            # reshape to correct shape
            shape = [1]*y.ndim
            shape[axis] = d.shape[0]
            d = d.reshape(shape)
        else:
            d = diff(x, axis=axis)
    nd = y.ndim
    slice1 = [slice(None)]*nd
    slice2 = [slice(None)]*nd
    slice1[axis] = slice(1, None)
    slice2[axis] = slice(None, -1)
    try:
        # MODIFIED HERE
        #ret = (d * (y[tuple(slice1)] + y[tuple(slice2)]) / 2.0).sum(axis)
        ret = d * (y[tuple(slice1)] + y[tuple(slice2)]) / 2.0
    except ValueError:
        # Operations didn't work, cast to ndarray
        d = np.asarray(d)
        y = np.asarray(y)
        # MODIFIED HERE
        #ret = add.reduce(d * (y[tuple(slice1)]+y[tuple(slice2)])/2.0, axis)
        ret = d * (y[tuple(slice1)]+y[tuple(slice2)])/2.0

    return ret

我们还需要以下库，位于 file/script 的顶部：

from numpy import diff
from numpy import asanyarray

让我们看看输出：

>>>trapz_modified(y, x=x)
array([0.0835, 0.198 , 0.3015, 0.364 , 0.3945, 0.3945, 0.364 , 0.3015,
       0.198 , 0.067 ])

Scipy 对多个列表进行数值积分

Scipy with numerical integration on multiple lists

python

python-3.x

scipy

numerical-integration