Re: 3M multiplication for t_COMPLEX

[Bill Allombert <allomber@math.u-bordeaux.fr>]

On Wed, Apr 13, 2005 at 08:09:49PM +0200, Karim Belabas wrote:
> Hi pari-dev,
> 
>   I have been sent the following example:
> 
>   \p36
>   x = 3.32307372544194849414053517840951697E-19 - 1.12813451938986614095048572644641324E-38*I;
>   y = 3.05042759366319004212246229138501820E36 + 103557517877613079.124097118142914269*I;
> 
>   xr = real(x); xi = imag(x);
>   yr = real(y); yi = imag(y);
> 
> ? imag(x * y)
> %1 = -3.388131790 E-21
> 
> ? xr*yi + xi*yr
> %2 = 1.563040145 E-37
> 
> The problem comes from the 3M formula :                                         
>   Im(y * x) := (xr + xi) * (yr + yi) - (xr * yr + xi * yi)            
>   Re(y * x) := (xr * yr) - (xi * yi)  \\ no new multiplication here       
>                                                                                 
> Which has the side effect of combining the imaginary and real accuracy.         
That might be completly irrelevant, but a standard trick here is to pick
the 3M formula that does not cause cancellation amongst the 4 possible
formulas below:

Im = (xr + xi) * (yr + yi) - (xr * yr + xi * yi), Re accordingly
Im = (xr - xi) * (yi - yr) + (xr * yr + xi * yi), Re accordingly
Re = (xr + xi) * (yr - yi) + (xr * yi - xi * yr), Im accordingly
Re = (xr - xi) * (yr + yi) + (xi * yr - xr * yi), Im accordingly

In the case at hand, the 2 last formulas will give the good result.

What do mpc here ?

Cheers,
Bill.