3

I am using the method of lines with 'solve_ivp` to solve a nonlinear PDE:

@njit(fastmath=True,error_model="numpy",cache=True)
def thinFilmEq(t,h,dx,Ma,phiFun,tempFun):
    phi = phiFun(h)
    temperature = tempFun(h)
    hxx = (np.roll(h,1) - 2*h + np.roll(h,-1))/dx**2  # use np.roll as I'm implementing periodic BC
    p = phi - hxx
    px = (np.roll(p,-1) - np.roll(p,1))/(2*dx)
    Tx = (np.roll(temperature,-1) - np.roll(temperature,1))/(2*dx)

    flux = h**3*px/3 + Ma*h**2*Tx/2
    dhdt = (np.roll(flux,-1) - np.roll(flux,1))/(2*dx)

    return dhdt

I get the following error: TypingError: non-precise type pyobject [1] During: typing of argument at C:/Users/yhcha/method_of_lines/test_01_thinFilmEq.py (28) I suspect it is due to phiFun and tempFun. They are functions which I supply at the time of calling. I make the functions arguments to the dhdt function just to keep things more general. When I try to remove phiFun and tempFun and explicitly give the function form inside thinFilmEq, the error goes away.

Then, I see the following error TypingError: Use of unsupported NumPy function 'numpy.roll' or unsupported use of the function. I thought maybe np.roll is not supported although it is included in the official website. I tried to 'enlarge' the array to somehow manually apply the same thing as np.roll when dealing with the finite difference for periodic BC:

    def augment(x):
        x2 = np.empty(len(x)+2)
        x2[1:-1] = x
        x2[0] = x[-1]
        x2[-1] = x[0]
        return x2

    H = augment(x)
    hx = (H[2:]-[H:-2])/dx   # use this instead of hx=(roll(h,-1)-roll(h,1))/dx

My questions are:

  1. It seems that I can get numba to work, at the expense of making the codes less generally (cannot supply an arbitrary function like phiFun and elegant (e.g. cannot use a one-liner with np.roll). Are there ways to get around it or is it just the price I need to pay when using numba to 'compile' the code?

  2. The original version without numba is close to 10x slower than the Matlab version I coded, and the numba version is still around 3-4 times slower than Matlab. I don't really expect scipy to outperform Matlab, but are there other ways to speedup the code to bridge the gap?

Physicist
  • 2,848
  • 8
  • 33
  • 62

0 Answers0