Python：计算数组的积分

Question

import numpy

import matplotlib.pyplot as plt

from scipy import integrate

def f(x,y):
    return x*y + x**2 

def integral(x,y):
    I = integrate.quad(f, 0, x, args=(y,))[0]
    return I

def gau(x,y):
    return (1+x)*integral(x,y)


xlist = numpy.linspace(-3.0, 3.0, 100)
ylist = numpy.linspace(-3.0, 3.0, 100)
X, Y = numpy.meshgrid(xlist, ylist)
Z = gau(2, Y)

print(Z)

我不断收到错误消息 "Supplied function does not return a valid float." ，我认为问题在于我尝试将数组传递给 quad 函数。我想过用类似的东西评估数组的每个条目的积分：

yi=numpy.linspace(-3.0,3.0,100)
for i, item in enumerate(yi):
    return integral[i]=integrate.quad(f,0,x,args=(yi,))[0]

它不起作用，但它是正确的方法吗？有什么other/better建议吗？

Answer 1

如果您的 f() 在传递数组时没有提供有效的浮点数，而不是 scipy.integral 本身；

为什么要将数组传递给 f() ？

Answer 2

您可以使用通用函数（参见 https://docs.scipy.org/doc/numpy/reference/ufuncs.html) which operates on arrays element-by-element. You can create these universal functions from any function using the frompyfunc function (https://docs.scipy.org/doc/numpy/reference/generated/numpy.frompyfunc.html）：

ugau = numpy.frompyfunc(gau,2,1)
Z=ugau(X,Y)

Answer 3

你可以使用quadpy（我的一个项目）。 quadpy 在函数范围和域的维数方面完全矢量化，因此您可以插入一个 returns 向量的函数，并一次在多个区间上对该函数进行积分。您只需确保输入函数正确处理矢量化输入。在你的情况下，那将是

import numpy
import quadpy


def f(x, y):
    return numpy.multiply.outer(y, x) + numpy.multiply.outer(numpy.ones_like(y), x ** 2)


def integral(x, y):
    scheme = quadpy.line_segment.gauss_legendre(5)
    intervals = numpy.array([numpy.zeros_like(x), x])
    out = scheme.integrate(lambda t: f(t, y), intervals)
    return out


def gau(x, y):
    return (1 + x) * integral(x, y)


xlist = numpy.linspace(-3.0, 3.0, 100)
ylist = numpy.linspace(-3.0, 3.0, 100)
Z = gau(2, ylist)

print(Z)

您还可以在此处插入 xlist 而不是 2 以一次计算全部。

Python：计算数组的积分

Python: Evaluating Integral for array

python

integration

scipy

quad