-1

So I've been trying to use the general integration method (quad) from Scipy on numpy on a specific expression; however, I get the following error

TypeError: unsupported operand type(s) for *: 'function' and 'float'

Here is the function I'm trying to integrate (seems like theres no mathjax here):

integral from 0 to a of t*y(t) * J0(u_{i} * t/a) dt where a is the length of y(t), J0 is a zeroth-order Bessel Function and u_{i} is the roots of J0(u) = 0

My code:

import numpy as np
import scipy.integrate as integrate
from scipy.special import *
y = np.load('dataset.npy')  # 1000 pts long
a = 1000
roots = jn_zeros(0, 1000)
func = lambda t: t * jv(0, t * (roots/a) )
result = integrate.quad(func * y , 0, a)
Sam
  • 145
  • 1
  • 3
  • 12

1 Answers1

0

The first problem is you are not calling your function, what you meant was

func(y)*y

otherwise you are just multiplying a function name by a number/array, which is the error you are seeing. Beyond that, the initial argument to integrate.quad is a callable, which is the function you want to integrate. It does not accept a vector, so your initial y vector is completely ignored. If anything, you want this:

result = integrate.quad(lambda y: func(y) * y , 0, a)

but this ignored the dataset completely. To integrate on samples use trapz, or simps, which accept an x array, and a y array. You probably need to generate your y array using np.linspace or something similar.

kabanus
  • 24,623
  • 6
  • 41
  • 74
  • I noticed the `arr` typo, but let's make it a pun intended. – kabanus Aug 13 '18 at 06:20
  • Thanks! but now it seems like I got the following error: ValueError: operands could not be broadcast together with shapes (5,) (1000,) – Sam Aug 13 '18 at 06:26
  • @Sam Looks like `func` is returning a matrix. You did not post enough info to debug this ([mcve]). I suggest you write a new question for your new problem, but only after you looked around at what `func` is returning and if it makes sense to you. – kabanus Aug 13 '18 at 06:29