gcc/ada/libgnat/s-valrea.adb

   1 ------------------------------------------------------------------------------
   2 --                                                                          --
   3 --                         GNAT COMPILER COMPONENTS                         --
   4 --                                                                          --
   5 --                      S Y S T E M . V A L _ R E A L                       --
   6 --                                                                          --
   7 --                                 B o d y                                  --
   8 --                                                                          --
   9 --          Copyright (C) 1992-2024, Free Software Foundation, Inc.         --
  10 --                                                                          --
  11 -- GNAT is free software;  you can  redistribute it  and/or modify it under --
  12 -- terms of the  GNU General Public License as published  by the Free Soft- --
  13 -- ware  Foundation;  either version 3,  or (at your option) any later ver- --
  14 -- sion.  GNAT is distributed in the hope that it will be useful, but WITH- --
  15 -- OUT ANY WARRANTY;  without even the  implied warranty of MERCHANTABILITY --
  16 -- or FITNESS FOR A PARTICULAR PURPOSE.                                     --
  17 --                                                                          --
  18 -- As a special exception under Section 7 of GPL version 3, you are granted --
  19 -- additional permissions described in the GCC Runtime Library Exception,   --
  20 -- version 3.1, as published by the Free Software Foundation.               --
  21 --                                                                          --
  22 -- You should have received a copy of the GNU General Public License and    --
  23 -- a copy of the GCC Runtime Library Exception along with this program;     --
  24 -- see the files COPYING3 and COPYING.RUNTIME respectively.  If not, see    --
  25 -- <http://www.gnu.org/licenses/>.                                          --
  26 --                                                                          --
  27 -- GNAT was originally developed  by the GNAT team at  New York University. --
  28 -- Extensive contributions were provided by Ada Core Technologies Inc.      --
  29 --                                                                          --
  30 ------------------------------------------------------------------------------
  31
  32 with System.Double_Real;
  33 with System.Float_Control;
  34 with System.Unsigned_Types; use System.Unsigned_Types;
  35 with System.Val_Util;       use System.Val_Util;
  36 with System.Value_R;
  37
  38 pragma Warnings (Off, "non-static constant in preelaborated unit");
  39 --  Every constant is static given our instantiation model
  40
  41 package body System.Val_Real is
  42
  43    pragma Assert (Num'Machine_Mantissa <= Uns'Size);
  44    --  We need an unsigned type large enough to represent the mantissa
  45
  46    Is_Large_Type : constant Boolean := Num'Machine_Mantissa >= 53;
  47    --  True if the floating-point type is at least IEEE Double
  48
  49    Precision_Limit : constant Uns := 2**Num'Machine_Mantissa - 1;
  50    --  See below for the rationale
  51
  52    package Impl is new Value_R (Uns, 2, Precision_Limit, Round => False);
  53
  54    subtype Base_T is Unsigned range 2 .. 16;
  55
  56    --  The following tables compute the maximum exponent of the base that can
  57    --  fit in the given floating-point format, that is to say the element at
  58    --  index N is the largest K such that N**K <= Num'Last.
  59
  60    Maxexp32 : constant array (Base_T) of Positive :=
  61      [2  => 127, 3 => 80,  4 => 63,  5 => 55,  6 => 49,
  62       7  => 45,  8 => 42,  9 => 40, 10 => 55, 11 => 37,
  63       12 => 35, 13 => 34, 14 => 33, 15 => 32, 16 => 31];
  64    --  The actual value for 10 is 38 but we also use scaling for 10
  65
  66    Maxexp64 : constant array (Base_T) of Positive :=
  67      [2  => 1023, 3 => 646,  4 => 511,  5 => 441,  6 => 396,
  68       7  => 364,  8 => 341,  9 => 323, 10 => 441, 11 => 296,
  69       12 => 285, 13 => 276, 14 => 268, 15 => 262, 16 => 255];
  70    --  The actual value for 10 is 308 but we also use scaling for 10
  71
  72    Maxexp80 : constant array (Base_T) of Positive :=
  73      [2  => 16383, 3 => 10337, 4 => 8191,  5 => 7056,  6 => 6338,
  74       7  => 5836,  8 => 5461,  9 => 5168, 10 => 7056, 11 => 4736,
  75       12 => 4570, 13 => 4427, 14 => 4303, 15 => 4193, 16 => 4095];
  76    --  The actual value for 10 is 4932 but we also use scaling for 10
  77
  78    package Double_Real is new System.Double_Real (Num);
  79    use type Double_Real.Double_T;
  80
  81    subtype Double_T is Double_Real.Double_T;
  82    --  The double floating-point type
  83
  84    function Exact_Log2 (N : Unsigned) return Positive is
  85      (case N is
  86         when  2     => 1,
  87         when  4     => 2,
  88         when  8     => 3,
  89         when 16     => 4,
  90         when others => raise Program_Error);
  91    --  Return the exponent of a power of 2
  92
  93    function Integer_to_Real
  94      (Str   : String;
  95       Val   : Impl.Value_Array;
  96       Base  : Unsigned;
  97       Scale : Impl.Scale_Array;
  98       Minus : Boolean) return Num;
  99    --  Convert the real value from integer to real representation
 100
 101    function Large_Powfive (Exp : Natural) return Double_T;
 102    --  Return 5.0**Exp as a double number, where Exp > Maxpow
 103
 104    function Large_Powfive (Exp : Natural; S : out Natural) return Double_T;
 105    --  Return Num'Scaling (5.0**Exp, -S) as a double number where Exp > Maxexp
 106
 107    ---------------------
 108    -- Integer_to_Real --
 109    ---------------------
 110
 111    function Integer_to_Real
 112      (Str   : String;
 113       Val   : Impl.Value_Array;
 114       Base  : Unsigned;
 115       Scale : Impl.Scale_Array;
 116       Minus : Boolean) return Num
 117    is
 118       pragma Assert (Base in 2 .. 16);
 119
 120       pragma Assert (Num'Machine_Radix = 2);
 121
 122       pragma Unsuppress (Range_Check);
 123
 124       Maxexp : constant Positive :=
 125                  (if    Num'Size = 32             then Maxexp32 (Base)
 126                   elsif Num'Size = 64             then Maxexp64 (Base)
 127                   elsif Num'Machine_Mantissa = 64 then Maxexp80 (Base)
 128                   else  raise Program_Error);
 129       --  Maximum exponent of the base that can fit in Num
 130
 131       D_Val : Double_T;
 132       R_Val : Num;
 133       S     : Integer;
 134
 135    begin
 136       --  We call the floating-point processor reset routine so we can be sure
 137       --  that the x87 FPU is properly set for conversions. This is especially
 138       --  needed on Windows, where calls to the operating system randomly reset
 139       --  the processor into 64-bit mode.
 140
 141       if Num'Machine_Mantissa = 64 then
 142          System.Float_Control.Reset;
 143       end if;
 144
 145       --  First convert the integer mantissa into a double real. The conversion
 146       --  of each part is exact, given the precision limit we used above. Then,
 147       --  if the contribution of the low part might be nonnull, scale the high
 148       --  part appropriately and add the low part to the result.
 149
 150       if Val (2) = 0 then
 151          D_Val := Double_Real.To_Double (Num (Val (1)));
 152          S := Scale (1);
 153
 154       else
 155          declare
 156             V1 : constant Num := Num (Val (1));
 157             V2 : constant Num := Num (Val (2));
 158
 159             DS : Positive;
 160
 161          begin
 162             DS := Scale (1) - Scale (2);
 163
 164             case Base is
 165                --  If the base is a power of two, we use the efficient Scaling
 166                --  attribute up to an amount worth a double mantissa.
 167
 168                when 2 | 4 | 8 | 16 =>
 169                   declare
 170                      L : constant Positive := Exact_Log2 (Base);
 171
 172                   begin
 173                      if DS <= 2 * Num'Machine_Mantissa / L then
 174                         DS := DS * L;
 175                         D_Val :=
 176                           Double_Real.Quick_Two_Sum (Num'Scaling (V1, DS), V2);
 177                         S := Scale (2);
 178
 179                      else
 180                         D_Val := Double_Real.To_Double (V1);
 181                         S := Scale (1);
 182                      end if;
 183                   end;
 184
 185                --  If the base is 10, we also scale up to an amount worth a
 186                --  double mantissa.
 187
 188                when 10 =>
 189                   declare
 190                      Powfive : constant array (0 .. Maxpow) of Double_T;
 191                      pragma Import (Ada, Powfive);
 192                      for Powfive'Address use Powfive_Address;
 193
 194                   begin
 195                      if DS <= Maxpow then
 196                         D_Val := Powfive (DS) * Num'Scaling (V1, DS) + V2;
 197                         S := Scale (2);
 198
 199                      else
 200                         D_Val := Double_Real.To_Double (V1);
 201                         S := Scale (1);
 202                      end if;
 203                   end;
 204
 205                --  Inaccurate implementation for other bases
 206
 207                when others =>
 208                   D_Val := Double_Real.To_Double (V1);
 209                   S := Scale (1);
 210             end case;
 211          end;
 212       end if;
 213
 214       --  Compute the final value by applying the scaling, if any
 215
 216       if (Val (1) = 0 and then Val (2) = 0) or else S = 0 then
 217          R_Val := Double_Real.To_Single (D_Val);
 218
 219       else
 220          case Base is
 221             --  If the base is a power of two, we use the efficient Scaling
 222             --  attribute with an overflow check, if it is not 2, to catch
 223             --  ludicrous exponents that would result in an infinity or zero.
 224
 225             when 2 | 4 | 8 | 16 =>
 226                declare
 227                   L : constant Positive := Exact_Log2 (Base);
 228
 229                begin
 230                   if Integer'First / L <= S and then S <= Integer'Last / L then
 231                      S := S * L;
 232                   end if;
 233
 234                   R_Val := Num'Scaling (Double_Real.To_Single (D_Val), S);
 235                end;
 236
 237             --  If the base is 10, we use a double implementation for the sake
 238             --  of accuracy combining powers of 5 and scaling attribute. Using
 239             --  this combination is better than using powers of 10 only because
 240             --  the Large_Powfive function may overflow only if the final value
 241             --  will also either overflow or underflow, thus making it possible
 242             --  to use a single division for the case of negative powers of 10.
 243
 244             when 10 =>
 245                declare
 246                   Powfive : constant array (0 .. Maxpow) of Double_T;
 247                   pragma Import (Ada, Powfive);
 248                   for Powfive'Address use Powfive_Address;
 249
 250                   RS : Natural;
 251
 252                begin
 253                   if S > 0 then
 254                      if S <= Maxpow then
 255                         D_Val := D_Val * Powfive (S);
 256                      else
 257                         D_Val := D_Val * Large_Powfive (S);
 258                      end if;
 259
 260                   else
 261                      if S >= -Maxpow then
 262                         D_Val := D_Val / Powfive (-S);
 263
 264                      --  For small types, typically IEEE Single, the trick
 265                      --  described above does not fully work.
 266
 267                      elsif not Is_Large_Type and then S < -Maxexp then
 268                         D_Val := D_Val / Large_Powfive (-S, RS);
 269                         S := S - RS;
 270
 271                      else
 272                         D_Val := D_Val / Large_Powfive (-S);
 273                      end if;
 274                   end if;
 275
 276                   R_Val := Num'Scaling (Double_Real.To_Single (D_Val), S);
 277                end;
 278
 279             --  Implementation for other bases with exponentiation
 280
 281             --  When the exponent is positive, we can do the computation
 282             --  directly because, if the exponentiation overflows, then
 283             --  the final value overflows as well. But when the exponent
 284             --  is negative, we may need to do it in two steps to avoid
 285             --  an artificial underflow.
 286
 287             when others =>
 288                declare
 289                   B : constant Num := Num (Base);
 290
 291                begin
 292                   R_Val := Double_Real.To_Single (D_Val);
 293
 294                   if S > 0 then
 295                      R_Val := R_Val * B ** S;
 296
 297                   else
 298                      if S < -Maxexp then
 299                         R_Val := R_Val / B ** Maxexp;
 300                         S := S + Maxexp;
 301                      end if;
 302
 303                      R_Val := R_Val / B ** (-S);
 304                   end if;
 305                end;
 306          end case;
 307       end if;
 308
 309       --  Finally deal with initial minus sign, note that this processing is
 310       --  done even if Uval is zero, so that -0.0 is correctly interpreted.
 311
 312       return (if Minus then -R_Val else R_Val);
 313
 314    exception
 315       when Constraint_Error => Bad_Value (Str);
 316    end Integer_to_Real;
 317
 318    -------------------
 319    -- Large_Powfive --
 320    -------------------
 321
 322    function Large_Powfive (Exp : Natural) return Double_T is
 323       Powfive : constant array (0 .. Maxpow) of Double_T;
 324       pragma Import (Ada, Powfive);
 325       for Powfive'Address use Powfive_Address;
 326
 327       Powfive_100 : constant Double_T;
 328       pragma Import (Ada, Powfive_100);
 329       for Powfive_100'Address use Powfive_100_Address;
 330
 331       Powfive_200 : constant Double_T;
 332       pragma Import (Ada, Powfive_200);
 333       for Powfive_200'Address use Powfive_200_Address;
 334
 335       Powfive_300 : constant Double_T;
 336       pragma Import (Ada, Powfive_300);
 337       for Powfive_300'Address use Powfive_300_Address;
 338
 339       R : Double_T;
 340       E : Natural;
 341
 342    begin
 343       pragma Assert (Exp > Maxpow);
 344
 345       if Is_Large_Type and then Exp >= 300 then
 346          R := Powfive_300;
 347          E := Exp - 300;
 348
 349       elsif Is_Large_Type and then Exp >= 200 then
 350          R := Powfive_200;
 351          E := Exp - 200;
 352
 353       elsif Is_Large_Type and then Exp >= 100 then
 354          R := Powfive_100;
 355          E := Exp - 100;
 356
 357       else
 358          R := Powfive (Maxpow);
 359          E := Exp - Maxpow;
 360       end if;
 361
 362       while E > Maxpow loop
 363          R := R * Powfive (Maxpow);
 364          E := E - Maxpow;
 365       end loop;
 366
 367       R := R * Powfive (E);
 368
 369       return R;
 370    end Large_Powfive;
 371
 372    function Large_Powfive (Exp : Natural; S : out Natural) return Double_T is
 373       Maxexp : constant Positive :=
 374         (if    Num'Size = 32             then Maxexp32 (5)
 375          elsif Num'Size = 64             then Maxexp64 (5)
 376          elsif Num'Machine_Mantissa = 64 then Maxexp80 (5)
 377          else  raise Program_Error);
 378       --  Maximum exponent of 5 that can fit in Num
 379
 380       Powfive : constant array (0 .. Maxpow) of Double_T;
 381       pragma Import (Ada, Powfive);
 382       for Powfive'Address use Powfive_Address;
 383
 384       R : Double_T;
 385       E : Natural;
 386
 387    begin
 388       pragma Assert (Exp > Maxexp);
 389
 390       pragma Warnings (Off, "-gnatw.a");
 391       pragma Assert (not Is_Large_Type);
 392       pragma Warnings (On, "-gnatw.a");
 393
 394       R := Powfive (Maxpow);
 395       E := Exp - Maxpow;
 396
 397       --  If the exponent is not too large, then scale down the result so that
 398       --  its final value does not overflow but, if it's too large, then do not
 399       --  bother doing it since overflow is just fine. The scaling factor is -3
 400       --  for every power of 5 above the maximum, in other words division by 8.
 401
 402       if Exp - Maxexp <= Maxpow then
 403          S := 3 * (Exp - Maxexp);
 404          R.Hi := Num'Scaling (R.Hi, -S);
 405          R.Lo := Num'Scaling (R.Lo, -S);
 406       else
 407          S := 0;
 408       end if;
 409
 410       while E > Maxpow loop
 411          R := R * Powfive (Maxpow);
 412          E := E - Maxpow;
 413       end loop;
 414
 415       R := R * Powfive (E);
 416
 417       return R;
 418    end Large_Powfive;
 419
 420    ---------------
 421    -- Scan_Real --
 422    ---------------
 423
 424    function Scan_Real
 425      (Str : String;
 426       Ptr : not null access Integer;
 427       Max : Integer) return Num
 428    is
 429       Base  : Unsigned;
 430       Scale : Impl.Scale_Array;
 431       Extra : Unsigned;
 432       Minus : Boolean;
 433       Val   : Impl.Value_Array;
 434
 435    begin
 436       Val := Impl.Scan_Raw_Real (Str, Ptr, Max, Base, Scale, Extra, Minus);
 437
 438       return Integer_to_Real (Str, Val, Base, Scale, Minus);
 439    end Scan_Real;
 440
 441    ----------------
 442    -- Value_Real --
 443    ----------------
 444
 445    function Value_Real (Str : String) return Num is
 446       Base  : Unsigned;
 447       Scale : Impl.Scale_Array;
 448       Extra : Unsigned;
 449       Minus : Boolean;
 450       Val   : Impl.Value_Array;
 451
 452    begin
 453       Val := Impl.Value_Raw_Real (Str, Base, Scale, Extra, Minus);
 454
 455       return Integer_to_Real (Str, Val, Base, Scale, Minus);
 456    end Value_Real;
 457
 458 end System.Val_Real;