calculate_prolate_0_0_pi.m

protable = Table[{r, N[Re[SpheroidalS1[0, 0, \[Pi], r]]]}, {r, -5, 5, 0.01}]
Export["SpheroidalS1_0_0_pi.csv", protable, "CSV"]

compare_prolate.py

from matplotlib import pyplot as plt
from scipy.special import pro_rad1
import numpy as np

mathematica_x, mathematica_y = np.genfromtxt('SpheroidalS1_0_0_pi.csv', delimiter=',', unpack=True)
y = pro_rad1(0, 0, np.pi, mathematica_x)

plt.subplot(211)
plt.plot(mathematica_x, y[0], label='SciPy')
plt.plot(mathematica_x, mathematica_y, label='Mathematica')
plt.title('SciPy and Mathematica Implementations')
plt.xticks(np.arange(-5, 5 + 1))
plt.legend()

plt.subplot(212)
plt.title('Difference between SciPy and Mathematica implementations')
plt.plot(mathematica_x, (y[0] - mathematica_y) / mathematica_y, label='(SciPy - Mathematica) / Mathematica')
plt.xticks(np.arange(-5, 5 + 1))
plt.legend(loc='lower right')
plt.show()

specfun_wrappers.c.diff

diff --git a/scipy/special/specfun_wrappers.c b/scipy/special/specfun_wrappers.c
index 8772dd1..cf2b71d 100644
--- a/scipy/special/specfun_wrappers.c
+++ b/scipy/special/specfun_wrappers.c
@@ -898,7 +898,7 @@ double prolate_radial1_nocv_wrap(double m, double n, double c, double x, double
   double r2f, r2d, r1f, cv, *eg;
   int int_m, int_n;

-  if ((x <=1.0) || (m<0) || (n<m) || \
+  if ((m<0) || (n<m) || \
      (m!=floor(m)) || (n!=floor(n)) || ((n-m)>198)) {
     sf_error("prolate_radial1_nocv", SF_ERROR_DOMAIN, NULL);
     *r1d = NPY_NAN;
@@ -950,7 +950,7 @@ int prolate_radial1_wrap(double m, double n, double c, double cv, double x, doub
   double r2f, r2d;
   int int_m, int_n;

-  if ((x <= 1.0) || (m<0) || (n<m) || \
+  if ((m<0) || (n<m) || \
      (m!=floor(m)) || (n!=floor(n))) {
     sf_error("prolate_radial1", SF_ERROR_DOMAIN, NULL);
     *r1f = NPY_NAN;