blob: eb13faeea48c2b2703363741eb93a498c128cc6b [file] [log] [blame]
------------------------------------------------------------------------------
-- --
-- GNAT COMPILER COMPONENTS --
-- --
-- S Y S T E M . V A L _ R E A L --
-- --
-- B o d y --
-- --
-- Copyright (C) 1992-2023, Free Software Foundation, Inc. --
-- --
-- GNAT is free software; you can redistribute it and/or modify it under --
-- terms of the GNU General Public License as published by the Free Soft- --
-- ware Foundation; either version 3, or (at your option) any later ver- --
-- sion. GNAT is distributed in the hope that it will be useful, but WITH- --
-- OUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY --
-- or FITNESS FOR A PARTICULAR PURPOSE. --
-- --
-- As a special exception under Section 7 of GPL version 3, you are granted --
-- additional permissions described in the GCC Runtime Library Exception, --
-- version 3.1, as published by the Free Software Foundation. --
-- --
-- You should have received a copy of the GNU General Public License and --
-- a copy of the GCC Runtime Library Exception along with this program; --
-- see the files COPYING3 and COPYING.RUNTIME respectively. If not, see --
-- <http://www.gnu.org/licenses/>. --
-- --
-- GNAT was originally developed by the GNAT team at New York University. --
-- Extensive contributions were provided by Ada Core Technologies Inc. --
-- --
------------------------------------------------------------------------------
with System.Double_Real;
with System.Float_Control;
with System.Unsigned_Types; use System.Unsigned_Types;
with System.Val_Util; use System.Val_Util;
with System.Value_R;
pragma Warnings (Off, "non-static constant in preelaborated unit");
-- Every constant is static given our instantiation model
package body System.Val_Real is
pragma Assert (Num'Machine_Mantissa <= Uns'Size);
-- We need an unsigned type large enough to represent the mantissa
Is_Large_Type : constant Boolean := Num'Machine_Mantissa >= 53;
-- True if the floating-point type is at least IEEE Double
Precision_Limit : constant Uns := 2**Num'Machine_Mantissa - 1;
-- See below for the rationale
package Impl is new Value_R (Uns, 2, Precision_Limit, Round => False);
subtype Base_T is Unsigned range 2 .. 16;
-- The following tables compute the maximum exponent of the base that can
-- fit in the given floating-point format, that is to say the element at
-- index N is the largest K such that N**K <= Num'Last.
Maxexp32 : constant array (Base_T) of Positive :=
[2 => 127, 3 => 80, 4 => 63, 5 => 55, 6 => 49,
7 => 45, 8 => 42, 9 => 40, 10 => 55, 11 => 37,
12 => 35, 13 => 34, 14 => 33, 15 => 32, 16 => 31];
-- The actual value for 10 is 38 but we also use scaling for 10
Maxexp64 : constant array (Base_T) of Positive :=
[2 => 1023, 3 => 646, 4 => 511, 5 => 441, 6 => 396,
7 => 364, 8 => 341, 9 => 323, 10 => 441, 11 => 296,
12 => 285, 13 => 276, 14 => 268, 15 => 262, 16 => 255];
-- The actual value for 10 is 308 but we also use scaling for 10
Maxexp80 : constant array (Base_T) of Positive :=
[2 => 16383, 3 => 10337, 4 => 8191, 5 => 7056, 6 => 6338,
7 => 5836, 8 => 5461, 9 => 5168, 10 => 7056, 11 => 4736,
12 => 4570, 13 => 4427, 14 => 4303, 15 => 4193, 16 => 4095];
-- The actual value for 10 is 4932 but we also use scaling for 10
package Double_Real is new System.Double_Real (Num);
use type Double_Real.Double_T;
subtype Double_T is Double_Real.Double_T;
-- The double floating-point type
function Exact_Log2 (N : Unsigned) return Positive is
(case N is
when 2 => 1,
when 4 => 2,
when 8 => 3,
when 16 => 4,
when others => raise Program_Error);
-- Return the exponent of a power of 2
function Integer_to_Real
(Str : String;
Val : Impl.Value_Array;
Base : Unsigned;
Scale : Impl.Scale_Array;
Minus : Boolean) return Num;
-- Convert the real value from integer to real representation
function Large_Powfive (Exp : Natural) return Double_T;
-- Return 5.0**Exp as a double number, where Exp > Maxpow
function Large_Powfive (Exp : Natural; S : out Natural) return Double_T;
-- Return Num'Scaling (5.0**Exp, -S) as a double number where Exp > Maxexp
---------------------
-- Integer_to_Real --
---------------------
function Integer_to_Real
(Str : String;
Val : Impl.Value_Array;
Base : Unsigned;
Scale : Impl.Scale_Array;
Minus : Boolean) return Num
is
pragma Assert (Base in 2 .. 16);
pragma Assert (Num'Machine_Radix = 2);
pragma Unsuppress (Range_Check);
Maxexp : constant Positive :=
(if Num'Size = 32 then Maxexp32 (Base)
elsif Num'Size = 64 then Maxexp64 (Base)
elsif Num'Machine_Mantissa = 64 then Maxexp80 (Base)
else raise Program_Error);
-- Maximum exponent of the base that can fit in Num
D_Val : Double_T;
R_Val : Num;
S : Integer;
begin
-- We call the floating-point processor reset routine so we can be sure
-- that the x87 FPU is properly set for conversions. This is especially
-- needed on Windows, where calls to the operating system randomly reset
-- the processor into 64-bit mode.
if Num'Machine_Mantissa = 64 then
System.Float_Control.Reset;
end if;
-- First convert the integer mantissa into a double real. The conversion
-- of each part is exact, given the precision limit we used above. Then,
-- if the contribution of the low part might be nonnull, scale the high
-- part appropriately and add the low part to the result.
if Val (2) = 0 then
D_Val := Double_Real.To_Double (Num (Val (1)));
S := Scale (1);
else
declare
V1 : constant Num := Num (Val (1));
V2 : constant Num := Num (Val (2));
DS : Positive;
begin
DS := Scale (1) - Scale (2);
case Base is
-- If the base is a power of two, we use the efficient Scaling
-- attribute up to an amount worth a double mantissa.
when 2 | 4 | 8 | 16 =>
declare
L : constant Positive := Exact_Log2 (Base);
begin
if DS <= 2 * Num'Machine_Mantissa / L then
DS := DS * L;
D_Val :=
Double_Real.Quick_Two_Sum (Num'Scaling (V1, DS), V2);
S := Scale (2);
else
D_Val := Double_Real.To_Double (V1);
S := Scale (1);
end if;
end;
-- If the base is 10, we also scale up to an amount worth a
-- double mantissa.
when 10 =>
declare
Powfive : constant array (0 .. Maxpow) of Double_T;
pragma Import (Ada, Powfive);
for Powfive'Address use Powfive_Address;
begin
if DS <= Maxpow then
D_Val := Powfive (DS) * Num'Scaling (V1, DS) + V2;
S := Scale (2);
else
D_Val := Double_Real.To_Double (V1);
S := Scale (1);
end if;
end;
-- Inaccurate implementation for other bases
when others =>
D_Val := Double_Real.To_Double (V1);
S := Scale (1);
end case;
end;
end if;
-- Compute the final value by applying the scaling, if any
if (Val (1) = 0 and then Val (2) = 0) or else S = 0 then
R_Val := Double_Real.To_Single (D_Val);
else
case Base is
-- If the base is a power of two, we use the efficient Scaling
-- attribute with an overflow check, if it is not 2, to catch
-- ludicrous exponents that would result in an infinity or zero.
when 2 | 4 | 8 | 16 =>
declare
L : constant Positive := Exact_Log2 (Base);
begin
if Integer'First / L <= S and then S <= Integer'Last / L then
S := S * L;
end if;
R_Val := Num'Scaling (Double_Real.To_Single (D_Val), S);
end;
-- If the base is 10, we use a double implementation for the sake
-- of accuracy combining powers of 5 and scaling attribute. Using
-- this combination is better than using powers of 10 only because
-- the Large_Powfive function may overflow only if the final value
-- will also either overflow or underflow, thus making it possible
-- to use a single division for the case of negative powers of 10.
when 10 =>
declare
Powfive : constant array (0 .. Maxpow) of Double_T;
pragma Import (Ada, Powfive);
for Powfive'Address use Powfive_Address;
RS : Natural;
begin
if S > 0 then
if S <= Maxpow then
D_Val := D_Val * Powfive (S);
else
D_Val := D_Val * Large_Powfive (S);
end if;
else
if S >= -Maxpow then
D_Val := D_Val / Powfive (-S);
-- For small types, typically IEEE Single, the trick
-- described above does not fully work.
elsif not Is_Large_Type and then S < -Maxexp then
D_Val := D_Val / Large_Powfive (-S, RS);
S := S - RS;
else
D_Val := D_Val / Large_Powfive (-S);
end if;
end if;
R_Val := Num'Scaling (Double_Real.To_Single (D_Val), S);
end;
-- Implementation for other bases with exponentiation
-- When the exponent is positive, we can do the computation
-- directly because, if the exponentiation overflows, then
-- the final value overflows as well. But when the exponent
-- is negative, we may need to do it in two steps to avoid
-- an artificial underflow.
when others =>
declare
B : constant Num := Num (Base);
begin
R_Val := Double_Real.To_Single (D_Val);
if S > 0 then
R_Val := R_Val * B ** S;
else
if S < -Maxexp then
R_Val := R_Val / B ** Maxexp;
S := S + Maxexp;
end if;
R_Val := R_Val / B ** (-S);
end if;
end;
end case;
end if;
-- Finally deal with initial minus sign, note that this processing is
-- done even if Uval is zero, so that -0.0 is correctly interpreted.
return (if Minus then -R_Val else R_Val);
exception
when Constraint_Error => Bad_Value (Str);
end Integer_to_Real;
-------------------
-- Large_Powfive --
-------------------
function Large_Powfive (Exp : Natural) return Double_T is
Powfive : constant array (0 .. Maxpow) of Double_T;
pragma Import (Ada, Powfive);
for Powfive'Address use Powfive_Address;
Powfive_100 : constant Double_T;
pragma Import (Ada, Powfive_100);
for Powfive_100'Address use Powfive_100_Address;
Powfive_200 : constant Double_T;
pragma Import (Ada, Powfive_200);
for Powfive_200'Address use Powfive_200_Address;
Powfive_300 : constant Double_T;
pragma Import (Ada, Powfive_300);
for Powfive_300'Address use Powfive_300_Address;
R : Double_T;
E : Natural;
begin
pragma Assert (Exp > Maxpow);
if Is_Large_Type and then Exp >= 300 then
R := Powfive_300;
E := Exp - 300;
elsif Is_Large_Type and then Exp >= 200 then
R := Powfive_200;
E := Exp - 200;
elsif Is_Large_Type and then Exp >= 100 then
R := Powfive_100;
E := Exp - 100;
else
R := Powfive (Maxpow);
E := Exp - Maxpow;
end if;
while E > Maxpow loop
R := R * Powfive (Maxpow);
E := E - Maxpow;
end loop;
R := R * Powfive (E);
return R;
end Large_Powfive;
function Large_Powfive (Exp : Natural; S : out Natural) return Double_T is
Maxexp : constant Positive :=
(if Num'Size = 32 then Maxexp32 (5)
elsif Num'Size = 64 then Maxexp64 (5)
elsif Num'Machine_Mantissa = 64 then Maxexp80 (5)
else raise Program_Error);
-- Maximum exponent of 5 that can fit in Num
Powfive : constant array (0 .. Maxpow) of Double_T;
pragma Import (Ada, Powfive);
for Powfive'Address use Powfive_Address;
R : Double_T;
E : Natural;
begin
pragma Assert (Exp > Maxexp);
pragma Warnings (Off, "-gnatw.a");
pragma Assert (not Is_Large_Type);
pragma Warnings (On, "-gnatw.a");
R := Powfive (Maxpow);
E := Exp - Maxpow;
-- If the exponent is not too large, then scale down the result so that
-- its final value does not overflow but, if it's too large, then do not
-- bother doing it since overflow is just fine. The scaling factor is -3
-- for every power of 5 above the maximum, in other words division by 8.
if Exp - Maxexp <= Maxpow then
S := 3 * (Exp - Maxexp);
R.Hi := Num'Scaling (R.Hi, -S);
R.Lo := Num'Scaling (R.Lo, -S);
else
S := 0;
end if;
while E > Maxpow loop
R := R * Powfive (Maxpow);
E := E - Maxpow;
end loop;
R := R * Powfive (E);
return R;
end Large_Powfive;
---------------
-- Scan_Real --
---------------
function Scan_Real
(Str : String;
Ptr : not null access Integer;
Max : Integer) return Num
is
Base : Unsigned;
Scale : Impl.Scale_Array;
Extra : Unsigned;
Minus : Boolean;
Val : Impl.Value_Array;
begin
Val := Impl.Scan_Raw_Real (Str, Ptr, Max, Base, Scale, Extra, Minus);
return Integer_to_Real (Str, Val, Base, Scale, Minus);
end Scan_Real;
----------------
-- Value_Real --
----------------
function Value_Real (Str : String) return Num is
Base : Unsigned;
Scale : Impl.Scale_Array;
Extra : Unsigned;
Minus : Boolean;
Val : Impl.Value_Array;
begin
Val := Impl.Value_Raw_Real (Str, Base, Scale, Extra, Minus);
return Integer_to_Real (Str, Val, Base, Scale, Minus);
end Value_Real;
end System.Val_Real;