Python Scipy integrate.quad()

Question

I'm trying to solve an integral inside a for loop where the value of the integral changes with every loop.

Ni_1D = np.array([(x/2)*(x-1), x-x**2, (x/2)*(x+1)])
Tn= [10.0, 10.0, 10.0]
Tt=[0, 0, 0]
y_coor = [0.0, 1.25, 2.5]
x_coor = [4.0, 4.0, 4.0]
dN1D = np.array([x-.5, x-2*x, x+.5])

for i in range(3):
    def f(x):
         return Ni_1D[i]*(Tn[i]*y_coor[i]+Tt[i]*x_coor[i])* dN1D[i]
    fx[i] = integrate.quad(f,-1,1)

But I keep getting the error "cannot convert expression to float." Obviously something is wrong with my first input to quad.

I can't figure out how to get the correct input for f for quad(). When I paste in the function manually everything works fine. For example for i=1 f=-x*(-12.5x**2 + 12.5x). When I define the function direction as -x*(-12.5x**2 + 12.5x) everything works fine.

Yes sorry I should have included that. x is symbolic and defined using sympy: x = sympy.symbols('x') — Lane, Apr 29 '23 at 00:51
That's what I suspected. Why do you think that will work in a numeric program? — hpaulj, Apr 29 '23 at 01:23
Well I suppose I thought it would work because when I bypass the for loop and ask it to evaluate the integral directly it works fine. For example: def f(x): return -x*(-12.5x**2 + 12.5x) fx = integrate.quad(f,-1,1) works fine. If my integrand changes with each iteration I need to keep "x" as symbolic. Is there a way to evaluate the integral with a symbolic program? — Lane, Apr 29 '23 at 01:28
With `f(x)`, `quad` calls it supplying a number as the integration variable. The `x` in function argument is not the same as `x` used a sympy expression. DO NOT TRY USE SYMPY WITH SCIPY. — hpaulj, Apr 29 '23 at 02:15
How would you recommend I evaluate the integral then? I need to use sympy to create the integrand since it changes each loop. — Lane, Apr 29 '23 at 02:25
Don't confuse sympy symbols with python variables. Or sympy expressions with python functions. You are integrating 3 different functions. Sympy can't be used to hide that. — hpaulj, Apr 29 '23 at 14:32

hpaulj · Answer 1 · 2023-04-29T15:50:35.397

When code doesn't work as expected, look at the pieces. Don't assume anything.

For specific i value:

In [30]: i=1
In [31]: Ni_1D[i]*(Tn[i]*y_coor[i]+Tt[i]*x_coor[i])* dN1D[i]
Out[31]: -x*(-12.5*x**2 + 12.5*x)

Or wrapping that in the function:

In [32]: def f(x):
    ...:     return Ni_1D[i]*(Tn[i]*y_coor[i]+Tt[i]*x_coor[i])* dN1D[i]
    ...:
In [33]: f(1)
Out[33]: -x*(-12.5*x**2 + 12.5*x)
In [34]: type(_)
Out[34]: sympy.core.mul.Mul

Regardless of what argument you supply, the result is still a sympy mul object, not some number. The x in the sympy expression is not the function argument.

To supply a number for a symbol in a sympy expression, you have to use subs. And quad is not going to do that for you.

In [35]: f(1).subs(x,1)
Out[35]: 0
In [36]: f(1).subs(x,10)
Out[36]: 11250.0000000000

If you really must use sympy, use lambdify, for example as in

In [41]: def f(z):
    ...:     fn = sympy.lambdify('x',Ni_1D[i]*(Tn[i]*y_coor[i]+Tt[i]*x_coor[i])* dN1D[i]
    ...: )
    ...:     return fn(z)
    ...:
In [42]: f(10)
Out[42]: 11250.0

Doing the sympy lambdify every time quad wants a value is going to slow down the evaluation considerably. sympy is not fast.

Python Scipy integrate.quad()

1 Answers1