Re: Factor 2 error in Inverse Laplace Transform
- To: mathgroup at smc.vnet.net
- Subject: [mg51381] Re: Factor 2 error in Inverse Laplace Transform
- From: p-valko at tamu.edu (Peter Valko)
- Date: Fri, 15 Oct 2004 02:47:50 -0400 (EDT)
- References: <ck87ta$9mf$1@smc.vnet.net> <ckajqc$m6v$1@smc.vnet.net> <ckfsgg$ius$1@smc.vnet.net> <cklmnu$f3d$1@smc.vnet.net>
- Sender: owner-wri-mathgroup at wolfram.com
Steven/Maxim,
I agree.
Peter
steve at smc.vnet.net (Steven M. Christensen) wrote in message news:<cklmnu$f3d$1 at smc.vnet.net>...
> I do not question your suggestion to replace the interval of
> integration [0,Infinity] with [-eps,Infinity]. Curiously, Mathematica
> does exactly that for the direct transform; if we pass it the
> (indeterminate?) expression DiracDelta[t]*Log[t], we can see the
> evaluation path:
>
> In[1]:=
> LaplaceTransform[DiracDelta[t]*Log[t], t, p]
>
> Out[1]=
> Limit[Integrate[(DiracDelta[t]*Log[t])/E^(p*t),
> {t, System`LaplaceTransformDump`low$38, Infinity},
> Assumptions -> p > 0 && System`LaplaceTransformDump`low$38 < 0,
> GenerateConditions -> False, PrincipalValue -> False],
> System`LaplaceTransformDump`low$38 -> 0, Direction -> 1]
>
> The low$38 parameter is your -eps. Taking the limit is a more accurate
> way than making the substitution eps->0. Consider the next integral:
>
> In[2]:=
> Assuming[p>0 && eps>0,
> Integrate[Log[t]*E^(-p*t), {t, -eps, Infinity}]]
>
> Out[2]=
> (Gamma[0, (-eps)*p] + E^(eps*p)*(I*Pi + Log[eps]))/p
>
> It would be incorrect to substitute eps->0, since this expression is
> undefined at 0. However, Limit[%2,eps->0] gives the correct value for
> LaplaceTransform[Log[t],t,p].
>
> The substitution seems to be safe for DiracDelta[t], because then the
> integral just samples the value at 0 and eps simply will not appear in
> the output. Therefore, you have either to take the limit or to
> separate the terms containing DiracDelta from the rest of the
> integrand before making the substitution eps->0, because the integrand
> may contain both generalized and ordinary functions, e.g.
> DiracDelta[t]+Log[t].
>
> Maxim Rytin
> m.r at inbox.ru
>
> p-valko at tamu.edu (Peter Valko) wrote in message news:<ckfsgg$ius$1 at smc.vnet.net>...
> > Maxim,
> >
> > As far as your Log example is considered, it really shows
> > inconsistency.
> > The other examples are also interesting, but you make very general
> > comments and they bother me a bit. Somehow I feel they are unfair.
> >
> > For instance, in your second example, first you consider the
> > input-output pair
> > In[1]=
> > Integrate[DiracDelta[x - a], {x, 0, 1}]
> > Out[1]=
> > UnitStep[1 - a]*UnitStep[a]
> > then you are substituting a->1 into the result and Mathematica gives 0.
> >
> > You compare this to the pair
> > In[2]=
> > Integrate[DiracDelta[x], {x, 0, 1}]
> > Out[2]=
> > 1/2
> >
> > and draw the conclusion that something is wrong.
> > I respectfully disaggree with you.
> >
> > In the Help of Integrate, Mathematica explicitely warns us:
> > "When an integrand depends on a parameter, the indefinite integral
> > should be considered valid for "generic" values of the parameter. For
> > certain values the reported integral may be meaningless,..."
> >
> > In general, it is extremely easy to ask Mathematica to "derive" a parametric
> > result, and then to substitute something very special into the
> > parameter to obtain contradiction.
> >
> > Personally, I am not really interested in that.
> > My examples try to mimic the situation when a user wants to solve a
> > specific problem, arrives at a Laplace transform and tries to find the
> > original using the "InverseLaplaceTransform" service of Mathematica.
> > In such a situation a factor of 2 error is dangerous and I think with
> > a little modification Mathematica could avoid it.
> >
> > For instance here is Problem 3:
> > InverseLaplaceTransform[s/(1+s)^(3/2),s,t]
> > It gives an answer that is wrong by a factor of two, while
> > Problem 3a:
> > InverseLaplaceTransform[Apart[s/(1+s)^(3/2)],s,t]
> > gives the right answer because then Mathematica finds other rules to apply and
> > they are correct.
> >
> > I suggest to put a little "safeguard" into InverseLaplaceTransform and
> > force it to avoid this trap. This can be done without redefining
> > DiracDelta and UnitStep. (Those are used extensively by Mathematica and
> > therefore any change would impact the whole system in too many ways.)
> >
> > Regards
> > Peter
> >
> >
> >
> > ab_def at prontomail.com (Maxim) wrote in message news:<ckajqc$m6v$1 at smc.vnet.net>...
> > > p-valko at tamu.edu (Peter Valko) wrote in message news:<ck87ta$9mf$1 at smc.vnet.net>...
> > > > Hi,
> > > >
> > > > InverseLaplaceTransform is an extremely useful part of Mathematica
> > > > (since v 4.1).
> > > > However, in the following simple problem it gives the wrong answer:
> > > > Problem 1:
> > > > InverseLaplaceTransform[s/(s+1),s,t]
> > > > -1/(2*E^t)+DiracDelta[t]
> > > > where the factor 2 is completely wrong.
> > > >
> > > > To see that I slightly rewrite Problem 1 into
> > > > Problem 1a:
> > > > InverseLaplaceTransform[Apart[s/(s+1)],s,t]
> > > > and then I get the correct answer:
> > > > -E^(-t)+DiracDelta[t]
> > > >
> > > > Of course one can "Unprotect" InverseLaplaceTransform and teach it to
> > > > give the correct answer but that is not the point.
> > > > (Also one can start a long debate about the meaning of DiracDelta in
> > > > Mathematica, but that is also not the point here. )
> > > >
> > > > There are several similar simple examples when the wrong factor of two
> > > > shows up, for instance
> > > > Problem 2:
> > > > InverseLaplaceTransform[ s ArcTan[1/s],s,t]
> > > >
> > > > Using the Trace one can find out that all these "factor 2" errors have
> > > > a common origin.
> > > > Solving Problem 1 Mathematica calculates the convolution integral
> > > >
> > > > Integrate[E^(-t+x)*Derivative[1][DiracDelta][x],{x,0,t}]
> > > >
> > > > and because the lower limit is exactly zero,the factor 2 shows up in
> > > > -1/(2*E^t), that is Mathematica "halves" the Dirac delta and all its
> > > > derivatives at the origin.
> > > >
> > > > I think the InverseLaplaceTransform function could be much improved if
> > > > the above convolution integral would be evaluated more carefully.
> > > >
> > > > For instance, doing it in two steps:
> > > > res1=Integrate[E^(-t+x)*Derivative[1][DiracDelta][x],{x,-eps,t},
> > > > Assumptions -> eps>0];
> > > > res2=res1/.eps -> 0
> > > > would give the right result.
> > > > (This caution is necessary only, if generalized functions are involved
> > > > in the integration.)
> > > >
> > > > I wonder if further examples/suggestions are welcome in this group
> > > > regarding InverseLaplaceTransform???
> > > >
> > > > Peter
> > >
> > > I'd say that this is two messes mixed together. One is a rather poor
> > > implementation of integral transforms, especially when there are
> > > generalized functions involved. For example:
> > >
> > > In[1]:=
> > > LaplaceTransform[InverseLaplaceTransform[Log[p], p, t], t, p] - Log[p]
> > > LaplaceTransform[InverseLaplaceTransform[PolyGamma[p], p, t], t, p] -
> > > PolyGamma[p]
> > >
> > > Out[1]=
> > > EulerGamma
> > >
> > > Out[2]=
> > > EulerGamma
> > >
> > > We can see that LaplaceTransform/InverseLaplaceTransform aren't
> > > consistenly defined for these functions (here Mathematica doesn't
> > > internally take integrals of distributions). Another issue is the
> > > question of how the integral of DiracDelta on [0,a] should be
> > > interpreted:
> > >
> > > In[3]:=
> > > Integrate[DiracDelta[x], {x, 0, 1}]
> > > Integrate[DiracDelta[x - a], {x, 0, 1}]
> > > Integrate[DiracDelta[x - 1], {x, 1, 2}]
> > >
> > > Out[3]=
> > > 1/2
> > >
> > > Out[4]=
> > > UnitStep[1 - a]*UnitStep[a]
> > >
> > > Out[5]=
> > > 0
> > >
> > > The value of the first integral is by convention taken to be 1/2;
> > > however, substituting a=0 into Out[4] we obtain 1, and making the
> > > change of variables x->x-1 (the third integral) we get 0. Similarly
> > > for the integrals of the DiracDelta derivatives:
> > >
> > > In[6]:=
> > > Integrate[DiracDelta'[x]*phi[x], {x, -eps, Infinity}]
> > > Integrate[DiracDelta'[x]*phi[x], {x, 0, Infinity}]
> > > Integrate[DiracDelta'[x]*x, {x, 0, Infinity}]
> > >
> > > Out[6]=
> > > (-DiracDelta[eps])*phi[0] - phi'[0]
> > >
> > > Out[7]=
> > > 0
> > >
> > > Out[8]=
> > > -(1/2)
> > >
> > > The value of the first integral for eps<0 is incorrect in any case,
> > > and Out[7] and Out[8] are not consistent with each other.
> > >
> > > Maxim Rytin
> > > m.r at inbox.ru