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 )
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]])
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. )
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).
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
has a parameter vec_func
which gave me hope. I was expecting something like:
g_int = quadrature( g, x1, x2 )
to get the fully integrated matrix, but it gives the ValueError:
matrix must be 2-dimensional
EDIT: What I am trying to do can apparently be done in Matlab using quadv
and has already been discussed for SciPy
The real case has been made available here.
To run it you will need:
- numpy
- scipy
- matplotlib
- sympy
Just run: "python curved_beam_mrs.py"
.
You will see that the procedure is already slow, mainly because of the integration, indicated by the TODO
in file curved_beam.py
.
It will go much slower if you remove the comment indicated after the TODO
in file curved_beam_mrs.py
.
The matrix of functions which is integrated is showed in the print.txt
file.
Thank you!
See Question&Answers more detail:
os