[Scipy-svn] r3043 - in trunk/Lib/special: cephes tests

[scipy-svn at scipy.org]

Author: cookedm
Date: 2007-05-25 13:31:38 -0500 (Fri, 25 May 2007)
New Revision: 3043

Modified:
   trunk/Lib/special/cephes/hyp2f1.c
   trunk/Lib/special/tests/test_basic.py
Log:
#424: hyp2f1 gives wrong value for specific values. Plug a singularity in a transformation by averaging over it

Modified: trunk/Lib/special/cephes/hyp2f1.c
===================================================================
--- trunk/Lib/special/cephes/hyp2f1.c	2007-05-25 16:12:00 UTC (rev 3042)
+++ trunk/Lib/special/cephes/hyp2f1.c	2007-05-25 18:31:38 UTC (rev 3043)
@@ -181,8 +181,10 @@
     double t1;
     t1 = fabs(b - a);
     if (fabs(t1 - round(t1)) < EPS) {
-        /* this transformation has a pole for b-a= +-integer */
-        goto hypdiv;
+        /* this transformation has a pole for b-a= +-integer,
+           so we average around it.
+         */
+        return 0.5*(hyp2f1(a, b*(1+1e-9), c, x) + hyp2f1(a, b*(1-1e-9), c, x));
     }
     p = hyp2f1(a, 1-c+a, 1-b+a, 1.0/x);
     q = hyp2f1(b, 1-c+b, 1-a+b, 1.0/x);

Modified: trunk/Lib/special/tests/test_basic.py
===================================================================
--- trunk/Lib/special/tests/test_basic.py	2007-05-25 16:12:00 UTC (rev 3042)
+++ trunk/Lib/special/tests/test_basic.py	2007-05-25 18:31:38 UTC (rev 3043)
@@ -1171,6 +1171,9 @@
                       gamma(1+5-2)/gamma(1+0.5*5-2)/gamma(0.5+0.5*5)],
                   [4, 0.5+4, 1.5-2*4, -1./3, (8./9)**(-2*4)*gamma(4./3)*
                       gamma(1.5-2*4)/gamma(3./2)/gamma(4./3-2*4)],
+                  # and some others
+                  # ticket #424
+                  [1.5, -0.5, 1.0, -10.0, 4.1300097765277476484],
                   ]
         for i, (a, b, c, x, v) in enumerate(values):
             cv = hyp2f1(a, b, c, x)