Re: Ada.Numerics, Accuracy of trigonometric functions

["Jeffrey R. Carter" <spam.jrcarter.not@spam.not.acm.org>]

On 10/07/2016 03:13 AM, Markus Schöpflin wrote:
>
> Size of SHORT_FLOAT =  32
> Maximum relative error of COS =  2.38419E-07
> Angle threshold =  12
> X = 4294967296.00000; COS(X) =  1.00000      <-- WHAT?

For GNAT (I have 4.9.3) Ada.Numerics.Generic_Elementary_Functions is defined in 
terms of Ada.Numerics.Aux, which defines

type Double is new Long_Long_Float;

Cos, after some checks for zero, converts its argument to Double and calls 
Aux.Cos. Aux.Cos, for abs X > Pi/4, reduces abs X and then invokes the FPU 
instruction fcos or fsin on the reduced value, depending on the quadrant.

Reduce is defined as

    procedure Reduce (X : in out Double; Q : out Natural);
    --  Implements reduction of X by Pi/2. Q is the quadrant of the final
    --  result in the range 0 .. 3. The absolute value of X is at most Pi.

and contains the following

       K  : Double := X * Two_Over_Pi;
    begin
       --  For X < 2.0**32, all products below are computed exactly.
       --  Due to cancellation effects all subtractions are exact as well.
       --  As no double extended floating-point number has more than 75
       --  zeros after the binary point, the result will be the correctly
       --  rounded result of X - K * (Pi / 2.0).

Since X = 2**32, this doesn't apply. It looks to me as if Reduce is converting X 
to zero, which is then passed to fcos, giving 1.

This result is independent of the actual type used to instantiate 
Generic_Elementary_Functions.

The C cosf function probably just passes its argument to fcos.

Investigation of what GNAT does for Cos with Cycle is left as an exercise for 
the reader.

-- 
Jeff Carter
"It's symbolic of his struggle against reality."
Monty Python's Life of Brian
78