Mathematica 9 is now available
Services & Resources / Wolfram Forums
-----
 /
MathGroup Archive
2004
*January
*February
*March
*April
*May
*June
*July
*August
*September
*October
*November
*December
*Archive Index
*Ask about this page
*Print this page
*Give us feedback
*Sign up for the Wolfram Insider

MathGroup Archive 2004

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

Search the Archive

Re: Accuracy problem in Mathematica

  • To: mathgroup at smc.vnet.net
  • Subject: [mg49108] Re: [mg49061] Accuracy problem in Mathematica
  • From: Garry Helzer <gah at math.umd.edu>
  • Date: Thu, 1 Jul 2004 05:26:06 -0400 (EDT)
  • References: <200406300934.FAA05317@smc.vnet.net>
  • Sender: owner-wri-mathgroup at wolfram.com

On Jun 30, 2004, at 5:34 AM, aaaa wrote:

> Hello,
>
> I'm having a problem with a calculation in Mathematica which I can't
> solve. I have an expression which I know to be (from analytical  
> reasons)
> always between 0 and 1. It's a function of a and n ( n being natural  
> and
> a rational) and it looks like this:
>
> 1/(1-a^2)^n +
>
> Sum[((2*n - k - 1)!/((n - 1)!*(n - k)!*2^(2*n - k)))*(1/(1 + a)^k -  
> 1/(1
> - a)^k), {k, 1, n}]
>
>
>
> Let's say a=0.5.
>
> Now, when I try to calculate for small n, it's ok. When calculating for
> large n's (around 400 and above) I'm starting to get wrong results (the
> number not being between 0 and 1). The problem is that the first term
> (the first line before the sum) is VERY VERY close to the negative of
> the second term (the sum), and it's getting closer as n grows. When
> using large n's, Mathematica says they are the same number or even that
> the last term is bigger (which means the whole expression becomes
> negative) - which is wrong. It's a matter of accuracy, and I'm not sure
> how I can fix it.
>
> Can anybody help me?
>
>
>
> Itamar
>
>

Well, you can increase accuracy. In fact you can calculate it exactly  
and then convert to floating point.

g[a_, n_] := 1/(1 - a^2)^n +

Sum[((2*n - k - 1)!/((n - 1)!*(n - k)!*2^(2*n - k)))*(1/(1 + a)^k - 1/(1
- a)^k), {k, 1, n}]

In[25]:=
g[5/10,400]

Out[25]=
516307148039893425547168143532290205071415791918265630465198976168080842 
750613\
809499416727419043626362456623184319361397387028360858469908614633295744 
633328\
862859160709473304622880050348452004683108502066519849464065629946092128 
637050\
746846813671905860596256092371406885602718086010232583371639746378672299 
278469\
876235477885047244577634657496497715524728140641132606408239769827006443 
281375\
57151280261324604136115868279422184239/\
918815807062947172002982616630447315694409668429175232919346394540127635 
379951\
214984232652958021474172881030243332247949592807000680057999868079932529 
077985\
374292871633644370258227968982931003670645733650408779245993508719468999 
377001\
957500228588384000620222104017057282938170614727767933541288546863279277 
349439\
745536223805363981702198726988983034665773318277629283750950553689943505 
174547\
935790431818538430777692511101882728448

In[26]:=
%//N

Out[26]=
0.0561927

Interestingly, Mathematica 5.0 gives a closed form for this function

In[28]:=
f[a_,n_]=1/(1-a^2)^n +

Sum[((2*n - k - 1)!/((n - 1)!*(n - k)!*2^(2*n - k)))*(1/(1 + a)^k - 1/(1
- a)^k), {k, 1, n}]


Out[28]=
\!\(\((1 - a\^2)\)\^\(-n\) - \((2\^\(1 - 2\ n\)\ ¹\ Csc[2\
         n\ ¹]\ \((HypergeometricPFQRegularized[{1, 1 - n}, {
             2 - 2\ n}, \(-\(2\/\(\(-1\) + a\)\)\)] + a\ \
HypergeometricPFQRegularized[{1, 1 - n}, {2 - 2\
                   n}, \(-\(2\/\(\(-1\) +
                      a\)\)\)] - HypergeometricPFQRegularized[{1, 1 -  
n}, {2 -
                         2\ n}, 2\/\(1 +
               a\)] + a\ HypergeometricPFQRegularized[{1, 1 - n}, {2 -  
2\ n}, \
2\/\(1 + a\)])\))\)/\((\((\(-1\) +
                       a)\)\ \((1 + a)\)\ \(\((\(-1\) + n)\)!\)\  
Gamma[n])\)\)

but it has problems:

In[27]:=
f[5/10,400]

 From In[27]:=
\!\(\*
   RowBox[{\(°::"indet"\), \(\(:\)\(\ \)\), "\<\"Indeterminate
       expression \\!\\(0\\\\ ¹\\\\
     ComplexInfinity\\) encountered. \\!\\(\\*ButtonBox[\\\"MoreÉ\\\", \
ButtonStyle->\\\"RefGuideLinkText\\\", ButtonFrame->None, \
ButtonData:>\\\"General::indet\\\"]\\)\"\>"}]\)

Out[27]=
Indeterminate


Garry Helzer
gah at math.umd.edu


  • Prev by Date: Re: MATHEMATICA AND FINANCIAL CALCULATIONS?
  • Next by Date: Joining 2D arrays
  • Previous by thread: Re: Accuracy problem in Mathematica
  • Next by thread: Re: Accuracy problem in Mathematica