Francesco Biscani added the comment:

Hi Mark,

the original code is C++, and the inf + nanj result can be reproduced by the 
following snippet:

"""
#include <complex>
#include <iostream>

int main(int argc, char * argv[]) {
  std::cout << std::complex<double>(3,0) / 0. << '\n';
  return 0;
}
"""

Here is the C99 version:

"""
#include <complex.h>
#include <math.h>
#include <stdio.h>

int main(int argc, char * argv[]) {
  double complex a = 3.0 + 0.0*I;
  printf("%f + i%f\n", creal(a/0.), cimag(a/0.));
  return 0;
}
"""

This is on Linux x86_64 with GCC 4.9. Clang gives the same result. Adding the 
"-ffast-math" compilation flag changes the result for the C version but 
apparently not for the C++ one.

The original code came from an implementation of a special function f(z) which 
has a pole near the origin. When computing f(0), the C++ code returns 
inf+nan*j, which is fine (in the sense that it is a complex infinity of some 
kind). I then need to use this result in a larger piece of code which at one 
point needs to compute 1/f(0), with the expected result being 0 mathematically. 
This works if I implement the calculation all within C/C++, but if I switch to 
Python when computing 1/f(0) I get nan + nan*j as a result.

Of course I could do all sorts of stuff to improve this specific calculation 
and avoid the problems (possibly improving the numerical behaviour near the 
pole via a Taylor expansion, etc. etc. etc.).

Beware that implementing C99 behaviour will result in a noticeable overhead for 
complex operations. In compiled C/C++ code one usually has the option to switch 
on/off the -ffast-math flag to improve performance at the expense of 
standard-conformance, but I imagine this is not feasible in Python.

----------

_______________________________________
Python tracker <rep...@bugs.python.org>
<http://bugs.python.org/issue25453>
_______________________________________
_______________________________________________
Python-bugs-list mailing list
Unsubscribe: 
https://mail.python.org/mailman/options/python-bugs-list/archive%40mail-archive.com

Reply via email to