[Cfp-interest 1476] Re: complex multiplication and division

[paul zimmermann]

sorry Damian, I indeed didn't checked that my inputs numbers were representable
in the binary64 format. Here is a new example:

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

int
main()
{
  double complex x, y, z;

  /* generated with ./mpcheck-double -seed 1473128 */
  x = -4.8798814420289904e-22 + I * -6.4569205261627488e+209;
  y = -1.0918976936190148e+147 + I * -8.2942999263695497e-85;
  printf ("x=(%.16e,%.16e)\n", creal (x), cimag (x));
  printf ("y=(%.16e,%.16e)\n", creal (y), cimag (y));
  z = x * y;
  printf ("x*y=(%.16e,%.16e)\n", creal (z), cimag (z));
  /* MPC gives       (-2.7232252786885612e+123,inf) */
  /* GCC 9.2.1 gives (-2.7232252786885195e+123,inf) */
  /* ICC 19.0.4.243: (-2.7232252786885195e+123,inf) */

  printf ("\n");

  x = -1.0689346237980500e+252 + I * 7.8846096489686452e-145;
  y = 5.1152592576527620e-318 + I * 6.7327111521288819e+78;
  printf ("x=(%.16e,%.16e)\n", creal (x), cimag (x));
  printf ("y=(%.16e,%.16e)\n", creal (y), cimag (y));
  z = x / y;
  printf ("x/y=(%.16e,%.16e)\n", creal (z), cimag (z));
  /* MPC gives       (-3.5164330596706266e-225,1.5876733750267203e+173) */
  /* GCC 9.2.1 gives (1.1710898434244477e-223,1.5876733750267203e+173) */
  /* ICC 19.0.4.243: (-3.5164330596706266e-225,1.5876733750267203e+173) */
}

The GCC results are with -O0 (GCC uses constant folding with -O1 or larger,
but strangely this does not change the x * y result).

I find it strange that there is no requirement for correct rounding for
multiplication or division. How can we then build accurate mathematical
functions if we can't rely on multiplication/division?

Best regards,
Paul