Best way to take a function (lambda) as input, and also have the variable within the computation itself?

Question

I have a function that should compute an integral, taking in some function as input. I'd like the code to compute a definite integral of: <some function, in terms of x. e.g., 3*x or 3*x*(1-x), etc.> * np.sin(np.pi * x)). I'm using scipy for this:

import scipy.integrate as integrate
def calculate(a):
    test = integrate.quad(a*np.sin(np.pi * x), 0, 1)
    return test

a = lambda x: 3*x
calculate(a)

Now this implementation will fail because of the discrepancy between a and x. I tried defining x as x = lambda x: x, but that won't work because I get an error of multiplying a float by a function.

Any suggestions?

Answer 1

Since you are trying to combine two symbolic expressions before computing the definite integral numerically, I think this might be a good application for sympy's symbolic manipulation tools.

from sympy import symbols, Integral, sin, pi

def calculate(a_exp):
    test = Integral(a_exp * sin(pi * x), (x, 0, 1)).evalf()
    return test

x = symbols('x')
a_exp = 3*x 
print(calculate(a_exp))

# 0.954929658551372

Note: I changed the name of a to a_exp to make it clear that this is an expression rather than a function.

If you decide to use sympy then note that you might also be able to compute the expression for the integral symbolically as well.

Update: Importing Sympy might be overkill for this

If computation speed is more important than precision, you can easily calculate the integral approximately using some simple discretized method.

For example, the functions below calculate the integral approximately with increasingly sophisticated methods. The accuracy of the first two will improve as n is increased and also depends on the nature of a_func etc.

import numpy as np
from scipy.integrate import trapz, quad

def calculate2(a_func, n=100):
    dx = 1/n
    x = np.linspace(0, 1-dx, n)
    y = a_func(x) * np.sin(np.pi*x)
    return np.sum(y) * dx

def calculate3(a_func, n=100):
    x = np.linspace(0, 1, n+1)
    y = a_func(x) * np.sin(np.pi*x)
    return trapz(y, x)

def calculate4(a_func):
    f = lambda x: a_func(x) * np.sin(np.pi*x)
    return quad(f, 0, 1)

a_func = lambda x: 3*x
print(calculate2(a_func))
# 0.9548511174430737

print(calculate3(a_func))
# 0.9548511174430737

print(calculate4(a_func)[0])
# 0.954929658551372

I'm not an expert on numerical integration so there may be better ways to do this than these.