laplace.py

@nb.njit
def z(theta):  # funny bromwich contour
    return (
            .1446 + 3.0232 * theta ** 2 / (theta ** 2 - 3.0767 * np.pi ** 2)
            + 0.2339 * 1j * theta
    )


@nb.njit
def z_dash(theta, h=.001):
    return (z(theta + h) - z(theta - h)) / (2 * h)


terms = 200
k = np.arange(1, terms)

h = 2 * np.pi / terms
theta = -np.pi + (k - .5) * h
z_computed = z(theta)
z_dash_computed = z_dash(theta)

t1 = np.exp(z_computed * terms) * z_dash_computed / 1j


@nb.njit(fastmath=True, cache=True)
def inverse_laplace_transform_(F, t):
    return np.sum(t1 * F(terms / t * z_computed) / t, axis=-1)


def inverse_laplace_transform(F):
    F = nb.njit(F, fastmath=True, cache=True)
    return lambda t: inverse_laplace_transform_(F, t)