OpenCores
URL https://opencores.org/ocsvn/openrisc/openrisc/trunk

Subversion Repositories openrisc

[/] [openrisc/] [trunk/] [gnu-dev/] [or1k-gcc/] [gcc/] [ada/] [s-fatgen.adb] - Blame information for rev 706

Details | Compare with Previous | View Log

Line No. Rev Author Line
1 706 jeremybenn
------------------------------------------------------------------------------
2
--                                                                          --
3
--                         GNAT COMPILER COMPONENTS                         --
4
--                                                                          --
5
--                       S Y S T E M . F A T _ G E N                        --
6
--                                                                          --
7
--                                 B o d y                                  --
8
--                                                                          --
9
--          Copyright (C) 1992-2010, 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
--  The implementation here is portable to any IEEE implementation. It does
33
--  not handle non-binary radix, and also assumes that model numbers and
34
--  machine numbers are basically identical, which is not true of all possible
35
--  floating-point implementations. On a non-IEEE machine, this body must be
36
--  specialized appropriately, or better still, its generic instantiations
37
--  should be replaced by efficient machine-specific code.
38
 
39
with Ada.Unchecked_Conversion;
40
with System;
41
package body System.Fat_Gen is
42
 
43
   Float_Radix        : constant T := T (T'Machine_Radix);
44
   Radix_To_M_Minus_1 : constant T := Float_Radix ** (T'Machine_Mantissa - 1);
45
 
46
   pragma Assert (T'Machine_Radix = 2);
47
   --  This version does not handle radix 16
48
 
49
   --  Constants for Decompose and Scaling
50
 
51
   Rad    : constant T := T (T'Machine_Radix);
52
   Invrad : constant T := 1.0 / Rad;
53
 
54
   subtype Expbits is Integer range 0 .. 6;
55
   --  2 ** (2 ** 7) might overflow.  How big can radix-16 exponents get?
56
 
57
   Log_Power : constant array (Expbits) of Integer := (1, 2, 4, 8, 16, 32, 64);
58
 
59
   R_Power : constant array (Expbits) of T :=
60
     (Rad **  1,
61
      Rad **  2,
62
      Rad **  4,
63
      Rad **  8,
64
      Rad ** 16,
65
      Rad ** 32,
66
      Rad ** 64);
67
 
68
   R_Neg_Power : constant array (Expbits) of T :=
69
     (Invrad **  1,
70
      Invrad **  2,
71
      Invrad **  4,
72
      Invrad **  8,
73
      Invrad ** 16,
74
      Invrad ** 32,
75
      Invrad ** 64);
76
 
77
   -----------------------
78
   -- Local Subprograms --
79
   -----------------------
80
 
81
   procedure Decompose (XX : T; Frac : out T; Expo : out UI);
82
   --  Decomposes a floating-point number into fraction and exponent parts.
83
   --  Both results are signed, with Frac having the sign of XX, and UI has
84
   --  the sign of the exponent. The absolute value of Frac is in the range
85
   --  0.0 <= Frac < 1.0. If Frac = 0.0 or -0.0, then Expo is always zero.
86
 
87
   function Gradual_Scaling  (Adjustment : UI) return T;
88
   --  Like Scaling with a first argument of 1.0, but returns the smallest
89
   --  denormal rather than zero when the adjustment is smaller than
90
   --  Machine_Emin. Used for Succ and Pred.
91
 
92
   --------------
93
   -- Adjacent --
94
   --------------
95
 
96
   function Adjacent (X, Towards : T) return T is
97
   begin
98
      if Towards = X then
99
         return X;
100
      elsif Towards > X then
101
         return Succ (X);
102
      else
103
         return Pred (X);
104
      end if;
105
   end Adjacent;
106
 
107
   -------------
108
   -- Ceiling --
109
   -------------
110
 
111
   function Ceiling (X : T) return T is
112
      XT : constant T := Truncation (X);
113
   begin
114
      if X <= 0.0 then
115
         return XT;
116
      elsif X = XT then
117
         return X;
118
      else
119
         return XT + 1.0;
120
      end if;
121
   end Ceiling;
122
 
123
   -------------
124
   -- Compose --
125
   -------------
126
 
127
   function Compose (Fraction : T; Exponent : UI) return T is
128
      Arg_Frac : T;
129
      Arg_Exp  : UI;
130
      pragma Unreferenced (Arg_Exp);
131
   begin
132
      Decompose (Fraction, Arg_Frac, Arg_Exp);
133
      return Scaling (Arg_Frac, Exponent);
134
   end Compose;
135
 
136
   ---------------
137
   -- Copy_Sign --
138
   ---------------
139
 
140
   function Copy_Sign (Value, Sign : T) return T is
141
      Result : T;
142
 
143
      function Is_Negative (V : T) return Boolean;
144
      pragma Import (Intrinsic, Is_Negative);
145
 
146
   begin
147
      Result := abs Value;
148
 
149
      if Is_Negative (Sign) then
150
         return -Result;
151
      else
152
         return Result;
153
      end if;
154
   end Copy_Sign;
155
 
156
   ---------------
157
   -- Decompose --
158
   ---------------
159
 
160
   procedure Decompose (XX : T; Frac : out T; Expo : out UI) is
161
      X : constant T := T'Machine (XX);
162
 
163
   begin
164
      if X = 0.0 then
165
 
166
         --  The normalized exponent of zero is zero, see RM A.5.2(15)
167
 
168
         Frac := X;
169
         Expo := 0;
170
 
171
      --  Check for infinities, transfinites, whatnot
172
 
173
      elsif X > T'Safe_Last then
174
         Frac := Invrad;
175
         Expo := T'Machine_Emax + 1;
176
 
177
      elsif X < T'Safe_First then
178
         Frac := -Invrad;
179
         Expo := T'Machine_Emax + 2;    -- how many extra negative values?
180
 
181
      else
182
         --  Case of nonzero finite x. Essentially, we just multiply
183
         --  by Rad ** (+-2**N) to reduce the range.
184
 
185
         declare
186
            Ax : T  := abs X;
187
            Ex : UI := 0;
188
 
189
         --  Ax * Rad ** Ex is invariant
190
 
191
         begin
192
            if Ax >= 1.0 then
193
               while Ax >= R_Power (Expbits'Last) loop
194
                  Ax := Ax * R_Neg_Power (Expbits'Last);
195
                  Ex := Ex + Log_Power (Expbits'Last);
196
               end loop;
197
 
198
               --  Ax < Rad ** 64
199
 
200
               for N in reverse Expbits'First .. Expbits'Last - 1 loop
201
                  if Ax >= R_Power (N) then
202
                     Ax := Ax * R_Neg_Power (N);
203
                     Ex := Ex + Log_Power (N);
204
                  end if;
205
 
206
                  --  Ax < R_Power (N)
207
 
208
               end loop;
209
 
210
               --  1 <= Ax < Rad
211
 
212
               Ax := Ax * Invrad;
213
               Ex := Ex + 1;
214
 
215
            else
216
               --  0 < ax < 1
217
 
218
               while Ax < R_Neg_Power (Expbits'Last) loop
219
                  Ax := Ax * R_Power (Expbits'Last);
220
                  Ex := Ex - Log_Power (Expbits'Last);
221
               end loop;
222
 
223
               --  Rad ** -64 <= Ax < 1
224
 
225
               for N in reverse Expbits'First .. Expbits'Last - 1 loop
226
                  if Ax < R_Neg_Power (N) then
227
                     Ax := Ax * R_Power (N);
228
                     Ex := Ex - Log_Power (N);
229
                  end if;
230
 
231
                  --  R_Neg_Power (N) <= Ax < 1
232
 
233
               end loop;
234
            end if;
235
 
236
            Frac := (if X > 0.0 then Ax else -Ax);
237
            Expo := Ex;
238
         end;
239
      end if;
240
   end Decompose;
241
 
242
   --------------
243
   -- Exponent --
244
   --------------
245
 
246
   function Exponent (X : T) return UI is
247
      X_Frac : T;
248
      X_Exp  : UI;
249
      pragma Unreferenced (X_Frac);
250
   begin
251
      Decompose (X, X_Frac, X_Exp);
252
      return X_Exp;
253
   end Exponent;
254
 
255
   -----------
256
   -- Floor --
257
   -----------
258
 
259
   function Floor (X : T) return T is
260
      XT : constant T := Truncation (X);
261
   begin
262
      if X >= 0.0 then
263
         return XT;
264
      elsif XT = X then
265
         return X;
266
      else
267
         return XT - 1.0;
268
      end if;
269
   end Floor;
270
 
271
   --------------
272
   -- Fraction --
273
   --------------
274
 
275
   function Fraction (X : T) return T is
276
      X_Frac : T;
277
      X_Exp  : UI;
278
      pragma Unreferenced (X_Exp);
279
   begin
280
      Decompose (X, X_Frac, X_Exp);
281
      return X_Frac;
282
   end Fraction;
283
 
284
   ---------------------
285
   -- Gradual_Scaling --
286
   ---------------------
287
 
288
   function Gradual_Scaling  (Adjustment : UI) return T is
289
      Y  : T;
290
      Y1 : T;
291
      Ex : UI := Adjustment;
292
 
293
   begin
294
      if Adjustment < T'Machine_Emin - 1 then
295
         Y  := 2.0 ** T'Machine_Emin;
296
         Y1 := Y;
297
         Ex := Ex - T'Machine_Emin;
298
         while Ex < 0 loop
299
            Y := T'Machine (Y / 2.0);
300
 
301
            if Y = 0.0 then
302
               return Y1;
303
            end if;
304
 
305
            Ex := Ex + 1;
306
            Y1 := Y;
307
         end loop;
308
 
309
         return Y1;
310
 
311
      else
312
         return Scaling (1.0, Adjustment);
313
      end if;
314
   end Gradual_Scaling;
315
 
316
   ------------------
317
   -- Leading_Part --
318
   ------------------
319
 
320
   function Leading_Part (X : T; Radix_Digits : UI) return T is
321
      L    : UI;
322
      Y, Z : T;
323
 
324
   begin
325
      if Radix_Digits >= T'Machine_Mantissa then
326
         return X;
327
 
328
      elsif Radix_Digits <= 0 then
329
         raise Constraint_Error;
330
 
331
      else
332
         L := Exponent (X) - Radix_Digits;
333
         Y := Truncation (Scaling (X, -L));
334
         Z := Scaling (Y, L);
335
         return Z;
336
      end if;
337
   end Leading_Part;
338
 
339
   -------------
340
   -- Machine --
341
   -------------
342
 
343
   --  The trick with Machine is to force the compiler to store the result
344
   --  in memory so that we do not have extra precision used. The compiler
345
   --  is clever, so we have to outwit its possible optimizations! We do
346
   --  this by using an intermediate pragma Volatile location.
347
 
348
   function Machine (X : T) return T is
349
      Temp : T;
350
      pragma Volatile (Temp);
351
   begin
352
      Temp := X;
353
      return Temp;
354
   end Machine;
355
 
356
   ----------------------
357
   -- Machine_Rounding --
358
   ----------------------
359
 
360
   --  For now, the implementation is identical to that of Rounding, which is
361
   --  a permissible behavior, but is not the most efficient possible approach.
362
 
363
   function Machine_Rounding (X : T) return T is
364
      Result : T;
365
      Tail   : T;
366
 
367
   begin
368
      Result := Truncation (abs X);
369
      Tail   := abs X - Result;
370
 
371
      if Tail >= 0.5  then
372
         Result := Result + 1.0;
373
      end if;
374
 
375
      if X > 0.0 then
376
         return Result;
377
 
378
      elsif X < 0.0 then
379
         return -Result;
380
 
381
      --  For zero case, make sure sign of zero is preserved
382
 
383
      else
384
         return X;
385
      end if;
386
   end Machine_Rounding;
387
 
388
   -----------
389
   -- Model --
390
   -----------
391
 
392
   --  We treat Model as identical to Machine. This is true of IEEE and other
393
   --  nice floating-point systems, but not necessarily true of all systems.
394
 
395
   function Model (X : T) return T is
396
   begin
397
      return Machine (X);
398
   end Model;
399
 
400
   ----------
401
   -- Pred --
402
   ----------
403
 
404
   --  Subtract from the given number a number equivalent to the value of its
405
   --  least significant bit. Given that the most significant bit represents
406
   --  a value of 1.0 * radix ** (exp - 1), the value we want is obtained by
407
   --  shifting this by (mantissa-1) bits to the right, i.e. decreasing the
408
   --  exponent by that amount.
409
 
410
   --  Zero has to be treated specially, since its exponent is zero
411
 
412
   function Pred (X : T) return T is
413
      X_Frac : T;
414
      X_Exp  : UI;
415
 
416
   begin
417
      if X = 0.0 then
418
         return -Succ (X);
419
 
420
      else
421
         Decompose (X, X_Frac, X_Exp);
422
 
423
         --  A special case, if the number we had was a positive power of
424
         --  two, then we want to subtract half of what we would otherwise
425
         --  subtract, since the exponent is going to be reduced.
426
 
427
         --  Note that X_Frac has the same sign as X, so if X_Frac is 0.5,
428
         --  then we know that we have a positive number (and hence a
429
         --  positive power of 2).
430
 
431
         if X_Frac = 0.5 then
432
            return X - Gradual_Scaling (X_Exp - T'Machine_Mantissa - 1);
433
 
434
         --  Otherwise the exponent is unchanged
435
 
436
         else
437
            return X - Gradual_Scaling (X_Exp - T'Machine_Mantissa);
438
         end if;
439
      end if;
440
   end Pred;
441
 
442
   ---------------
443
   -- Remainder --
444
   ---------------
445
 
446
   function Remainder (X, Y : T) return T is
447
      A        : T;
448
      B        : T;
449
      Arg      : T;
450
      P        : T;
451
      P_Frac   : T;
452
      Sign_X   : T;
453
      IEEE_Rem : T;
454
      Arg_Exp  : UI;
455
      P_Exp    : UI;
456
      K        : UI;
457
      P_Even   : Boolean;
458
 
459
      Arg_Frac : T;
460
      pragma Unreferenced (Arg_Frac);
461
 
462
   begin
463
      if Y = 0.0 then
464
         raise Constraint_Error;
465
      end if;
466
 
467
      if X > 0.0 then
468
         Sign_X :=  1.0;
469
         Arg := X;
470
      else
471
         Sign_X := -1.0;
472
         Arg := -X;
473
      end if;
474
 
475
      P := abs Y;
476
 
477
      if Arg < P then
478
         P_Even := True;
479
         IEEE_Rem := Arg;
480
         P_Exp := Exponent (P);
481
 
482
      else
483
         Decompose (Arg, Arg_Frac, Arg_Exp);
484
         Decompose (P,   P_Frac,   P_Exp);
485
 
486
         P := Compose (P_Frac, Arg_Exp);
487
         K := Arg_Exp - P_Exp;
488
         P_Even := True;
489
         IEEE_Rem := Arg;
490
 
491
         for Cnt in reverse 0 .. K loop
492
            if IEEE_Rem >= P then
493
               P_Even := False;
494
               IEEE_Rem := IEEE_Rem - P;
495
            else
496
               P_Even := True;
497
            end if;
498
 
499
            P := P * 0.5;
500
         end loop;
501
      end if;
502
 
503
      --  That completes the calculation of modulus remainder. The final
504
      --  step is get the IEEE remainder. Here we need to compare Rem with
505
      --  (abs Y) / 2. We must be careful of unrepresentable Y/2 value
506
      --  caused by subnormal numbers
507
 
508
      if P_Exp >= 0 then
509
         A := IEEE_Rem;
510
         B := abs Y * 0.5;
511
 
512
      else
513
         A := IEEE_Rem * 2.0;
514
         B := abs Y;
515
      end if;
516
 
517
      if A > B or else (A = B and then not P_Even) then
518
         IEEE_Rem := IEEE_Rem - abs Y;
519
      end if;
520
 
521
      return Sign_X * IEEE_Rem;
522
   end Remainder;
523
 
524
   --------------
525
   -- Rounding --
526
   --------------
527
 
528
   function Rounding (X : T) return T is
529
      Result : T;
530
      Tail   : T;
531
 
532
   begin
533
      Result := Truncation (abs X);
534
      Tail   := abs X - Result;
535
 
536
      if Tail >= 0.5  then
537
         Result := Result + 1.0;
538
      end if;
539
 
540
      if X > 0.0 then
541
         return Result;
542
 
543
      elsif X < 0.0 then
544
         return -Result;
545
 
546
      --  For zero case, make sure sign of zero is preserved
547
 
548
      else
549
         return X;
550
      end if;
551
   end Rounding;
552
 
553
   -------------
554
   -- Scaling --
555
   -------------
556
 
557
   --  Return x * rad ** adjustment quickly, or quietly underflow to zero,
558
   --  or overflow naturally.
559
 
560
   function Scaling (X : T; Adjustment : UI) return T is
561
   begin
562
      if X = 0.0 or else Adjustment = 0 then
563
         return X;
564
      end if;
565
 
566
      --  Nonzero x essentially, just multiply repeatedly by Rad ** (+-2**n)
567
 
568
      declare
569
         Y  : T  := X;
570
         Ex : UI := Adjustment;
571
 
572
      --  Y * Rad ** Ex is invariant
573
 
574
      begin
575
         if Ex < 0 then
576
            while Ex <= -Log_Power (Expbits'Last) loop
577
               Y := Y * R_Neg_Power (Expbits'Last);
578
               Ex := Ex + Log_Power (Expbits'Last);
579
            end loop;
580
 
581
            --  -64 < Ex <= 0
582
 
583
            for N in reverse Expbits'First .. Expbits'Last - 1 loop
584
               if Ex <= -Log_Power (N) then
585
                  Y := Y * R_Neg_Power (N);
586
                  Ex := Ex + Log_Power (N);
587
               end if;
588
 
589
               --  -Log_Power (N) < Ex <= 0
590
 
591
            end loop;
592
 
593
            --  Ex = 0
594
 
595
         else
596
            --  Ex >= 0
597
 
598
            while Ex >= Log_Power (Expbits'Last) loop
599
               Y := Y * R_Power (Expbits'Last);
600
               Ex := Ex - Log_Power (Expbits'Last);
601
            end loop;
602
 
603
            --  0 <= Ex < 64
604
 
605
            for N in reverse Expbits'First .. Expbits'Last - 1 loop
606
               if Ex >= Log_Power (N) then
607
                  Y := Y * R_Power (N);
608
                  Ex := Ex - Log_Power (N);
609
               end if;
610
 
611
               --  0 <= Ex < Log_Power (N)
612
 
613
            end loop;
614
 
615
            --  Ex = 0
616
 
617
         end if;
618
 
619
         return Y;
620
      end;
621
   end Scaling;
622
 
623
   ----------
624
   -- Succ --
625
   ----------
626
 
627
   --  Similar computation to that of Pred: find value of least significant
628
   --  bit of given number, and add. Zero has to be treated specially since
629
   --  the exponent can be zero, and also we want the smallest denormal if
630
   --  denormals are supported.
631
 
632
   function Succ (X : T) return T is
633
      X_Frac : T;
634
      X_Exp  : UI;
635
      X1, X2 : T;
636
 
637
   begin
638
      if X = 0.0 then
639
         X1 := 2.0 ** T'Machine_Emin;
640
 
641
         --  Following loop generates smallest denormal
642
 
643
         loop
644
            X2 := T'Machine (X1 / 2.0);
645
            exit when X2 = 0.0;
646
            X1 := X2;
647
         end loop;
648
 
649
         return X1;
650
 
651
      else
652
         Decompose (X, X_Frac, X_Exp);
653
 
654
         --  A special case, if the number we had was a negative power of two,
655
         --  then we want to add half of what we would otherwise add, since the
656
         --  exponent is going to be reduced.
657
 
658
         --  Note that X_Frac has the same sign as X, so if X_Frac is -0.5,
659
         --  then we know that we have a negative number (and hence a negative
660
         --  power of 2).
661
 
662
         if X_Frac = -0.5 then
663
            return X + Gradual_Scaling (X_Exp - T'Machine_Mantissa - 1);
664
 
665
         --  Otherwise the exponent is unchanged
666
 
667
         else
668
            return X + Gradual_Scaling (X_Exp - T'Machine_Mantissa);
669
         end if;
670
      end if;
671
   end Succ;
672
 
673
   ----------------
674
   -- Truncation --
675
   ----------------
676
 
677
   --  The basic approach is to compute
678
 
679
   --    T'Machine (RM1 + N) - RM1
680
 
681
   --  where N >= 0.0 and RM1 = radix ** (mantissa - 1)
682
 
683
   --  This works provided that the intermediate result (RM1 + N) does not
684
   --  have extra precision (which is why we call Machine). When we compute
685
   --  RM1 + N, the exponent of N will be normalized and the mantissa shifted
686
   --  shifted appropriately so the lower order bits, which cannot contribute
687
   --  to the integer part of N, fall off on the right. When we subtract RM1
688
   --  again, the significant bits of N are shifted to the left, and what we
689
   --  have is an integer, because only the first e bits are different from
690
   --  zero (assuming binary radix here).
691
 
692
   function Truncation (X : T) return T is
693
      Result : T;
694
 
695
   begin
696
      Result := abs X;
697
 
698
      if Result >= Radix_To_M_Minus_1 then
699
         return Machine (X);
700
 
701
      else
702
         Result := Machine (Radix_To_M_Minus_1 + Result) - Radix_To_M_Minus_1;
703
 
704
         if Result > abs X  then
705
            Result := Result - 1.0;
706
         end if;
707
 
708
         if X > 0.0 then
709
            return  Result;
710
 
711
         elsif X < 0.0 then
712
            return -Result;
713
 
714
         --  For zero case, make sure sign of zero is preserved
715
 
716
         else
717
            return X;
718
         end if;
719
      end if;
720
   end Truncation;
721
 
722
   -----------------------
723
   -- Unbiased_Rounding --
724
   -----------------------
725
 
726
   function Unbiased_Rounding (X : T) return T is
727
      Abs_X  : constant T := abs X;
728
      Result : T;
729
      Tail   : T;
730
 
731
   begin
732
      Result := Truncation (Abs_X);
733
      Tail   := Abs_X - Result;
734
 
735
      if Tail > 0.5  then
736
         Result := Result + 1.0;
737
 
738
      elsif Tail = 0.5 then
739
         Result := 2.0 * Truncation ((Result / 2.0) + 0.5);
740
      end if;
741
 
742
      if X > 0.0 then
743
         return Result;
744
 
745
      elsif X < 0.0 then
746
         return -Result;
747
 
748
      --  For zero case, make sure sign of zero is preserved
749
 
750
      else
751
         return X;
752
      end if;
753
   end Unbiased_Rounding;
754
 
755
   -----------
756
   -- Valid --
757
   -----------
758
 
759
   --  Note: this routine does not work for VAX float. We compensate for this
760
   --  in Exp_Attr by using the Valid functions in Vax_Float_Operations rather
761
   --  than the corresponding instantiation of this function.
762
 
763
   function Valid (X : not null access T) return Boolean is
764
 
765
      IEEE_Emin : constant Integer := T'Machine_Emin - 1;
766
      IEEE_Emax : constant Integer := T'Machine_Emax - 1;
767
 
768
      IEEE_Bias : constant Integer := -(IEEE_Emin - 1);
769
 
770
      subtype IEEE_Exponent_Range is
771
        Integer range IEEE_Emin - 1 .. IEEE_Emax + 1;
772
 
773
      --  The implementation of this floating point attribute uses a
774
      --  representation type Float_Rep that allows direct access to the
775
      --  exponent and mantissa parts of a floating point number.
776
 
777
      --  The Float_Rep type is an array of Float_Word elements. This
778
      --  representation is chosen to make it possible to size the type based
779
      --  on a generic parameter. Since the array size is known at compile
780
      --  time, efficient code can still be generated. The size of Float_Word
781
      --  elements should be large enough to allow accessing the exponent in
782
      --  one read, but small enough so that all floating point object sizes
783
      --  are a multiple of the Float_Word'Size.
784
 
785
      --  The following conditions must be met for all possible instantiations
786
      --  of the attributes package:
787
 
788
      --    - T'Size is an integral multiple of Float_Word'Size
789
 
790
      --    - The exponent and sign are completely contained in a single
791
      --      component of Float_Rep, named Most_Significant_Word (MSW).
792
 
793
      --    - The sign occupies the most significant bit of the MSW and the
794
      --      exponent is in the following bits. Unused bits (if any) are in
795
      --      the least significant part.
796
 
797
      type Float_Word is mod 2**Positive'Min (System.Word_Size, 32);
798
      type Rep_Index is range 0 .. 7;
799
 
800
      Rep_Words : constant Positive :=
801
                    (T'Size + Float_Word'Size - 1) / Float_Word'Size;
802
      Rep_Last  : constant Rep_Index :=
803
                    Rep_Index'Min
804
                      (Rep_Index (Rep_Words - 1),
805
                       (T'Mantissa + 16) / Float_Word'Size);
806
      --  Determine the number of Float_Words needed for representing the
807
      --  entire floating-point value. Do not take into account excessive
808
      --  padding, as occurs on IA-64 where 80 bits floats get padded to 128
809
      --  bits. In general, the exponent field cannot be larger than 15 bits,
810
      --  even for 128-bit floating-point types, so the final format size
811
      --  won't be larger than T'Mantissa + 16.
812
 
813
      type Float_Rep is
814
         array (Rep_Index range 0 .. Rep_Index (Rep_Words - 1)) of Float_Word;
815
 
816
      pragma Suppress_Initialization (Float_Rep);
817
      --  This pragma suppresses the generation of an initialization procedure
818
      --  for type Float_Rep when operating in Initialize/Normalize_Scalars
819
      --  mode. This is not just a matter of efficiency, but of functionality,
820
      --  since Valid has a pragma Inline_Always, which is not permitted if
821
      --  there are nested subprograms present.
822
 
823
      Most_Significant_Word : constant Rep_Index :=
824
                                Rep_Last * Standard'Default_Bit_Order;
825
      --  Finding the location of the Exponent_Word is a bit tricky. In general
826
      --  we assume Word_Order = Bit_Order. This expression needs to be refined
827
      --  for VMS.
828
 
829
      Exponent_Factor : constant Float_Word :=
830
                          2**(Float_Word'Size - 1) /
831
                            Float_Word (IEEE_Emax - IEEE_Emin + 3) *
832
                              Boolean'Pos (Most_Significant_Word /= 2) +
833
                                Boolean'Pos (Most_Significant_Word = 2);
834
      --  Factor that the extracted exponent needs to be divided by to be in
835
      --  range 0 .. IEEE_Emax - IEEE_Emin + 2. Special kludge: Exponent_Factor
836
      --  is 1 for x86/IA64 double extended as GCC adds unused bits to the
837
      --  type.
838
 
839
      Exponent_Mask : constant Float_Word :=
840
                        Float_Word (IEEE_Emax - IEEE_Emin + 2) *
841
                          Exponent_Factor;
842
      --  Value needed to mask out the exponent field. This assumes that the
843
      --  range IEEE_Emin - 1 .. IEEE_Emax + contains 2**N values, for some N
844
      --  in Natural.
845
 
846
      function To_Float is new Ada.Unchecked_Conversion (Float_Rep, T);
847
 
848
      type Float_Access is access all T;
849
      function To_Address is
850
         new Ada.Unchecked_Conversion (Float_Access, System.Address);
851
 
852
      XA : constant System.Address := To_Address (Float_Access (X));
853
 
854
      R : Float_Rep;
855
      pragma Import (Ada, R);
856
      for R'Address use XA;
857
      --  R is a view of the input floating-point parameter. Note that we
858
      --  must avoid copying the actual bits of this parameter in float
859
      --  form (since it may be a signalling NaN.
860
 
861
      E  : constant IEEE_Exponent_Range :=
862
             Integer ((R (Most_Significant_Word) and Exponent_Mask) /
863
                                                        Exponent_Factor)
864
               - IEEE_Bias;
865
      --  Mask/Shift T to only get bits from the exponent. Then convert biased
866
      --  value to integer value.
867
 
868
      SR : Float_Rep;
869
      --  Float_Rep representation of significant of X.all
870
 
871
   begin
872
      if T'Denorm then
873
 
874
         --  All denormalized numbers are valid, so the only invalid numbers
875
         --  are overflows and NaNs, both with exponent = Emax + 1.
876
 
877
         return E /= IEEE_Emax + 1;
878
 
879
      end if;
880
 
881
      --  All denormalized numbers except 0.0 are invalid
882
 
883
      --  Set exponent of X to zero, so we end up with the significand, which
884
      --  definitely is a valid number and can be converted back to a float.
885
 
886
      SR := R;
887
      SR (Most_Significant_Word) :=
888
           (SR (Most_Significant_Word)
889
             and not Exponent_Mask) + Float_Word (IEEE_Bias) * Exponent_Factor;
890
 
891
      return (E in IEEE_Emin .. IEEE_Emax) or else
892
         ((E = IEEE_Emin - 1) and then abs To_Float (SR) = 1.0);
893
   end Valid;
894
 
895
   ---------------------
896
   -- Unaligned_Valid --
897
   ---------------------
898
 
899
   function Unaligned_Valid (A : System.Address) return Boolean is
900
      subtype FS is String (1 .. T'Size / Character'Size);
901
      type FSP is access FS;
902
 
903
      function To_FSP is new Ada.Unchecked_Conversion (Address, FSP);
904
 
905
      Local_T : aliased T;
906
 
907
   begin
908
      --  Note that we have to be sure that we do not load the value into a
909
      --  floating-point register, since a signalling NaN may cause a trap.
910
      --  The following assignment is what does the actual alignment, since
911
      --  we know that the target Local_T is aligned.
912
 
913
      To_FSP (Local_T'Address).all := To_FSP (A).all;
914
 
915
      --  Now that we have an aligned value, we can use the normal aligned
916
      --  version of Valid to obtain the required result.
917
 
918
      return Valid (Local_T'Access);
919
   end Unaligned_Valid;
920
 
921
end System.Fat_Gen;

powered by: WebSVN 2.1.0

© copyright 1999-2024 OpenCores.org, equivalent to Oliscience, all rights reserved. OpenCores®, registered trademark.