0

I am trying to implement the Galerkin discontinous method in Python to solve Partial Differential Equations. This method involves solving a large number of systems of ordinary differential equations with the form

enter image description here

where the matrix $M^i$ is definied by the $L^2$ product, i.e. each entry of the matrix is an integral of the product of two functions $N_k$, and $N_l$, that is

enter image description here

where $I_i=[-1,1]$.

To do these integrals, I have defined the function DotProduct(G, H) as follows

from scipy.integrate import quad

def DotProduct(G, H):
    F = lambda x:G(x)*H(x)
    w = quad(F, -1.0, 1.0)[0]
return w

The use of this function is as follows

f01 = lambda x: 0.5*(1-x)
f11 = lambda x: 0.5*(1+x)
Integral = DotProduct(f01, f11)
print(Integral)

My problem is that 90% of the execution time of the code is dedicated to calculating these integrals.

I've read about speeding up the computation time for integrals using numba, but have had no success implementing it in the function DotProduct(G, H).

Does anyone have any advice?

  • At first all functions have to be compileable by numba (like @Jérôme Richard showed in his answer) The second step is to create a low level callable. This avoids the function call overhead in the Python Interperter. Example: https://stackoverflow.com/a/61089817/4045774 – max9111 Feb 14 '22 at 08:33

1 Answers1

1

Almost all SciPy function are not yet supported by Numba. This includes quad. That being said, the function called by quad could be accelerated by Numba. The point is you cannot use lambda functions the way you do too because it is also not supported by Numba. You need to move f01 and f11 in the jitted function.

import numba as nb
from scipy.integrate import quad

@nb.njit
def F(x):
    G = lambda x: 0.5*(1-x)
    H = lambda x: 0.5*(1+x)
    return G(x) * H(x)

def DotProduct():
    w = quad(F, -1.0, 1.0)[0]
    return w

Integral = DotProduct()
print(Integral)

This is 36% faster on my machine but this is not a case where Numba is really useful... Numba can be faster at least if the F function is expensive. The initial code takes less than 6 us to complete. This is far too small.

Jérôme Richard
  • 41,678
  • 6
  • 29
  • 59