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?

UPDATE

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).

Answers


You can turn off warnings with numpy.seterr:

    numpy.seterr(invalid='ignore')

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.


Need Your Help

never ending job in mapreduce

java hadoop configuration mapreduce

I have set some MapReduce configuration in my main method as so

iOS I can't get my carrier name

ios core-telephony carrier

as explained here -> Retrieving Carrier Name from iPhone Programmatically