Repository: mpmath/mpmath (python)

Issue:
Hankel functions incorrect for large imaginary parts and standard precision
For example:
```
>>> hankel1(1,1+30j)
mpc(real='0.0', imag='0.0')
>>> hankel2(1,1-30j)
mpc(real='0.0', imag='0.0')
```
whereas Mathematica gives
```
In[1]:= HankelH1[1, 1. + 30 I]
Out[1]= -7.25495*10^-15 - 1.17346*10^-14 I
In[2]:= HankelH2[1, 1. - 30 I]
Out[2]= -7.25495*10^-15 + 1.17346*10^-14 I
```
Correct values can be found by increasing the precision, but this should be automatic as clearly Mathematica is doing.

The problem arises because the standard definitions of the functions involve two terms with many digits in common, and there need to be sufficient digits available when they are combined.

Pull request description:
Use sum_accurately() in hankel1/2()
Closes #637

Patch:
diff --git a/mpmath/functions/bessel.py b/mpmath/functions/bessel.py
index 836bed6..7c69afa 100644
--- a/mpmath/functions/bessel.py
+++ b/mpmath/functions/bessel.py
@@ -180,11 +180,17 @@ def h(n):
 
 @defun_wrapped
 def hankel1(ctx,n,x,**kwargs):
-    return ctx.besselj(n,x,**kwargs) + ctx.j*ctx.bessely(n,x,**kwargs)
+    def terms():
+        return [ctx.besselj(n,x,**kwargs),
+                ctx.j*ctx.bessely(n,x,**kwargs)]
+    return ctx.sum_accurately(terms)
 
 @defun_wrapped
 def hankel2(ctx,n,x,**kwargs):
-    return ctx.besselj(n,x,**kwargs) - ctx.j*ctx.bessely(n,x,**kwargs)
+    def terms():
+        return [ctx.besselj(n,x,**kwargs),
+                -ctx.j*ctx.bessely(n,x,**kwargs)]
+    return ctx.sum_accurately(terms)
 
 @defun_wrapped
 def whitm(ctx,k,m,z,**kwargs):


Test patch:
diff --git a/mpmath/tests/test_functions2.py b/mpmath/tests/test_functions2.py
index 188327f..610f7ca 100644
--- a/mpmath/tests/test_functions2.py
+++ b/mpmath/tests/test_functions2.py
@@ -2433,3 +2433,7 @@ def test_issue_634():
 
 def test_issue_908():
     assert mp.besselj(-10+0j, 0+0j) == 0
+
+def test_issue_637():
+    assert hankel1(1, 1 + 30j).ae(-7.25495e-15 - 1.17346e-14j)
+    assert hankel2(1, 1 - 30j).ae(-7.25495e-15 + 1.17346e-14j)