1------------------------------------------------------------------------------
2--                                                                          --
3--                         GNAT RUN-TIME COMPONENTS                         --
4--                                                                          --
5--                       S Y S T E M . E X N _ L L F                        --
6--                                                                          --
7--                                 B o d y                                  --
8--                                                                          --
9--          Copyright (C) 1992-2015, 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: the reason for treating exponents in the range 0 .. 4 specially is
33--  to ensure identical results to the static inline expansion in the case of
34--  a compile time known exponent in this range. The use of Float'Machine and
35--  Long_Float'Machine is to avoid unwanted extra precision in the results.
36
37package body System.Exn_LLF is
38
39   function Exp
40     (Left  : Long_Long_Float;
41      Right : Integer) return Long_Long_Float;
42   --  Common routine used if Right not in 0 .. 4
43
44   ---------------
45   -- Exn_Float --
46   ---------------
47
48   function Exn_Float
49     (Left  : Float;
50      Right : Integer) return Float
51   is
52      Temp : Float;
53   begin
54      case Right is
55         when 0 =>
56            return 1.0;
57         when 1 =>
58            return Left;
59         when 2 =>
60            return Float'Machine (Left * Left);
61         when 3 =>
62            return Float'Machine (Left * Left * Left);
63         when 4 =>
64            Temp := Float'Machine (Left * Left);
65            return Float'Machine (Temp * Temp);
66         when others =>
67            return
68              Float'Machine
69                (Float (Exp (Long_Long_Float (Left), Right)));
70      end case;
71   end Exn_Float;
72
73   --------------------
74   -- Exn_Long_Float --
75   --------------------
76
77   function Exn_Long_Float
78     (Left  : Long_Float;
79      Right : Integer) return Long_Float
80   is
81      Temp : Long_Float;
82   begin
83      case Right is
84         when 0 =>
85            return 1.0;
86         when 1 =>
87            return Left;
88         when 2 =>
89            return Long_Float'Machine (Left * Left);
90         when 3 =>
91            return Long_Float'Machine (Left * Left * Left);
92         when 4 =>
93            Temp := Long_Float'Machine (Left * Left);
94            return Long_Float'Machine (Temp * Temp);
95         when others =>
96            return
97              Long_Float'Machine
98                (Long_Float (Exp (Long_Long_Float (Left), Right)));
99      end case;
100   end Exn_Long_Float;
101
102   -------------------------
103   -- Exn_Long_Long_Float --
104   -------------------------
105
106   function Exn_Long_Long_Float
107     (Left  : Long_Long_Float;
108      Right : Integer) return Long_Long_Float
109   is
110      Temp : Long_Long_Float;
111   begin
112      case Right is
113         when 0 =>
114            return 1.0;
115         when 1 =>
116            return Left;
117         when 2 =>
118            return Left * Left;
119         when 3 =>
120            return Left * Left * Left;
121         when 4 =>
122            Temp := Left * Left;
123            return Temp * Temp;
124         when others =>
125            return Exp (Left, Right);
126      end case;
127   end Exn_Long_Long_Float;
128
129   ---------
130   -- Exp --
131   ---------
132
133   function Exp
134     (Left  : Long_Long_Float;
135      Right : Integer) return Long_Long_Float
136   is
137      Result : Long_Long_Float := 1.0;
138      Factor : Long_Long_Float := Left;
139      Exp    : Integer := Right;
140
141   begin
142      --  We use the standard logarithmic approach, Exp gets shifted right
143      --  testing successive low order bits and Factor is the value of the
144      --  base raised to the next power of 2. If the low order bit or Exp is
145      --  set, multiply the result by this factor. For negative exponents,
146      --  invert result upon return.
147
148      if Exp >= 0 then
149         loop
150            if Exp rem 2 /= 0 then
151               Result := Result * Factor;
152            end if;
153
154            Exp := Exp / 2;
155            exit when Exp = 0;
156            Factor := Factor * Factor;
157         end loop;
158
159         return Result;
160
161      --  Here we have a negative exponent, and we compute the result as:
162
163      --     1.0 / (Left ** (-Right))
164
165      --  Note that the case of Left being zero is not special, it will
166      --  simply result in a division by zero at the end, yielding a
167      --  correctly signed infinity, or possibly generating an overflow.
168
169      --  Note on overflow: The coding of this routine assumes that the
170      --  target generates infinities with standard IEEE semantics. If this
171      --  is not the case, then the code below may raise Constraint_Error.
172      --  This follows the implementation permission given in RM 4.5.6(12).
173
174      else
175         begin
176            loop
177               if Exp rem 2 /= 0 then
178                  Result := Result * Factor;
179               end if;
180
181               Exp := Exp / 2;
182               exit when Exp = 0;
183               Factor := Factor * Factor;
184            end loop;
185
186            return 1.0 / Result;
187         end;
188      end if;
189   end Exp;
190
191end System.Exn_LLF;
192