MathGroup Archive 2005

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

Search the Archive

Re: [Mathematica 5.1] Bug Report - Two numerical values for a same variable

  • To: mathgroup at smc.vnet.net
  • Subject: [mg54195] Re: [Mathematica 5.1] Bug Report - Two numerical values for a same variable
  • From: "David W. Cantrell" <DWCantrell at sigmaxi.org>
  • Date: Sun, 13 Feb 2005 00:21:26 -0500 (EST)
  • References: <cuhr4q$978$1@smc.vnet.net> <cukc65$lup$1@smc.vnet.net>
  • Sender: owner-wri-mathgroup at wolfram.com

Maxim <ab_def at prontomail.com> wrote:
> On Fri, 11 Feb 2005 08:42:02 +0000 (UTC), Ling Li <ling at caltech.edu>
> wrote:
>
> > Hi,
> >
> > Mathematica v5.1 is much better than v5.0 on the accuracy of
> > integration. I am not talking about numerical accuracy---v5.0 sometimes
> > gives out two answers that are off by a constant for the same input.
> >
> > However, I still get two different numerical values for a same variable
> > in v5.1
> >
> > a=Exp[Integrate[(Cos[w]-Exp[-w^2/2]*Cos[23/100*w])/w, {w,0,Infinity}]]
> > N[a]
> > N[a,20]
> >
> > Thanks for looking into this problem.
> > --Ling

I apologize for not having checked the result given in my previous post
via numerical integration. Had I done so, I would have realized that the
result was incorrect. Your result of approx 0.5439146 indeed appears to be
correct.

Problems may be compounded by the fact that I'm still using version 5.0.

> The expression which is evaluated incorrectly is N[a, 20] (I get
> Divide::infy warnings and 0.52983935469483822112 as the output). This is
> not so unusual when Hypergeometric and MeijerG are involved; for example,
> see http://forums.wolfram.com/mathgroup/archive/2004/Jul/msg00185.html .
> If we start with arbitrary-precision numbers in the input, we obtain the
> correct result without any problems:

But wasn't Ling _already_ using arbitrary-precision numbers in the input?

David

> In[1]:=
> f[w_] = (Cos[w] - E^(-w^2/2)*Cos[23/100*w])/w;
> a = E^Integrate[f[w], {w, 0, Infinity}];
> a /. -529/20000 -> -529/20000`20
>
> Out[3]=
> 0.54391462270509428
>
> Here are two ways to verify the result. First, you can use the lim
> function from
> http://forums.wolfram.com/mathgroup/archive/2003/Aug/msg00233.html to
> simplify the derivative of Hypergeometric1F1:
>
> In[7]:=
> a /. Derivative[1, 0, 0][Hypergeometric1F1][a_, b_, z_] :>
>    lim[(Hypergeometric1F1[a + eps, b, z] -
>          Hypergeometric1F1[a, b, z])/eps,
>      eps -> 0]
> N[%, 20]
>
> Out[7]=
> E^((-EulerGamma + (529*HypergeometricPFQ[{1, 1}, {3/2, 2},
> -529/20000])/10000 - Log[2])/2)
>
> Out[8]=
> 0.54391462270509428214
>
> Another way is to evaluate the integral numerically (well, almost). To do
> that, we will split the integration range into [0, 1] and [1, Infinity)
> and then integrate the term with E^(-w^2/2) using Method ->
> DoubleExponential and integrate Cos[w]/w symbolically:
>
> In[9]:=
> E^(NIntegrate[f[w], {w, 0, 1}, WorkingPrecision -> 30] +
>      (If[FreeQ[#, E],
>         Integrate[#, {w, 1, Infinity}],
>         NIntegrate[#, {w, 1, Infinity},
>           Method -> DoubleExponential, WorkingPrecision -> 30]
>       ]& /@ Expand[f[w]]))
>
> Out[9]=
> 0.5439146227050942821365022967
>
> We could also do it in a purely numerical way, integrating the
> oscillatory terms on [1, Infinity) separately by using Method ->
> Oscillatory and increasing the value of MaxRecursion, but I wasn't able
> to get the result with more than 15 digits of precision then.
>
> Mathematica has trouble with high-precision evaluation of MeijerG as
> well: N[MeijerG[{{0, 1}, {}}, {{1/2}, {}}, 1/100], 20] seems unable to
> find even one significant digit with any setting of $MaxExtraPrecision,
> and gives an incorrect result if we increase $MinPrecision (MeijerG[{{0,
> 1}, {}}, {{1/2}, {}}, z] is O[Sqrt[z]] when z->0). Perhaps there are no
> general enough methods to evaluate MeijerG in such cases, but the claim
> that Mathematica "can evaluate special functions with any parameters to
> any precision" seems a little too strong anyway.
>
> Maxim Rytin
> m.r at inbox.ru


  • Prev by Date: Re: Re: Re: how to have a blind factorization of a polinomial?
  • Next by Date: Re: Fourier Transfer and a game?!?!
  • Previous by thread: Re: [Mathematica 5.1] Bug Report - Two numerical values for a same variable
  • Next by thread: Re: Re: [Mathematica 5.1] Bug Report - Two numerical values for a same variable