使用SciPy集成返回矩阵或数组的函数 [英] using SciPy to integrate a function that returns a matrix or array

问题描述

我有一个符号数组，可以表示为:

I have a symbolic array that can be expressed as:

from sympy import lambdify, Matrix

g_sympy = Matrix([[   x,  2*x,  3*x,  4*x,  5*x,  6*x,  7*x,  8*x,   9*x,  10*x],
                  [x**2, x**3, x**4, x**5, x**6, x**7, x**8, x**9, x**10, x**11]])

g = lambdify( (x), g_sympy )

因此对于每个x，我得到一个不同的矩阵:

So that for each x I get a different matrix:

g(1.) # matrix([[  1.,   2.,   3.,   4.,   5.,   6.,   7.,   8.,   9.,  10.],
      #         [  1.,   1.,   1.,   1.,   1.,   1.,   1.,   1.,   1.,   1.]])
g(2.) # matrix([[  2.00e+00,   4.00e+00,   6.00e+00,   8.00e+00,   1.00e+01, 1.20e+01,   1.40e+01,   1.60e+01,   1.80e+01,   2.00e+01],
      #         [  4.00e+00,   8.00e+00,   1.60e+01,   3.20e+01,   6.40e+01, 1.28e+02,   2.56e+02,   5.12e+02,   1.02e+03,   2.05e+03]])

依此类推...


我需要对x上的g进行数值积分，比如说from 0. to 100.(在实际情况下，积分没有确切的解决方案)，在我目前的方法中，我必须对g中的每个元素进行lambdify并积分它单独.我正在使用quad进行像元素一样的集成:

and so on...


I need to numerically integrate g over x, say from 0. to 100. (in the real case the integral does not have an exact solution) and in my current approach I have to lambdify each element in g and integrate it individually. I am using quad to do an element-wise integration like:

ans = np.zeros( g_sympy.shape )
for (i,j), func_sympy in ndenumerate(g_sympy):
    func = lambdify( (x), func_sympy)
    ans[i,j] = quad( func, 0., 100. )

这里有两个问题: 1)lambdify使用了多次和 2)进行循环；而且我相信第一个是瓶颈，因为g_sympy矩阵最多具有10000个项(对于for循环来说这没什么大不了的.)

There are two problems here: 1) lambdify used many times and 2) for loop; and I believe the first one is the bottleneck, because the g_sympy matrix has at most 10000 terms (which is not a big deal to a for loop).

如上所示，lambdify允许评估整个矩阵，所以我想:有没有办法整合整个矩阵?"

As shown above lambdify allows the evaluation of the whole matrix, so I thought: "Is there a way to integrate the whole matrix?"

scipy.integrate.quadrature具有参数vec_func，这给了我希望.我期待着这样的事情:

scipy.integrate.quadrature has a parameter vec_func which gave me hope. I was expecting something like:

g_int = quadrature( g, x1, x2 )

获得完全积分的矩阵，但给出ValueError:矩阵必须是二维的

to get the fully integrated matrix, but it gives the ValueError: matrix must be 2-dimensional

我正在尝试做的显然可以在Matlab中完成使用quadv和已针对SciPy进行了讨论

What I am trying to do can apparently be done in Matlab using quadv and has already been discussed for SciPy

真实案例已在此处提供.

要运行它，您将需要:

• numpy
• scipy
• matplotlib
• 象征

只需运行:"python curved_beam_mrs.py".

您将看到该过程已经很慢，主要是因为集成，由文件curved_beam.py中的TODO表示.

You will see that the procedure is already slow, mainly because of the integration, indicated by the TODO in file curved_beam.py.

如果您删除文件curved_beam_mrs.py中TODO后指示的注释，它将大大降低速度.

It will go much slower if you remove the comment indicated after the TODO in file curved_beam_mrs.py.

集成的功能矩阵显示在print.txt文件中.

The matrix of functions which is integrated is showed in the print.txt file.

谢谢！

推荐答案

quad或quadrature的第一个参数必须是可调用的. quadrature的vec_func自变量指的是此可调用项的自变量是否是(可能是多维的)向量.从技术上讲，您可以vectorize quad本身:

The first argument to either quad or quadrature must be a callable. The vec_func argument of the quadrature refers to whether the argument of this callable is a (possibly multidimensional) vector. Technically, you can vectorize the quad itself:

>>> from math import sin, cos, pi
>>> from scipy.integrate import quad
>>> from numpy import vectorize
>>> a = [sin, cos]
>>> vectorize(quad)(a, 0, pi)
(array([  2.00000000e+00,   4.92255263e-17]), array([  2.22044605e-14,   2.21022394e-14]))

但是，这等效于显式循环a的元素.具体来说，如果您要这样做，它不会给您带来任何性能提升.因此，总而言之，问题是为什么要在这里实现以及到底要实现什么.

But that's just equivalent to explicit looping over the elements of a. Specifically, it'll not give you any performance gains, if that's what you're after. So, all in all, the question is why and what exactly you are trying to achieve here.

这篇关于使用SciPy集成返回矩阵或数组的函数的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！