Re: Bug Report - Two numerical values for a same variable
- To: mathgroup at smc.vnet.net
- Subject: [mg54255] Re: Bug Report - Two numerical values for a same variable
- From: "Carl K. Woll" <carlw at u.washington.edu>
- Date: Mon, 14 Feb 2005 21:50:50 -0500 (EST)
- Organization: University of Washington
- References: <cuhr4q$978$1@smc.vnet.net> <cukc65$lup$1@smc.vnet.net> <200502130521.AAA03676@smc.vnet.net> <opsl4914q7iz9bcq@monster> <002301c51217$ad7438b0$b5803e44@Dell> <cup6bc$e11$1@smc.vnet.net>
- Sender: owner-wri-mathgroup at wolfram.com
Bob, Objects with head Rational are atoms, so you don't have access to the integers inside a rational object. It's the samething as complex numbers. For example, Complex[6,7] /. 7->1 6 + 7 I Carl Woll "DrBob" <drbob at bigfoot.com> wrote in message news:cup6bc$e11$1 at smc.vnet.net... > I'm not sure why Maxim's method doesn't encounter the bug. > > But this seems even more odd, to me: > > f[w_] = (Cos[w] - Cos[(23/100)*w]/E^(w^2/2))/w; > a = E^Integrate[f[w], {w, 0, Infinity}]; > a /. 20000 -> 20000.`20. > > E^((1/2)*(-EulerGamma - Log[2] - > Derivative[1, 0, 0][ > Hypergeometric1F1][0, 1/2, > -(529/20000)])) > > Cases[a, 20000, Infinity] > {} > > FullForm[a] > FullForm[E^(Rational[1, 2]* > (-EulerGamma - Log[2] - > Derivative[1, 0, 0][ > Hypergeometric1F1][0, > Rational[1, 2], Rational[ > -529, 20000]]))] > > If ReplaceAll can find -529/20000 in a, why can't it (or Cases) find > 20000? > > Bobby > > On Sun, 13 Feb 2005 22:02:14 -0000, David W. Cantrell > <DWCantrell at sigmaxi.org> wrote: > >> ----- Original Message ----- >> From: "DrBob" <drbob at bigfoot.com> To: mathgroup at smc.vnet.net >> To: "David W. Cantrell" <DWCantrell at sigmaxi.org>; >> <mathgroup at smc.vnet.net> >> Sent: Sunday, February 13, 2005 16:43 >> Subject: [mg54255] Re: Bug Report - Two numerical values for a same variable >> >> >>> >> But wasn't Ling _already_ using arbitrary-precision numbers in the >>> >> input? >>> >>> No, he was using exact numbers, and got an answer that was still exact. >>> The substitution >>> >>> -529/20000 -> -529/20000`20 >>> >>> made the expression approximate, but more than machine-precision. >>> Numbers like 20000`20 are called "arbitrary precision" for some >>> unfathomable reason. >> >> >> Hmm. Thanks for the info. For some reason (apparently, not a good one), I >> must have thought that "exact" and "arbitrary-precision" were synonymous. >> >> In any event, this thread has spurred me to upgrade to version 5.1. >> Having done so, I now see that a bug which I had diagnosed in version 5.0 >> is still present in version 5.1. >> >> In[1]:= >> Exp[Integrate[(Cos[w] - Exp[-(w^2/2)]*Cos[(23/100)*w])/w, {w, 0, >> Infinity}]] >> Out[1]= >> E^((1/2)*(-EulerGamma - Log[2] - Derivative[1, 0, >> 0][Hypergeometric1F1][0, 1/2, -(529/20000)])) >> Version 5.0's exact result was incorrect because it was >> missing -EulerGamma, BTW. >> >> In[2]:= N[%] >> Out[2]= 0.5439146227050919 >> >> which is correct. >> >> In[3]:= N[%%, 20] >> Divide :: infy : Infinite expression 1/0. encountered. ... >> Out[3]= >> 0.52983935469483822112357538077323822033`20.000000000000007 >> >> Note that this is incorrect for exactly the same reason which I diagnosed >> in my original response in this thread: Mathematica treats the derivative >> of the hypergeometric function as if it were 0. It's easy to see this >> because Out[6] is the same as Out[3]. So this is a bug present in both >> versions 5.0 and 5.1. >> >> In[6]:= >> N[E^((1/2)*(-EulerGamma - Log[2])), 20] >> Out[6]= >> 0.52983935469483822112357538077323822033`20.000000000000007 >> >> But note that, when asked about that derivative separately, Mathematica >> knows that it's not 0: >> >> In[8]:= >> N[Derivative[1, 0, 0][Hypergeometric1F1][0, 1/2, -(529/20000)], 20] >> Out[8]= >> -0.05243686946160788253943021003847243652`20. >> >> Finally, I'm rather surprised that Maxim's >> >>> >> 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 >> >> doesn't suffer from the same bug. Why does his method work correctly, >> while N[a, 20] encounters the bug? >> >> David Cantrell >> >> >> >>> On Sun, 13 Feb 2005 00:21:26 -0500 (EST), David W. Cantrell >>> <DWCantrell at sigmaxi.org> wrote: >>> >>> > 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. >> > > > > -- > DrBob at bigfoot.com > www.eclecticdreams.net >
- References:
- Re: [Mathematica 5.1] Bug Report - Two numerical values for a same variable
- From: "David W. Cantrell" <DWCantrell@sigmaxi.org>
- Re: [Mathematica 5.1] Bug Report - Two numerical values for a same variable