gcc/ada/libgnat/s-exponr.adb

   1 ------------------------------------------------------------------------------
   2 --                                                                          --
   3 --                         GNAT RUN-TIME COMPONENTS                         --
   4 --                                                                          --
   5 --                        S Y S T E M . E X P O N R                         --
   6 --                                                                          --
   7 --                                 B o d y                                  --
   8 --                                                                          --
   9 --            Copyright (C) 2021-2023, 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 --  Note that the reason for treating exponents in the range 0 .. 4 specially
  33 --  is to ensure identical results with the static expansion in the case of a
  34 --  compile-time known exponent in this range; similarly, the use 'Machine is
  35 --  to avoid unwanted extra precision in the results.
  36
  37 --  For a negative exponent, we compute the result as per RM 4.5.6(11/3):
  38
  39 --     Left ** Right = 1.0 / (Left ** (-Right))
  40
  41 --  Note that the case of Left being zero is not special, it will simply result
  42 --  in a division by zero at the end, yielding a correctly signed infinity, or
  43 --  possibly raising an overflow exception.
  44
  45 --  Note on overflow: this coding assumes that the target generates infinities
  46 --  with standard IEEE semantics. If this is not the case, then the code for
  47 --  negative exponents may raise Constraint_Error, which is in keeping with the
  48 --  implementation permission given in RM 4.5.6(12).
  49
  50 with System.Double_Real;
  51
  52 function System.Exponr (Left : Num; Right : Integer) return Num is
  53
  54    package Double_Real is new System.Double_Real (Num);
  55    use type Double_Real.Double_T;
  56
  57    subtype Double_T is Double_Real.Double_T;
  58    --  The double floating-point type
  59
  60    subtype Safe_Negative is Integer range Integer'First + 1 .. -1;
  61    --  The range of safe negative exponents
  62
  63    function Expon (Left : Num; Right : Natural) return Num;
  64    --  Routine used if Right is greater than 4
  65
  66    -----------
  67    -- Expon --
  68    -----------
  69
  70    function Expon (Left : Num; Right : Natural) return Num is
  71       Result : Double_T := Double_Real.To_Double (1.0);
  72       Factor : Double_T := Double_Real.To_Double (Left);
  73       Exp    : Natural  := Right;
  74
  75    begin
  76       --  We use the standard logarithmic approach, Exp gets shifted right
  77       --  testing successive low order bits and Factor is the value of the
  78       --  base raised to the next power of 2. If the low order bit or Exp
  79       --  is set, multiply the result by this factor.
  80
  81       loop
  82          if Exp rem 2 /= 0 then
  83             Result := Result * Factor;
  84             exit when Exp = 1;
  85          end if;
  86
  87          Exp := Exp / 2;
  88          Factor := Double_Real.Sqr (Factor);
  89       end loop;
  90
  91       return Double_Real.To_Single (Result);
  92    end Expon;
  93
  94 begin
  95    case Right is
  96       when 0 =>
  97          return 1.0;
  98
  99       when 1 =>
 100          return Left;
 101
 102       when 2 =>
 103          return Num'Machine (Left * Left);
 104
 105       when 3 =>
 106          return Num'Machine (Left * Left * Left);
 107
 108       when 4 =>
 109          declare
 110             Sqr : constant Num := Num'Machine (Left * Left);
 111
 112          begin
 113             return Num'Machine (Sqr * Sqr);
 114          end;
 115
 116       when Safe_Negative =>
 117          return Num'Machine (1.0 / Exponr (Left, -Right));
 118
 119       when Integer'First =>
 120          return Num'Machine (1.0 / (Exponr (Left, Integer'Last) * Left));
 121
 122       when others =>
 123          return Num'Machine (Expon (Left, Right));
 124    end case;
 125 end System.Exponr;