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

Subversion Repositories openrisc

[/] [openrisc/] [trunk/] [gnu-dev/] [or1k-gcc/] [gcc/] [ada/] [eval_fat.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
--                             E V A L _ F A T                              --
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.  See the GNU General Public License --
17
-- for  more details.  You should have  received  a copy of the GNU General --
18
-- Public License  distributed with GNAT; see file COPYING3.  If not, go to --
19
-- http://www.gnu.org/licenses for a complete copy of the license.          --
20
--                                                                          --
21
-- GNAT was originally developed  by the GNAT team at  New York University. --
22
-- Extensive contributions were provided by Ada Core Technologies Inc.      --
23
--                                                                          --
24
------------------------------------------------------------------------------
25
 
26
with Einfo;    use Einfo;
27
with Errout;   use Errout;
28
with Targparm; use Targparm;
29
 
30
package body Eval_Fat is
31
 
32
   Radix : constant Int := 2;
33
   --  This code is currently only correct for the radix 2 case. We use the
34
   --  symbolic value Radix where possible to help in the unlikely case of
35
   --  anyone ever having to adjust this code for another value, and for
36
   --  documentation purposes.
37
 
38
   --  Another assumption is that the range of the floating-point type is
39
   --  symmetric around zero.
40
 
41
   type Radix_Power_Table is array (Int range 1 .. 4) of Int;
42
 
43
   Radix_Powers : constant Radix_Power_Table :=
44
                    (Radix ** 1, Radix ** 2, Radix ** 3, Radix ** 4);
45
 
46
   -----------------------
47
   -- Local Subprograms --
48
   -----------------------
49
 
50
   procedure Decompose
51
     (RT       : R;
52
      X        : T;
53
      Fraction : out T;
54
      Exponent : out UI;
55
      Mode     : Rounding_Mode := Round);
56
   --  Decomposes a non-zero floating-point number into fraction and exponent
57
   --  parts. The fraction is in the interval 1.0 / Radix .. T'Pred (1.0) and
58
   --  uses Rbase = Radix. The result is rounded to a nearest machine number.
59
 
60
   procedure Decompose_Int
61
     (RT       : R;
62
      X        : T;
63
      Fraction : out UI;
64
      Exponent : out UI;
65
      Mode     : Rounding_Mode);
66
   --  This is similar to Decompose, except that the Fraction value returned
67
   --  is an integer representing the value Fraction * Scale, where Scale is
68
   --  the value (Machine_Radix_Value (RT) ** Machine_Mantissa_Value (RT)). The
69
   --  value is obtained by using biased rounding (halfway cases round away
70
   --  from zero), round to even, a floor operation or a ceiling operation
71
   --  depending on the setting of Mode (see corresponding descriptions in
72
   --  Urealp).
73
 
74
   --------------
75
   -- Adjacent --
76
   --------------
77
 
78
   function Adjacent (RT : R; X, Towards : T) return T is
79
   begin
80
      if Towards = X then
81
         return X;
82
      elsif Towards > X then
83
         return Succ (RT, X);
84
      else
85
         return Pred (RT, X);
86
      end if;
87
   end Adjacent;
88
 
89
   -------------
90
   -- Ceiling --
91
   -------------
92
 
93
   function Ceiling (RT : R; X : T) return T is
94
      XT : constant T := Truncation (RT, X);
95
   begin
96
      if UR_Is_Negative (X) then
97
         return XT;
98
      elsif X = XT then
99
         return X;
100
      else
101
         return XT + Ureal_1;
102
      end if;
103
   end Ceiling;
104
 
105
   -------------
106
   -- Compose --
107
   -------------
108
 
109
   function Compose (RT : R; Fraction : T; Exponent : UI) return T is
110
      Arg_Frac : T;
111
      Arg_Exp  : UI;
112
      pragma Warnings (Off, Arg_Exp);
113
   begin
114
      Decompose (RT, Fraction, Arg_Frac, Arg_Exp);
115
      return Scaling (RT, Arg_Frac, Exponent);
116
   end Compose;
117
 
118
   ---------------
119
   -- Copy_Sign --
120
   ---------------
121
 
122
   function Copy_Sign (RT : R; Value, Sign : T) return T is
123
      pragma Warnings (Off, RT);
124
      Result : T;
125
 
126
   begin
127
      Result := abs Value;
128
 
129
      if UR_Is_Negative (Sign) then
130
         return -Result;
131
      else
132
         return Result;
133
      end if;
134
   end Copy_Sign;
135
 
136
   ---------------
137
   -- Decompose --
138
   ---------------
139
 
140
   procedure Decompose
141
     (RT       : R;
142
      X        : T;
143
      Fraction : out T;
144
      Exponent : out UI;
145
      Mode     : Rounding_Mode := Round)
146
   is
147
      Int_F : UI;
148
 
149
   begin
150
      Decompose_Int (RT, abs X, Int_F, Exponent, Mode);
151
 
152
      Fraction := UR_From_Components
153
       (Num      => Int_F,
154
        Den      => Machine_Mantissa_Value (RT),
155
        Rbase    => Radix,
156
        Negative => False);
157
 
158
      if UR_Is_Negative (X) then
159
         Fraction := -Fraction;
160
      end if;
161
 
162
      return;
163
   end Decompose;
164
 
165
   -------------------
166
   -- Decompose_Int --
167
   -------------------
168
 
169
   --  This procedure should be modified with care, as there are many non-
170
   --  obvious details that may cause problems that are hard to detect. For
171
   --  zero arguments, Fraction and Exponent are set to zero. Note that sign
172
   --  of zero cannot be preserved.
173
 
174
   procedure Decompose_Int
175
     (RT       : R;
176
      X        : T;
177
      Fraction : out UI;
178
      Exponent : out UI;
179
      Mode     : Rounding_Mode)
180
   is
181
      Base : Int := Rbase (X);
182
      N    : UI  := abs Numerator (X);
183
      D    : UI  := Denominator (X);
184
 
185
      N_Times_Radix : UI;
186
 
187
      Even : Boolean;
188
      --  True iff Fraction is even
189
 
190
      Most_Significant_Digit : constant UI :=
191
                                 Radix ** (Machine_Mantissa_Value (RT) - 1);
192
 
193
      Uintp_Mark : Uintp.Save_Mark;
194
      --  The code is divided into blocks that systematically release
195
      --  intermediate values (this routine generates lots of junk!)
196
 
197
   begin
198
      if N = Uint_0 then
199
         Fraction := Uint_0;
200
         Exponent := Uint_0;
201
         return;
202
      end if;
203
 
204
      Calculate_D_And_Exponent_1 : begin
205
         Uintp_Mark := Mark;
206
         Exponent := Uint_0;
207
 
208
         --  In cases where Base > 1, the actual denominator is Base**D. For
209
         --  cases where Base is a power of Radix, use the value 1 for the
210
         --  Denominator and adjust the exponent.
211
 
212
         --  Note: Exponent has different sign from D, because D is a divisor
213
 
214
         for Power in 1 .. Radix_Powers'Last loop
215
            if Base = Radix_Powers (Power) then
216
               Exponent := -D * Power;
217
               Base := 0;
218
               D := Uint_1;
219
               exit;
220
            end if;
221
         end loop;
222
 
223
         Release_And_Save (Uintp_Mark, D, Exponent);
224
      end Calculate_D_And_Exponent_1;
225
 
226
      if Base > 0 then
227
         Calculate_Exponent : begin
228
            Uintp_Mark := Mark;
229
 
230
            --  For bases that are a multiple of the Radix, divide the base by
231
            --  Radix and adjust the Exponent. This will help because D will be
232
            --  much smaller and faster to process.
233
 
234
            --  This occurs for decimal bases on machines with binary floating-
235
            --  point for example. When calculating 1E40, with Radix = 2, N
236
            --  will be 93 bits instead of 133.
237
 
238
            --        N            E
239
            --      ------  * Radix
240
            --           D
241
            --       Base
242
 
243
            --                  N                        E
244
            --    =  --------------------------  *  Radix
245
            --                     D        D
246
            --         (Base/Radix)  * Radix
247
 
248
            --             N                  E-D
249
            --    =  ---------------  *  Radix
250
            --                    D
251
            --        (Base/Radix)
252
 
253
            --  This code is commented out, because it causes numerous
254
            --  failures in the regression suite. To be studied ???
255
 
256
            while False and then Base > 0 and then Base mod Radix = 0 loop
257
               Base := Base / Radix;
258
               Exponent := Exponent + D;
259
            end loop;
260
 
261
            Release_And_Save (Uintp_Mark, Exponent);
262
         end Calculate_Exponent;
263
 
264
         --  For remaining bases we must actually compute the exponentiation
265
 
266
         --  Because the exponentiation can be negative, and D must be integer,
267
         --  the numerator is corrected instead.
268
 
269
         Calculate_N_And_D : begin
270
            Uintp_Mark := Mark;
271
 
272
            if D < 0 then
273
               N := N * Base ** (-D);
274
               D := Uint_1;
275
            else
276
               D := Base ** D;
277
            end if;
278
 
279
            Release_And_Save (Uintp_Mark, N, D);
280
         end Calculate_N_And_D;
281
 
282
         Base := 0;
283
      end if;
284
 
285
      --  Now scale N and D so that N / D is a value in the interval [1.0 /
286
      --  Radix, 1.0) and adjust Exponent accordingly, so the value N / D *
287
      --  Radix ** Exponent remains unchanged.
288
 
289
      --  Step 1 - Adjust N so N / D >= 1 / Radix, or N = 0
290
 
291
      --  N and D are positive, so N / D >= 1 / Radix implies N * Radix >= D.
292
      --  As this scaling is not possible for N is Uint_0, zero is handled
293
      --  explicitly at the start of this subprogram.
294
 
295
      Calculate_N_And_Exponent : begin
296
         Uintp_Mark := Mark;
297
 
298
         N_Times_Radix := N * Radix;
299
         while not (N_Times_Radix >= D) loop
300
            N := N_Times_Radix;
301
            Exponent := Exponent - 1;
302
            N_Times_Radix := N * Radix;
303
         end loop;
304
 
305
         Release_And_Save (Uintp_Mark, N, Exponent);
306
      end Calculate_N_And_Exponent;
307
 
308
      --  Step 2 - Adjust D so N / D < 1
309
 
310
      --  Scale up D so N / D < 1, so N < D
311
 
312
      Calculate_D_And_Exponent_2 : begin
313
         Uintp_Mark := Mark;
314
 
315
         while not (N < D) loop
316
 
317
            --  As N / D >= 1, N / (D * Radix) will be at least 1 / Radix, so
318
            --  the result of Step 1 stays valid
319
 
320
            D := D * Radix;
321
            Exponent := Exponent + 1;
322
         end loop;
323
 
324
         Release_And_Save (Uintp_Mark, D, Exponent);
325
      end Calculate_D_And_Exponent_2;
326
 
327
      --  Here the value N / D is in the range [1.0 / Radix .. 1.0)
328
 
329
      --  Now find the fraction by doing a very simple-minded division until
330
      --  enough digits have been computed.
331
 
332
      --  This division works for all radices, but is only efficient for a
333
      --  binary radix. It is just like a manual division algorithm, but
334
      --  instead of moving the denominator one digit right, we move the
335
      --  numerator one digit left so the numerator and denominator remain
336
      --  integral.
337
 
338
      Fraction := Uint_0;
339
      Even := True;
340
 
341
      Calculate_Fraction_And_N : begin
342
         Uintp_Mark := Mark;
343
 
344
         loop
345
            while N >= D loop
346
               N := N - D;
347
               Fraction := Fraction + 1;
348
               Even := not Even;
349
            end loop;
350
 
351
            --  Stop when the result is in [1.0 / Radix, 1.0)
352
 
353
            exit when Fraction >= Most_Significant_Digit;
354
 
355
            N := N * Radix;
356
            Fraction := Fraction * Radix;
357
            Even := True;
358
         end loop;
359
 
360
         Release_And_Save (Uintp_Mark, Fraction, N);
361
      end Calculate_Fraction_And_N;
362
 
363
      Calculate_Fraction_And_Exponent : begin
364
         Uintp_Mark := Mark;
365
 
366
         --  Determine correct rounding based on the remainder which is in
367
         --  N and the divisor D. The rounding is performed on the absolute
368
         --  value of X, so Ceiling and Floor need to check for the sign of
369
         --  X explicitly.
370
 
371
         case Mode is
372
            when Round_Even =>
373
 
374
               --  This rounding mode should not be used for static
375
               --  expressions, but only for compile-time evaluation of
376
               --  non-static expressions.
377
 
378
               if (Even and then N * 2 > D)
379
                     or else
380
                  (not Even and then N * 2 >= D)
381
               then
382
                  Fraction := Fraction + 1;
383
               end if;
384
 
385
            when Round   =>
386
 
387
               --  Do not round to even as is done with IEEE arithmetic, but
388
               --  instead round away from zero when the result is exactly
389
               --  between two machine numbers. See RM 4.9(38).
390
 
391
               if N * 2 >= D then
392
                  Fraction := Fraction + 1;
393
               end if;
394
 
395
            when Ceiling =>
396
               if N > Uint_0 and then not UR_Is_Negative (X) then
397
                  Fraction := Fraction + 1;
398
               end if;
399
 
400
            when Floor   =>
401
               if N > Uint_0 and then UR_Is_Negative (X) then
402
                  Fraction := Fraction + 1;
403
               end if;
404
         end case;
405
 
406
         --  The result must be normalized to [1.0/Radix, 1.0), so adjust if
407
         --  the result is 1.0 because of rounding.
408
 
409
         if Fraction = Most_Significant_Digit * Radix then
410
            Fraction := Most_Significant_Digit;
411
            Exponent := Exponent + 1;
412
         end if;
413
 
414
         --  Put back sign after applying the rounding
415
 
416
         if UR_Is_Negative (X) then
417
            Fraction := -Fraction;
418
         end if;
419
 
420
         Release_And_Save (Uintp_Mark, Fraction, Exponent);
421
      end Calculate_Fraction_And_Exponent;
422
   end Decompose_Int;
423
 
424
   --------------
425
   -- Exponent --
426
   --------------
427
 
428
   function Exponent (RT : R; X : T) return UI is
429
      X_Frac : UI;
430
      X_Exp  : UI;
431
      pragma Warnings (Off, X_Frac);
432
   begin
433
      Decompose_Int (RT, X, X_Frac, X_Exp, Round_Even);
434
      return X_Exp;
435
   end Exponent;
436
 
437
   -----------
438
   -- Floor --
439
   -----------
440
 
441
   function Floor (RT : R; X : T) return T is
442
      XT : constant T := Truncation (RT, X);
443
 
444
   begin
445
      if UR_Is_Positive (X) then
446
         return XT;
447
 
448
      elsif XT = X then
449
         return X;
450
 
451
      else
452
         return XT - Ureal_1;
453
      end if;
454
   end Floor;
455
 
456
   --------------
457
   -- Fraction --
458
   --------------
459
 
460
   function Fraction (RT : R; X : T) return T is
461
      X_Frac : T;
462
      X_Exp  : UI;
463
      pragma Warnings (Off, X_Exp);
464
   begin
465
      Decompose (RT, X, X_Frac, X_Exp);
466
      return X_Frac;
467
   end Fraction;
468
 
469
   ------------------
470
   -- Leading_Part --
471
   ------------------
472
 
473
   function Leading_Part (RT : R; X : T; Radix_Digits : UI) return T is
474
      RD : constant UI := UI_Min (Radix_Digits, Machine_Mantissa_Value (RT));
475
      L  : UI;
476
      Y  : T;
477
   begin
478
      L := Exponent (RT, X) - RD;
479
      Y := UR_From_Uint (UR_Trunc (Scaling (RT, X, -L)));
480
      return Scaling (RT, Y, L);
481
   end Leading_Part;
482
 
483
   -------------
484
   -- Machine --
485
   -------------
486
 
487
   function Machine
488
     (RT    : R;
489
      X     : T;
490
      Mode  : Rounding_Mode;
491
      Enode : Node_Id) return T
492
   is
493
      X_Frac : T;
494
      X_Exp  : UI;
495
      Emin   : constant UI := Machine_Emin_Value (RT);
496
 
497
   begin
498
      Decompose (RT, X, X_Frac, X_Exp, Mode);
499
 
500
      --  Case of denormalized number or (gradual) underflow
501
 
502
      --  A denormalized number is one with the minimum exponent Emin, but that
503
      --  breaks the assumption that the first digit of the mantissa is a one.
504
      --  This allows the first non-zero digit to be in any of the remaining
505
      --  Mant - 1 spots. The gap between subsequent denormalized numbers is
506
      --  the same as for the smallest normalized numbers. However, the number
507
      --  of significant digits left decreases as a result of the mantissa now
508
      --  having leading seros.
509
 
510
      if X_Exp < Emin then
511
         declare
512
            Emin_Den : constant UI := Machine_Emin_Value (RT)
513
                                        - Machine_Mantissa_Value (RT) + Uint_1;
514
         begin
515
            if X_Exp < Emin_Den or not Denorm_On_Target then
516
               if UR_Is_Negative (X) then
517
                  Error_Msg_N
518
                    ("floating-point value underflows to -0.0?", Enode);
519
                  return Ureal_M_0;
520
 
521
               else
522
                  Error_Msg_N
523
                    ("floating-point value underflows to 0.0?", Enode);
524
                  return Ureal_0;
525
               end if;
526
 
527
            elsif Denorm_On_Target then
528
 
529
               --  Emin - Mant <= X_Exp < Emin, so result is denormal. Handle
530
               --  gradual underflow by first computing the number of
531
               --  significant bits still available for the mantissa and
532
               --  then truncating the fraction to this number of bits.
533
 
534
               --  If this value is different from the original fraction,
535
               --  precision is lost due to gradual underflow.
536
 
537
               --  We probably should round here and prevent double rounding as
538
               --  a result of first rounding to a model number and then to a
539
               --  machine number. However, this is an extremely rare case that
540
               --  is not worth the extra complexity. In any case, a warning is
541
               --  issued in cases where gradual underflow occurs.
542
 
543
               declare
544
                  Denorm_Sig_Bits : constant UI := X_Exp - Emin_Den + 1;
545
 
546
                  X_Frac_Denorm   : constant T := UR_From_Components
547
                    (UR_Trunc (Scaling (RT, abs X_Frac, Denorm_Sig_Bits)),
548
                     Denorm_Sig_Bits,
549
                     Radix,
550
                     UR_Is_Negative (X));
551
 
552
               begin
553
                  if X_Frac_Denorm /= X_Frac then
554
                     Error_Msg_N
555
                       ("gradual underflow causes loss of precision?",
556
                        Enode);
557
                     X_Frac := X_Frac_Denorm;
558
                  end if;
559
               end;
560
            end if;
561
         end;
562
      end if;
563
 
564
      return Scaling (RT, X_Frac, X_Exp);
565
   end Machine;
566
 
567
   -----------
568
   -- Model --
569
   -----------
570
 
571
   function Model (RT : R; X : T) return T is
572
      X_Frac : T;
573
      X_Exp  : UI;
574
   begin
575
      Decompose (RT, X, X_Frac, X_Exp);
576
      return Compose (RT, X_Frac, X_Exp);
577
   end Model;
578
 
579
   ----------
580
   -- Pred --
581
   ----------
582
 
583
   function Pred (RT : R; X : T) return T is
584
   begin
585
      return -Succ (RT, -X);
586
   end Pred;
587
 
588
   ---------------
589
   -- Remainder --
590
   ---------------
591
 
592
   function Remainder (RT : R; X, Y : T) return T is
593
      A        : T;
594
      B        : T;
595
      Arg      : T;
596
      P        : T;
597
      Arg_Frac : T;
598
      P_Frac   : T;
599
      Sign_X   : T;
600
      IEEE_Rem : T;
601
      Arg_Exp  : UI;
602
      P_Exp    : UI;
603
      K        : UI;
604
      P_Even   : Boolean;
605
 
606
      pragma Warnings (Off, Arg_Frac);
607
 
608
   begin
609
      if UR_Is_Positive (X) then
610
         Sign_X :=  Ureal_1;
611
      else
612
         Sign_X := -Ureal_1;
613
      end if;
614
 
615
      Arg := abs X;
616
      P   := abs Y;
617
 
618
      if Arg < P then
619
         P_Even := True;
620
         IEEE_Rem := Arg;
621
         P_Exp := Exponent (RT, P);
622
 
623
      else
624
         --  ??? what about zero cases?
625
         Decompose (RT, Arg, Arg_Frac, Arg_Exp);
626
         Decompose (RT, P,   P_Frac,   P_Exp);
627
 
628
         P := Compose (RT, P_Frac, Arg_Exp);
629
         K := Arg_Exp - P_Exp;
630
         P_Even := True;
631
         IEEE_Rem := Arg;
632
 
633
         for Cnt in reverse 0 .. UI_To_Int (K) loop
634
            if IEEE_Rem >= P then
635
               P_Even := False;
636
               IEEE_Rem := IEEE_Rem - P;
637
            else
638
               P_Even := True;
639
            end if;
640
 
641
            P := P * Ureal_Half;
642
         end loop;
643
      end if;
644
 
645
      --  That completes the calculation of modulus remainder. The final step
646
      --  is get the IEEE remainder. Here we compare Rem with (abs Y) / 2.
647
 
648
      if P_Exp >= 0 then
649
         A := IEEE_Rem;
650
         B := abs Y * Ureal_Half;
651
 
652
      else
653
         A := IEEE_Rem * Ureal_2;
654
         B := abs Y;
655
      end if;
656
 
657
      if A > B or else (A = B and then not P_Even) then
658
         IEEE_Rem := IEEE_Rem - abs Y;
659
      end if;
660
 
661
      return Sign_X * IEEE_Rem;
662
   end Remainder;
663
 
664
   --------------
665
   -- Rounding --
666
   --------------
667
 
668
   function Rounding (RT : R; X : T) return T is
669
      Result : T;
670
      Tail   : T;
671
 
672
   begin
673
      Result := Truncation (RT, abs X);
674
      Tail   := abs X - Result;
675
 
676
      if Tail >= Ureal_Half  then
677
         Result := Result + Ureal_1;
678
      end if;
679
 
680
      if UR_Is_Negative (X) then
681
         return -Result;
682
      else
683
         return Result;
684
      end if;
685
   end Rounding;
686
 
687
   -------------
688
   -- Scaling --
689
   -------------
690
 
691
   function Scaling (RT : R; X : T; Adjustment : UI) return T is
692
      pragma Warnings (Off, RT);
693
 
694
   begin
695
      if Rbase (X) = Radix then
696
         return UR_From_Components
697
           (Num      => Numerator (X),
698
            Den      => Denominator (X) - Adjustment,
699
            Rbase    => Radix,
700
            Negative => UR_Is_Negative (X));
701
 
702
      elsif Adjustment >= 0 then
703
         return X * Radix ** Adjustment;
704
      else
705
         return X / Radix ** (-Adjustment);
706
      end if;
707
   end Scaling;
708
 
709
   ----------
710
   -- Succ --
711
   ----------
712
 
713
   function Succ (RT : R; X : T) return T is
714
      Emin     : constant UI := Machine_Emin_Value (RT);
715
      Mantissa : constant UI := Machine_Mantissa_Value (RT);
716
      Exp      : UI := UI_Max (Emin, Exponent (RT, X));
717
      Frac     : T;
718
      New_Frac : T;
719
 
720
   begin
721
      if UR_Is_Zero (X) then
722
         Exp := Emin;
723
      end if;
724
 
725
      --  Set exponent such that the radix point will be directly following the
726
      --  mantissa after scaling.
727
 
728
      if Denorm_On_Target or Exp /= Emin then
729
         Exp := Exp - Mantissa;
730
      else
731
         Exp := Exp - 1;
732
      end if;
733
 
734
      Frac := Scaling (RT, X, -Exp);
735
      New_Frac := Ceiling (RT, Frac);
736
 
737
      if New_Frac = Frac then
738
         if New_Frac = Scaling (RT, -Ureal_1, Mantissa - 1) then
739
            New_Frac := New_Frac + Scaling (RT, Ureal_1, Uint_Minus_1);
740
         else
741
            New_Frac := New_Frac + Ureal_1;
742
         end if;
743
      end if;
744
 
745
      return Scaling (RT, New_Frac, Exp);
746
   end Succ;
747
 
748
   ----------------
749
   -- Truncation --
750
   ----------------
751
 
752
   function Truncation (RT : R; X : T) return T is
753
      pragma Warnings (Off, RT);
754
   begin
755
      return UR_From_Uint (UR_Trunc (X));
756
   end Truncation;
757
 
758
   -----------------------
759
   -- Unbiased_Rounding --
760
   -----------------------
761
 
762
   function Unbiased_Rounding (RT : R; X : T) return T is
763
      Abs_X  : constant T := abs X;
764
      Result : T;
765
      Tail   : T;
766
 
767
   begin
768
      Result := Truncation (RT, Abs_X);
769
      Tail   := Abs_X - Result;
770
 
771
      if Tail > Ureal_Half  then
772
         Result := Result + Ureal_1;
773
 
774
      elsif Tail = Ureal_Half then
775
         Result := Ureal_2 *
776
                     Truncation (RT, (Result / Ureal_2) + Ureal_Half);
777
      end if;
778
 
779
      if UR_Is_Negative (X) then
780
         return -Result;
781
      elsif UR_Is_Positive (X) then
782
         return Result;
783
 
784
      --  For zero case, make sure sign of zero is preserved
785
 
786
      else
787
         return X;
788
      end if;
789
   end Unbiased_Rounding;
790
 
791
end Eval_Fat;

powered by: WebSVN 2.1.0

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