Factor 2 error in Inverse Laplace Transform
- To: mathgroup at smc.vnet.net
- Subject: [mg52797] Factor 2 error in Inverse Laplace Transform
- From: ab_def at prontomail.com (Maxim)
- Date: Mon, 13 Dec 2004 04:23:11 -0500 (EST)
- Sender: owner-wri-mathgroup at wolfram.com
There was a discussion of a curious LaplaceTransform glitch a while ago in this group ( http://forums.wolfram.com/mathgroup/archive/2004/Oct/msg00248.html ). Here is a similar example in version 5.1: In[1]:= InverseLaplaceTransform[p*(PolyGamma[p + 2] - PolyGamma[p]), p, t] Out[1]= 2 - 1/(2*E^t) Since PolyGamma[p + 2] - PolyGamma[p] == 1/p + 1/(p + 1), the correct answer is 2*DiracDelta[t] - E^(-t). It is interesting to trace where the error comes from. Let f(t) = 1 + E^(-t), F(p) = L[f] = 1/p + 1/(p + 1). Mathematica evaluates the inverse transform of p*F(p) by using the identity (L^-1)[p*F(p)] = convolution((L^-1)[p], (L^-1)[F(p)]) = convolution(DiracDelta'(t), f(t)). Then it writes the convolution as In[2]:= Integrate[(1 + E^(-tau))*DiracDelta'[t - tau], {tau, low, t}, Assumptions -> low < 0] Out[2]= -1/(2*E^t) + DiracDelta[-low + t] + DiracDelta[-low + t]/E^low where the parameter low is used to 'include' the zero in the integration range, after which Mathematica takes the limit as low tends to 0 from the left (all this can be seen by doing Trace). The first error is that in fact it is needed to include the point tau=t, not tau=0; otherwise Mathematica multiplies the value of the function at tau=t by 1/2: in Mathematica Integrate[DiracDelta[t], {t, 0, Infinity}] == 1/2, and similarly in the above case we get -E^(-t)/2 instead of -E^(-t). The next error is encountered when taking the limit: In[3]:= Limit[-1/(2*E^t) + DiracDelta[-low + t] + DiracDelta[-low + t]/E^low, low -> 0, Direction -> 1] Out[3]= 2 - 1/(2*E^t) which is the answer returned by InverseLaplaceTransform; the DiracDelta terms somehow evaluated to 2. At the same time, taking the limit term by term gives In[4]:= Limit[#, low -> 0, Direction -> 1]& /@ (-1/(2*E^t) + DiracDelta[-low + t] + DiracDelta[-low + t]/E^low) Out[4]= -1/(2*E^t) + 2*DiracDelta[t, t] which is somewhat closer to the correct answer, but also incorrect: DiracDelta[t, t] is the same as DiracDelta[t]^2, which is indeterminate -- it contains a product of distributions with intersecting supports. In fact, rather than performing all those dubious transformations such as taking pointwise limits of generalized functions, we could evaluate the convolution in a straightforward way: it is simply In[5]:= Integrate[(1 + E^(-tau))*UnitStep[tau]*DiracDelta'[t - tau], {tau, -Infinity, Infinity}] Out[5]= -(t + Sqrt[t^2])/(2*E^t*Sqrt[t^2]) When evaluated correctly, with UnitStep viewed as a generalized function, this should give the required answer 2*DiracDelta[t] - E^(-t). Unfortunately, Mathematica is unable to handle this integral as well (although, curiously, there is a way to make it work -- by rewriting the integral as Integrate[(1 + E^(-tau))*DiracDelta'[t - tau], {tau, 0, Infinity}]). The second example from Peter's post still gives an incorrect result: In[6]:= InverseLaplaceTransform[p*ArcTan[1/p], p, t] Out[6]= (2*t^2 + t*Cos[t] - Sin[t])/(2*t^2) The correct value is DiracDelta[t] + (t*Cos[t] - Sin[t])/t^2. Besides, if we try to verify Out[6] by taking the direct transform, we find out that LaplaceTransform of (t*Cos[t] - Sin[t])/t^2 is incorrect too: In[7]:= LaplaceTransform[(t*Cos[t]-Sin[t])/t^2, t, p] // FullSimplify Out[7]= p*ArcCot[p] The correct value is p*ArcCot[p] - 1. The only conclusion that comes to mind is something along the lines "if it were a horse..." Maxim Rytin m.r at inbox.ru