scipy.special.erf raises a runtime warning with 0.j
>>> from scipy.special import erf >>> print (erf(0.j)) __main__:1: RuntimeWarning: invalid value encountered in erf 0j
This warning is only printed once (even if I do scipy.special.errprint(0)), but I don't understand why it's printed at all. Really, 0.j is the same number as 0. and it has no problem with that one.
I suppose there are two questions: 1) Is there any way to suppress this warning? 2) Is this warning a bug, or am I missing something?
I (think I) tracked down the error function in the scipy source tree. It is located in: scipy/special/specfun/specfun.f (subroutine CERROR). This function does not raise the warning (It works fine when called from a simple fortran program).
You can turn off warnings with numpy.seterr:
0.j is not the same as 0.. The former is a complex number, that latter is just a float.
>>> type(0.j) <class 'complex'> >>> type(0.) <class 'float'>
The complex erf and real erf use different algorithms, e.g.
>>> erf(complex(1)) (0.84270079294971512+0j) >>> erf(1) 0.84270079294971478
Since the real erf and complex erf use different algorithms, some warnings in the complex erf will not be present in the real erf. If we check the Fortran implementation, we'll find:
SUBROUTINE CERROR(Z,CER) C ... Z1=Z C ... CS=Z1 CR=Z1 DO 10 K=1,120 CR=CR*Z1*Z1/(K+0.5D0) CS=CS+CR IF (CDABS(CR/CS).LT.1.0D-15) GO TO 15 10 CONTINUE
In particular, Z = 0 + 0j, so Z1 = 0 + 0j, so CS = CR = 0 + 0j before the loop. At the first iteration of the loop, we get:
- CR ← CR × Z12 / (K + 0.5) = 0 + 0j
- CS ← CS + CR = 0 + 0j
and then the condition requires CR/CS, which is a 0/0, which is an invalid floating point operation, and thus the warning.
This is a minor issue which can be "fixed" easily by checking if Z == 0 at the beginning. You could report a bug if you find this behavior erratic.