python numpy.convolve to solve convolution integral with limits from 0 to t instead -t to t

Question

I have a convolution integral of the type:

To solve this integral numerically, I would like to use numpy.convolve(). Now, as you can see in the online help, the convolution is formally done from -infinity to +infinity meaning that the arrays are moved along each other completely for evaluation - which is not what I need. I obviously need to be sure to pick the correct part of the convolution - can you confirm that this is the right way to do it or alternatively tell me how to do it right and (maybe even more important) why?

res = np.convolve(J_t, dF, mode="full")[:len(dF)]

J_t is an analytical function and I can evaluate as many points as I need, dF are derivatives of measurement data. for this attempt I choose len(J_t) = len(dF) because from my understanding I do not need more.

Thank you for your thoughts, as always, I appreciate your help!

---

Background information (for those who might be interested)

These type of integrals can be used to evaluate viscoelastic behaviour of bodies (or the response of an electric circuit during change of voltage, if you feel more familiar on this topic). For viscoelasticity, J(t) is the creep compliance function and F(t) can be the deviatoric strains over time, then this integral would yield the deviatoric stresses. If you now e.g. have a J(t) of the form:

J_t = lambda p, t: p[0] + p[1]*N.exp(-t/p[2])

with p = [J_elastic, J_viscous, tau] this would be the "famous" standard linear solid. The integral limits are the start of the measurement t_0 = 0 and the moment of interest, t.

Answer 1

To get it right, I have chosen the following two functions:

a(t) = t
b(t) = t**2

It is easy to do the math and find that their "convolution" as defined in your case, takes on the values:

c(t) = t**4 / 12

So lets try them out:

>>> delta = 0.001
>>> t = np.arange(1000) * delta
>>> a = t
>>> b = t**2
>>> c = np.convolve(a, b) * delta
>>> d = t**4 / 12
>>> plt.plot(np.arange(len(c)) * delta, c)
[<matplotlib.lines.Line2D object at 0x00000000025C37B8>]
>>> plt.plot(t[::50], d[::50], 'o')
[<matplotlib.lines.Line2D object at 0x000000000637AB38>]
>>> plt.show()

So by doing the above, if both your a and b have n elements, you get the right convolution values in the first n elements of c.

Not sure if the following explanation will make any sense, but here it goes... If you think of convolution as mirroring one of the functions along the y-axis, then sliding it along the x axis and computing the integral of the product at each point, it is easy to see how, since outside of the area of definition numpy takes them as if padded with zeros, you are effectively setting an integration interval from 0 to t, since the first function is zero below zero, and the second is zero above t, since it originally was zero below zero, but has been mirrored and moved t to the right.