Learn how easy it is to sync an existing GitHub or Google Code repo to a SourceForge project! See Demo

Close

#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.