#2172 erf(1+5*%i) not accurate

closed
nobody
None
5
2012-08-17
2011-03-17
Raymond Toy
No

erf(1.0+5*%i) returns

- 2.7837671212125964e+9 %i - 1.0786562260474622e+9

But the true answer (as given by Wolfram/Alpha) is

- 2.7837702922089653e+9 %i - 1.0786931161985406e+9

Only about 4 digits are correct.

Discussion

  • Barton Willis

    Barton Willis - 2011-03-18

    By the way: using a 1F1 representation for erf, the values agree with Wolfram | Alpha. But the hypergeometric code has to switch to bigfloats to do this:

    (%i11) my_erf(x) := nfloat(2*x*hypergeometric([1/2],[3/2],-x^2)/sqrt(%pi),[]);
    (%o11) my_erf(x):=nfloat((2*x*hypergeometric([1/2],[3/2],-x^2))/sqrt(%pi),[])

    (%i12) load(hypergeometric)$

    (%i14) my_erf(1.0+5*%i);
    (%o14) -2.783777029220896b9*%i-1.078693116198541b9

     
  • Raymond Toy

    Raymond Toy - 2011-03-24

    I think the inaccuracy comes from gamma-incomplete, which is used to evaluate erf. For this particular argument, gamma-incomplete uses the series.

    I think the test of when to use the series could be relaxed a bit. When I set *gamma-imag* to 0.1, the continued fraction is used, which gives

    - 2.7837770292209086e+9 %i - 1.0786931161985383e+9

    This is better.

    I've done some experiments using the continued fraction http://functions.wolfram.com/06.06.10.0005.01. This seems to work a bit better.

     
  • Dieter Kaiser

    Dieter Kaiser - 2011-03-24

    Hello Ray,

    you are right. The problem is the function gamma_incomplete. The continued fraction expansion of the function tends to converge to a wrong value for some ranges of the argument.

    In the past, I had searched for some hints in the literature to solve this problem more general. But I had not found a simple solution. Therefore, I had done a heuristic approach and checked several ranges for the arguments to achieve a better approximation of the function gamma_incomplete.

    I have to work again on this problem.

    Dieter Kaiser

     
  • Raymond Toy

    Raymond Toy - 2011-03-24

    Here is an alternative method of computing gamma_incomplete:

    http://functions.wolfram.com/06.06.10.0007.01

    There appear to be no restrictions on the range of validity of this continued fraction, so it can probably be used for other areas. I've tested this a bit (with oct), and it converges quickly for things on or near the negative real axis. The few tests I've done indicate that accuracy is good.

     
  • Raymond Toy

    Raymond Toy - 2011-03-29

    Actually, the continued fraction

    http://functions.wolfram.com/06.06.10.0009.01

    appears to converge faster and to be more accurate. The previously mentioned fraction has some issues for points near the negative real axis where the accuracy is not as good as we might expect.

     
  • Raymond Toy

    Raymond Toy - 2012-08-17
    • status: open --> closed
     
  • Raymond Toy

    Raymond Toy - 2012-08-17

    Appears to have been fixed some time ago.

     

Log in to post a comment.

Get latest updates about Open Source Projects, Conferences and News.

Sign up for the SourceForge newsletter:





No, thanks