与数组参数积分时如何修复 "only length-1 arrays can be converted to Python scalars"

Question

我正在使用 scipy.integrate 中的四边形从对象中获取有限范围内的积分。假设目标物体在打击中：

∫expm(A*X).expm(B*X)dx

其中A和B都是numpy矩阵

为了解决这个问题，我使用了打击代码：

from scipy.integrate import quad
from scipy.linalg import expm
import numpy as np

def integrand(X, A, B):
    return np.dot(expm(A*X),expm(B*X))


A = np.array([[1, 2], [3, 4]])
B = np.array([[1, 2], [3, 4]])

I= quad(integrand, 0, 1, args=(A,B))

但是对于结果我得到这个错误：

TypeError: only length-1 arrays can be converted to Python scalars

我知道当函数需要单个值但您传递的是数组时，会引发错误 "only length-1 arrays can be converted to Python scalars"。但我的问题是基于数组。那我该如何解决呢。

Answer 1

正如评论中所指出的，quad 需要一个标量函数。您始终可以通过将索引添加为输出来将函数传递给标量：

def integrand(X, A, B, ix=None):
    """ pass ix=None to return the matrix, ix = 0,1,2,3 to return an element"""
    output = np.dot(expm(A*X),expm(B*X))
    if ix is None:
        return output
    i, j = ix//2, ix%2
    return output[i,j]
I= np.array([quad(integrand, 0, 1, args=(A,B, i))[0]
for i in range(4)]).reshape(2,2)
I
>>array([[1031.61668602, 1502.47836021],
       [2253.71754031, 3285.33422634]])

请注意，这是非常低效的，因为您要计算 4 次积分，只要这不会打扰您即可。

或者，使用 trapz:

x_i = np.linspace(0,1,60)
np.trapz([integrand(x, A, B) for x in x_i], x=x_i, axis=0)
>>array([[1034.46472361, 1506.62915374],
   [2259.94373062, 3294.40845422]])

Answer 2

quadpy 进行矢量化计算。 expm 仅适用于方阵（而不适用于方阵列表）这一事实需要对矩阵形状进行一些处理。

from quadpy import quad
import numpy as np
from scipy.linalg import expm


A = np.array([[1, 2], [3, 4]])
B = np.array([[1, 2], [3, 4]])


def integrand(X):
    expAX = np.array([expm(A * x) for x in X])
    expAX = np.moveaxis(expAX, 0, -1)
    #
    expBX = np.array([expm(B * x) for x in X])
    expBX = np.moveaxis(expBX, 0, -1)
    return np.einsum("ij...,jk...->ik...", expAX, expBX)


val, err = quad(integrand, 0, 1)
print(val)

[[1031.61668602 1502.47836021]
 [2253.71754031 3285.33422633]]

与数组参数积分时如何修复 "only length-1 arrays can be converted to Python scalars"

how to fix "only length-1 arrays can be converted to Python scalars" when getting integral with an array argument

python

numpy

integrate

scipy

quad