MathGroup Archive 2004

[Date Index] [Thread Index] [Author Index]

Search the Archive

Re: Factor 2 error in Inverse Laplace Transform

  • To: mathgroup at smc.vnet.net
  • Subject: [mg51373] Re: Factor 2 error in Inverse Laplace Transform
  • From: p-valko at tamu.edu (Peter Valko)
  • Date: Fri, 15 Oct 2004 02:47:16 -0400 (EDT)
  • References: <ck87ta$9mf$1@smc.vnet.net> <ckajqc$m6v$1@smc.vnet.net> <200410120557.BAA19233@smc.vnet.net> <cklmaj$f00$1@smc.vnet.net>
  • Sender: owner-wri-mathgroup at wolfram.com

Bob,


1) " If this works in a single example, that doesn't mean it's a
reliable strategy in general."

I am sorry if my previous message gave you the false impression that
it is enough safeguard to insert Apart before doing
InverseLaplaceTransform.
That was just an example to show that going through different paths,
the function InverseLaplaceTransform can come up with the wrong or
right answer.

In general, the safeguard has to be inserted into the code of
InverseLaplaceTransform itself. Whenever the convolution theorem is
used to obtain the result, the algorithm needs to check whether
DiracDelta (or its derivatives) are involved and need to do things
differently, if the answer is yes. Another fix should be inserted
regarding the inverse of Log[s] and or the transform (or rather
generalized transform) of 1/t.


2) " Everything I've seen on the group in three years says the Laplace
functions are unreliable."

If this statement is true it is very unfortunate. In fact the
LaplaceTransform and its inverse are both extremely powerful functions
and represent the quintessence of what a CAS is good for. Both these
functions rely heavily on Mathematica's integration package and I
think it is needless to tell this Group that Integrate is one of the
greatest things in Mathematica.

As far as I understand LaplaceTransform causes little trouble (in fact
it is quite good.) Most of the haevy usage and most of the problems
are related to InverseLaplaceTransform.

There are various reasons why InverseLaplaceTransform might have a bad
reputation.
 
i) Oftentimes the user knows a certain result but Mathematica rather
gives back the expression unevaluated. This is embarassing, but it can
be easily fixed by teaching Mathematica more and more Laplace
transform pairs. Andre Mallet posted some fixes on MathSource for
version 4.1 and it is my understanding that some of those fixes made
it into version 5.
Recently Urs Graf published an excellent book (Applied Laplace
Transforms and z-Transforms for Scientists and Engineers, 2004, 
Birkhauser) on the subject. In his package he "teaches" Mathematica a
lot of pairs (especially those involving Bessel functions and other
special functions either in the original or in the transform or both.)
Also he uses a more elaborate concept of "pseudo functions", showing
that if you really miss something you can always create it in
Mathematica.

ii) More embarrasing is (at least for me) if Mathematica comes up with
a "pseudo answer", which still contains an unevaluated "Integrate". I
think this happens inadvertently and should be eliminated.

iii) Finally in certain cases Mathematica comes up with an answer that
is just deadly wrong. This is the most dangerous situation and it is
almost always related to the interpretation of the integral of a
DiracDelta or one of its derivatives, or to the integral of the log
function. These problems could be easily fixed. (Unfortunately, I do
not see a process in place where reasonable suggestions are accepted
and handled to make the improvements.)
 

Additional remarks:

In my opinion teaching more and more pairs to Mathematica is not the
real solution. In fact the more and more the function
InverseLaplaceTransform is relying on general rules and calling the
Integrate function, the more and more powerful it will be. This is the
strength of Mathematica: applying very simple and general rules.
However, if there is a glitch in one of the general rules, it MUST be
corrected soon, because it effects more and more results and destroys
the reputation of an extreemly valuable tool.

In MathSource there is a function available called GWR. This function
can calculate the inverse of a Laplace transform at a given value of
time with any required accuracy. I use it often to check an "analytic"
result coming out from InverseLaplaceTransform and this check always
reveals if there is a problem (somewhat similarly to making a plot if
you are not sure about a limit...)
  
Regards
Peter




DrBob <drbob at bigfoot.com> wrote in message news:<cklmaj$f00$1 at smc.vnet.net>...
> >> I suggest to put a little "safeguard" into InverseLaplaceTransform and
> >> force it to avoid this trap.
> > Problem 3a:
> > InverseLaplaceTransform[Apart[s/(1+s)^(3/2)],s,t]
> 
> If this works in a single example, that doesn't mean it's a reliable strategy in general.
> 
> Everything I've seen on the group in three years says the Laplace functions are unreliable.
> 
> You'd best check every answer.
> 
> Bobby
> 
> On Tue, 12 Oct 2004 01:57:50 -0400 (EDT), Peter Valko <p-valko at tamu.edu> wrote:
> 
> > 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
> >
> >
> >
> >


  • Prev by Date: Re: Re: Re: Re: normal distribution random number generation
  • Next by Date: Integrate vs. NIntegrate
  • Previous by thread: Re: Re: Factor 2 error in Inverse Laplace Transform
  • Next by thread: Re: Re: Factor 2 error in Inverse Laplace Transform