BUG: special: Invalid arguments to ellip_harm can crash Python.

For example,

    In [1]: from scipy.special import ellip_harm

    In [2]: n = 4

    In [3]: p = 2*n + 2  # Too big!

    In [4]: ellip_harm(0.5, 2.0, n, p, 0.2)

    python(74331,0x7fff7c7ce310) malloc: *** error for object 0x120c1d110: pointer being freed was not allocated
    *** set a breakpoint in malloc_error_break to debug
    Abort trap: 6

The problem was that when the function `lame_coefficients` detected
invalid arguments, it returned without initializing bufferp[0].  In
normal use, bufferp[0] is malloc'ed by `lame_coefficients` and must
be freed by the caller.  If `lame_coefficients` failed to initialize
bufferp[0] with NULL, the value returned in bufferp[0] would be garbage,
and that would cause a crash when the caller attempted to free the memory.

scipy/special/_ellip_harm.pxd

@@ -49,6 +49,11 @@ cdef extern from "lapack_defs.h":
 cdef inline double* lame_coefficients(double h2, double k2, int n, int p,
                                       void **bufferp, double signm,
                                       double signn) nogil:
+
+    # Ensure that the caller can safely call free(*bufferp) even if an
+    # invalid argument is found in the following validation code.
+    bufferp[0] = NULL
+
     if n < 0:
         sf_error.error("ellip_harm", sf_error.ARG, "invalid value for n")
         return NULL

scipy/special/tests/test_ellip_harm.py

@@ -271,3 +271,11 @@ def ellip_harm_known(h2, k2, n, p, s):
     points = np.array(points)
     assert_func_equal(ellip_harm, ellip_harm_known, points, rtol=1e-12)
 
+
+def test_ellip_harm_invalid_p():
+    # Regression test. This should return nan.
+    n = 4
+    # Make p > 2*n + 1.
+    p = 2*n + 2
+    result = ellip_harm(0.5, 2.0, n, p, 0.2)
+    assert np.isnan(result)