Daily bump.
[official-gcc.git] / gcc / ada / eval_fat.adb
blob78dcea60165dd3f35dee5e391e1b85f57bc05082
1 ------------------------------------------------------------------------------
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-2007, 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 ------------------------------------------------------------------------------
26 with Einfo; use Einfo;
27 with Errout; use Errout;
28 with Sem_Util; use Sem_Util;
29 with Ttypef; use Ttypef;
30 with Targparm; use Targparm;
32 package body Eval_Fat is
34 Radix : constant Int := 2;
35 -- This code is currently only correct for the radix 2 case. We use the
36 -- symbolic value Radix where possible to help in the unlikely case of
37 -- anyone ever having to adjust this code for another value, and for
38 -- documentation purposes.
40 -- Another assumption is that the range of the floating-point type is
41 -- symmetric around zero.
43 type Radix_Power_Table is array (Int range 1 .. 4) of Int;
45 Radix_Powers : constant Radix_Power_Table :=
46 (Radix ** 1, Radix ** 2, Radix ** 3, Radix ** 4);
48 -----------------------
49 -- Local Subprograms --
50 -----------------------
52 procedure Decompose
53 (RT : R;
54 X : T;
55 Fraction : out T;
56 Exponent : out UI;
57 Mode : Rounding_Mode := Round);
58 -- Decomposes a non-zero floating-point number into fraction and exponent
59 -- parts. The fraction is in the interval 1.0 / Radix .. T'Pred (1.0) and
60 -- uses Rbase = Radix. The result is rounded to a nearest machine number.
62 procedure Decompose_Int
63 (RT : R;
64 X : T;
65 Fraction : out UI;
66 Exponent : out UI;
67 Mode : Rounding_Mode);
68 -- This is similar to Decompose, except that the Fraction value returned
69 -- is an integer representing the value Fraction * Scale, where Scale is
70 -- the value (Radix ** Machine_Mantissa (RT)). The value is obtained by
71 -- using biased rounding (halfway cases round away from zero), round to
72 -- even, a floor operation or a ceiling operation depending on the setting
73 -- of Mode (see corresponding descriptions in Urealp).
75 function Machine_Emin (RT : R) return Int;
76 -- Return value of the Machine_Emin attribute
78 --------------
79 -- Adjacent --
80 --------------
82 function Adjacent (RT : R; X, Towards : T) return T is
83 begin
84 if Towards = X then
85 return X;
86 elsif Towards > X then
87 return Succ (RT, X);
88 else
89 return Pred (RT, X);
90 end if;
91 end Adjacent;
93 -------------
94 -- Ceiling --
95 -------------
97 function Ceiling (RT : R; X : T) return T is
98 XT : constant T := Truncation (RT, X);
99 begin
100 if UR_Is_Negative (X) then
101 return XT;
102 elsif X = XT then
103 return X;
104 else
105 return XT + Ureal_1;
106 end if;
107 end Ceiling;
109 -------------
110 -- Compose --
111 -------------
113 function Compose (RT : R; Fraction : T; Exponent : UI) return T is
114 Arg_Frac : T;
115 Arg_Exp : UI;
116 pragma Warnings (Off, Arg_Exp);
117 begin
118 Decompose (RT, Fraction, Arg_Frac, Arg_Exp);
119 return Scaling (RT, Arg_Frac, Exponent);
120 end Compose;
122 ---------------
123 -- Copy_Sign --
124 ---------------
126 function Copy_Sign (RT : R; Value, Sign : T) return T is
127 pragma Warnings (Off, RT);
128 Result : T;
130 begin
131 Result := abs Value;
133 if UR_Is_Negative (Sign) then
134 return -Result;
135 else
136 return Result;
137 end if;
138 end Copy_Sign;
140 ---------------
141 -- Decompose --
142 ---------------
144 procedure Decompose
145 (RT : R;
146 X : T;
147 Fraction : out T;
148 Exponent : out UI;
149 Mode : Rounding_Mode := Round)
151 Int_F : UI;
153 begin
154 Decompose_Int (RT, abs X, Int_F, Exponent, Mode);
156 Fraction := UR_From_Components
157 (Num => Int_F,
158 Den => UI_From_Int (Machine_Mantissa (RT)),
159 Rbase => Radix,
160 Negative => False);
162 if UR_Is_Negative (X) then
163 Fraction := -Fraction;
164 end if;
166 return;
167 end Decompose;
169 -------------------
170 -- Decompose_Int --
171 -------------------
173 -- This procedure should be modified with care, as there are many non-
174 -- obvious details that may cause problems that are hard to detect. For
175 -- zero arguments, Fraction and Exponent are set to zero. Note that sign
176 -- of zero cannot be preserved.
178 procedure Decompose_Int
179 (RT : R;
180 X : T;
181 Fraction : out UI;
182 Exponent : out UI;
183 Mode : Rounding_Mode)
185 Base : Int := Rbase (X);
186 N : UI := abs Numerator (X);
187 D : UI := Denominator (X);
189 N_Times_Radix : UI;
191 Even : Boolean;
192 -- True iff Fraction is even
194 Most_Significant_Digit : constant UI :=
195 Radix ** (Machine_Mantissa (RT) - 1);
197 Uintp_Mark : Uintp.Save_Mark;
198 -- The code is divided into blocks that systematically release
199 -- intermediate values (this routine generates lots of junk!)
201 begin
202 if N = Uint_0 then
203 Fraction := Uint_0;
204 Exponent := Uint_0;
205 return;
206 end if;
208 Calculate_D_And_Exponent_1 : begin
209 Uintp_Mark := Mark;
210 Exponent := Uint_0;
212 -- In cases where Base > 1, the actual denominator is Base**D. For
213 -- cases where Base is a power of Radix, use the value 1 for the
214 -- Denominator and adjust the exponent.
216 -- Note: Exponent has different sign from D, because D is a divisor
218 for Power in 1 .. Radix_Powers'Last loop
219 if Base = Radix_Powers (Power) then
220 Exponent := -D * Power;
221 Base := 0;
222 D := Uint_1;
223 exit;
224 end if;
225 end loop;
227 Release_And_Save (Uintp_Mark, D, Exponent);
228 end Calculate_D_And_Exponent_1;
230 if Base > 0 then
231 Calculate_Exponent : begin
232 Uintp_Mark := Mark;
234 -- For bases that are a multiple of the Radix, divide the base by
235 -- Radix and adjust the Exponent. This will help because D will be
236 -- much smaller and faster to process.
238 -- This occurs for decimal bases on machines with binary floating-
239 -- point for example. When calculating 1E40, with Radix = 2, N
240 -- will be 93 bits instead of 133.
242 -- N E
243 -- ------ * Radix
244 -- D
245 -- Base
247 -- N E
248 -- = -------------------------- * Radix
249 -- D D
250 -- (Base/Radix) * Radix
252 -- N E-D
253 -- = --------------- * Radix
254 -- D
255 -- (Base/Radix)
257 -- This code is commented out, because it causes numerous
258 -- failures in the regression suite. To be studied ???
260 while False and then Base > 0 and then Base mod Radix = 0 loop
261 Base := Base / Radix;
262 Exponent := Exponent + D;
263 end loop;
265 Release_And_Save (Uintp_Mark, Exponent);
266 end Calculate_Exponent;
268 -- For remaining bases we must actually compute the exponentiation
270 -- Because the exponentiation can be negative, and D must be integer,
271 -- the numerator is corrected instead.
273 Calculate_N_And_D : begin
274 Uintp_Mark := Mark;
276 if D < 0 then
277 N := N * Base ** (-D);
278 D := Uint_1;
279 else
280 D := Base ** D;
281 end if;
283 Release_And_Save (Uintp_Mark, N, D);
284 end Calculate_N_And_D;
286 Base := 0;
287 end if;
289 -- Now scale N and D so that N / D is a value in the interval [1.0 /
290 -- Radix, 1.0) and adjust Exponent accordingly, so the value N / D *
291 -- Radix ** Exponent remains unchanged.
293 -- Step 1 - Adjust N so N / D >= 1 / Radix, or N = 0
295 -- N and D are positive, so N / D >= 1 / Radix implies N * Radix >= D.
296 -- As this scaling is not possible for N is Uint_0, zero is handled
297 -- explicitly at the start of this subprogram.
299 Calculate_N_And_Exponent : begin
300 Uintp_Mark := Mark;
302 N_Times_Radix := N * Radix;
303 while not (N_Times_Radix >= D) loop
304 N := N_Times_Radix;
305 Exponent := Exponent - 1;
306 N_Times_Radix := N * Radix;
307 end loop;
309 Release_And_Save (Uintp_Mark, N, Exponent);
310 end Calculate_N_And_Exponent;
312 -- Step 2 - Adjust D so N / D < 1
314 -- Scale up D so N / D < 1, so N < D
316 Calculate_D_And_Exponent_2 : begin
317 Uintp_Mark := Mark;
319 while not (N < D) loop
321 -- As N / D >= 1, N / (D * Radix) will be at least 1 / Radix, so
322 -- the result of Step 1 stays valid
324 D := D * Radix;
325 Exponent := Exponent + 1;
326 end loop;
328 Release_And_Save (Uintp_Mark, D, Exponent);
329 end Calculate_D_And_Exponent_2;
331 -- Here the value N / D is in the range [1.0 / Radix .. 1.0)
333 -- Now find the fraction by doing a very simple-minded division until
334 -- enough digits have been computed.
336 -- This division works for all radices, but is only efficient for a
337 -- binary radix. It is just like a manual division algorithm, but
338 -- instead of moving the denominator one digit right, we move the
339 -- numerator one digit left so the numerator and denominator remain
340 -- integral.
342 Fraction := Uint_0;
343 Even := True;
345 Calculate_Fraction_And_N : begin
346 Uintp_Mark := Mark;
348 loop
349 while N >= D loop
350 N := N - D;
351 Fraction := Fraction + 1;
352 Even := not Even;
353 end loop;
355 -- Stop when the result is in [1.0 / Radix, 1.0)
357 exit when Fraction >= Most_Significant_Digit;
359 N := N * Radix;
360 Fraction := Fraction * Radix;
361 Even := True;
362 end loop;
364 Release_And_Save (Uintp_Mark, Fraction, N);
365 end Calculate_Fraction_And_N;
367 Calculate_Fraction_And_Exponent : begin
368 Uintp_Mark := Mark;
370 -- Determine correct rounding based on the remainder which is in
371 -- N and the divisor D. The rounding is performed on the absolute
372 -- value of X, so Ceiling and Floor need to check for the sign of
373 -- X explicitly.
375 case Mode is
376 when Round_Even =>
378 -- This rounding mode should not be used for static
379 -- expressions, but only for compile-time evaluation of
380 -- non-static expressions.
382 if (Even and then N * 2 > D)
383 or else
384 (not Even and then N * 2 >= D)
385 then
386 Fraction := Fraction + 1;
387 end if;
389 when Round =>
391 -- Do not round to even as is done with IEEE arithmetic, but
392 -- instead round away from zero when the result is exactly
393 -- between two machine numbers. See RM 4.9(38).
395 if N * 2 >= D then
396 Fraction := Fraction + 1;
397 end if;
399 when Ceiling =>
400 if N > Uint_0 and then not UR_Is_Negative (X) then
401 Fraction := Fraction + 1;
402 end if;
404 when Floor =>
405 if N > Uint_0 and then UR_Is_Negative (X) then
406 Fraction := Fraction + 1;
407 end if;
408 end case;
410 -- The result must be normalized to [1.0/Radix, 1.0), so adjust if
411 -- the result is 1.0 because of rounding.
413 if Fraction = Most_Significant_Digit * Radix then
414 Fraction := Most_Significant_Digit;
415 Exponent := Exponent + 1;
416 end if;
418 -- Put back sign after applying the rounding
420 if UR_Is_Negative (X) then
421 Fraction := -Fraction;
422 end if;
424 Release_And_Save (Uintp_Mark, Fraction, Exponent);
425 end Calculate_Fraction_And_Exponent;
426 end Decompose_Int;
428 --------------
429 -- Exponent --
430 --------------
432 function Exponent (RT : R; X : T) return UI is
433 X_Frac : UI;
434 X_Exp : UI;
435 pragma Warnings (Off, X_Frac);
436 begin
437 Decompose_Int (RT, X, X_Frac, X_Exp, Round_Even);
438 return X_Exp;
439 end Exponent;
441 -----------
442 -- Floor --
443 -----------
445 function Floor (RT : R; X : T) return T is
446 XT : constant T := Truncation (RT, X);
448 begin
449 if UR_Is_Positive (X) then
450 return XT;
452 elsif XT = X then
453 return X;
455 else
456 return XT - Ureal_1;
457 end if;
458 end Floor;
460 --------------
461 -- Fraction --
462 --------------
464 function Fraction (RT : R; X : T) return T is
465 X_Frac : T;
466 X_Exp : UI;
467 pragma Warnings (Off, X_Exp);
468 begin
469 Decompose (RT, X, X_Frac, X_Exp);
470 return X_Frac;
471 end Fraction;
473 ------------------
474 -- Leading_Part --
475 ------------------
477 function Leading_Part (RT : R; X : T; Radix_Digits : UI) return T is
478 RD : constant UI := UI_Min (Radix_Digits, Machine_Mantissa (RT));
479 L : UI;
480 Y : T;
481 begin
482 L := Exponent (RT, X) - RD;
483 Y := UR_From_Uint (UR_Trunc (Scaling (RT, X, -L)));
484 return Scaling (RT, Y, L);
485 end Leading_Part;
487 -------------
488 -- Machine --
489 -------------
491 function Machine
492 (RT : R;
493 X : T;
494 Mode : Rounding_Mode;
495 Enode : Node_Id) return T
497 X_Frac : T;
498 X_Exp : UI;
499 Emin : constant UI := UI_From_Int (Machine_Emin (RT));
501 begin
502 Decompose (RT, X, X_Frac, X_Exp, Mode);
504 -- Case of denormalized number or (gradual) underflow
506 -- A denormalized number is one with the minimum exponent Emin, but that
507 -- breaks the assumption that the first digit of the mantissa is a one.
508 -- This allows the first non-zero digit to be in any of the remaining
509 -- Mant - 1 spots. The gap between subsequent denormalized numbers is
510 -- the same as for the smallest normalized numbers. However, the number
511 -- of significant digits left decreases as a result of the mantissa now
512 -- having leading seros.
514 if X_Exp < Emin then
515 declare
516 Emin_Den : constant UI :=
517 UI_From_Int
518 (Machine_Emin (RT) - Machine_Mantissa (RT) + 1);
519 begin
520 if X_Exp < Emin_Den or not Denorm_On_Target then
521 if UR_Is_Negative (X) then
522 Error_Msg_N
523 ("floating-point value underflows to -0.0?", Enode);
524 return Ureal_M_0;
526 else
527 Error_Msg_N
528 ("floating-point value underflows to 0.0?", Enode);
529 return Ureal_0;
530 end if;
532 elsif Denorm_On_Target then
534 -- Emin - Mant <= X_Exp < Emin, so result is denormal. Handle
535 -- gradual underflow by first computing the number of
536 -- significant bits still available for the mantissa and
537 -- then truncating the fraction to this number of bits.
539 -- If this value is different from the original fraction,
540 -- precision is lost due to gradual underflow.
542 -- We probably should round here and prevent double rounding as
543 -- a result of first rounding to a model number and then to a
544 -- machine number. However, this is an extremely rare case that
545 -- is not worth the extra complexity. In any case, a warning is
546 -- issued in cases where gradual underflow occurs.
548 declare
549 Denorm_Sig_Bits : constant UI := X_Exp - Emin_Den + 1;
551 X_Frac_Denorm : constant T := UR_From_Components
552 (UR_Trunc (Scaling (RT, abs X_Frac, Denorm_Sig_Bits)),
553 Denorm_Sig_Bits,
554 Radix,
555 UR_Is_Negative (X));
557 begin
558 if X_Frac_Denorm /= X_Frac then
559 Error_Msg_N
560 ("gradual underflow causes loss of precision?",
561 Enode);
562 X_Frac := X_Frac_Denorm;
563 end if;
564 end;
565 end if;
566 end;
567 end if;
569 return Scaling (RT, X_Frac, X_Exp);
570 end Machine;
572 ------------------
573 -- Machine_Emin --
574 ------------------
576 function Machine_Emin (RT : R) return Int is
577 Digs : constant UI := Digits_Value (RT);
578 Emin : Int;
580 begin
581 if Vax_Float (RT) then
582 if Digs = VAXFF_Digits then
583 Emin := VAXFF_Machine_Emin;
585 elsif Digs = VAXDF_Digits then
586 Emin := VAXDF_Machine_Emin;
588 else
589 pragma Assert (Digs = VAXGF_Digits);
590 Emin := VAXGF_Machine_Emin;
591 end if;
593 elsif Is_AAMP_Float (RT) then
594 if Digs = AAMPS_Digits then
595 Emin := AAMPS_Machine_Emin;
597 else
598 pragma Assert (Digs = AAMPL_Digits);
599 Emin := AAMPL_Machine_Emin;
600 end if;
602 else
603 if Digs = IEEES_Digits then
604 Emin := IEEES_Machine_Emin;
606 elsif Digs = IEEEL_Digits then
607 Emin := IEEEL_Machine_Emin;
609 else
610 pragma Assert (Digs = IEEEX_Digits);
611 Emin := IEEEX_Machine_Emin;
612 end if;
613 end if;
615 return Emin;
616 end Machine_Emin;
618 ----------------------
619 -- Machine_Mantissa --
620 ----------------------
622 function Machine_Mantissa (RT : R) return Nat is
623 Digs : constant UI := Digits_Value (RT);
624 Mant : Nat;
626 begin
627 if Vax_Float (RT) then
628 if Digs = VAXFF_Digits then
629 Mant := VAXFF_Machine_Mantissa;
631 elsif Digs = VAXDF_Digits then
632 Mant := VAXDF_Machine_Mantissa;
634 else
635 pragma Assert (Digs = VAXGF_Digits);
636 Mant := VAXGF_Machine_Mantissa;
637 end if;
639 elsif Is_AAMP_Float (RT) then
640 if Digs = AAMPS_Digits then
641 Mant := AAMPS_Machine_Mantissa;
643 else
644 pragma Assert (Digs = AAMPL_Digits);
645 Mant := AAMPL_Machine_Mantissa;
646 end if;
648 else
649 if Digs = IEEES_Digits then
650 Mant := IEEES_Machine_Mantissa;
652 elsif Digs = IEEEL_Digits then
653 Mant := IEEEL_Machine_Mantissa;
655 else
656 pragma Assert (Digs = IEEEX_Digits);
657 Mant := IEEEX_Machine_Mantissa;
658 end if;
659 end if;
661 return Mant;
662 end Machine_Mantissa;
664 -------------------
665 -- Machine_Radix --
666 -------------------
668 function Machine_Radix (RT : R) return Nat is
669 pragma Warnings (Off, RT);
670 begin
671 return Radix;
672 end Machine_Radix;
674 -----------
675 -- Model --
676 -----------
678 function Model (RT : R; X : T) return T is
679 X_Frac : T;
680 X_Exp : UI;
681 begin
682 Decompose (RT, X, X_Frac, X_Exp);
683 return Compose (RT, X_Frac, X_Exp);
684 end Model;
686 ----------
687 -- Pred --
688 ----------
690 function Pred (RT : R; X : T) return T is
691 begin
692 return -Succ (RT, -X);
693 end Pred;
695 ---------------
696 -- Remainder --
697 ---------------
699 function Remainder (RT : R; X, Y : T) return T is
700 A : T;
701 B : T;
702 Arg : T;
703 P : T;
704 Arg_Frac : T;
705 P_Frac : T;
706 Sign_X : T;
707 IEEE_Rem : T;
708 Arg_Exp : UI;
709 P_Exp : UI;
710 K : UI;
711 P_Even : Boolean;
713 pragma Warnings (Off, Arg_Frac);
715 begin
716 if UR_Is_Positive (X) then
717 Sign_X := Ureal_1;
718 else
719 Sign_X := -Ureal_1;
720 end if;
722 Arg := abs X;
723 P := abs Y;
725 if Arg < P then
726 P_Even := True;
727 IEEE_Rem := Arg;
728 P_Exp := Exponent (RT, P);
730 else
731 -- ??? what about zero cases?
732 Decompose (RT, Arg, Arg_Frac, Arg_Exp);
733 Decompose (RT, P, P_Frac, P_Exp);
735 P := Compose (RT, P_Frac, Arg_Exp);
736 K := Arg_Exp - P_Exp;
737 P_Even := True;
738 IEEE_Rem := Arg;
740 for Cnt in reverse 0 .. UI_To_Int (K) loop
741 if IEEE_Rem >= P then
742 P_Even := False;
743 IEEE_Rem := IEEE_Rem - P;
744 else
745 P_Even := True;
746 end if;
748 P := P * Ureal_Half;
749 end loop;
750 end if;
752 -- That completes the calculation of modulus remainder. The final step
753 -- is get the IEEE remainder. Here we compare Rem with (abs Y) / 2.
755 if P_Exp >= 0 then
756 A := IEEE_Rem;
757 B := abs Y * Ureal_Half;
759 else
760 A := IEEE_Rem * Ureal_2;
761 B := abs Y;
762 end if;
764 if A > B or else (A = B and then not P_Even) then
765 IEEE_Rem := IEEE_Rem - abs Y;
766 end if;
768 return Sign_X * IEEE_Rem;
769 end Remainder;
771 --------------
772 -- Rounding --
773 --------------
775 function Rounding (RT : R; X : T) return T is
776 Result : T;
777 Tail : T;
779 begin
780 Result := Truncation (RT, abs X);
781 Tail := abs X - Result;
783 if Tail >= Ureal_Half then
784 Result := Result + Ureal_1;
785 end if;
787 if UR_Is_Negative (X) then
788 return -Result;
789 else
790 return Result;
791 end if;
792 end Rounding;
794 -------------
795 -- Scaling --
796 -------------
798 function Scaling (RT : R; X : T; Adjustment : UI) return T is
799 pragma Warnings (Off, RT);
801 begin
802 if Rbase (X) = Radix then
803 return UR_From_Components
804 (Num => Numerator (X),
805 Den => Denominator (X) - Adjustment,
806 Rbase => Radix,
807 Negative => UR_Is_Negative (X));
809 elsif Adjustment >= 0 then
810 return X * Radix ** Adjustment;
811 else
812 return X / Radix ** (-Adjustment);
813 end if;
814 end Scaling;
816 ----------
817 -- Succ --
818 ----------
820 function Succ (RT : R; X : T) return T is
821 Emin : constant UI := UI_From_Int (Machine_Emin (RT));
822 Mantissa : constant UI := UI_From_Int (Machine_Mantissa (RT));
823 Exp : UI := UI_Max (Emin, Exponent (RT, X));
824 Frac : T;
825 New_Frac : T;
827 begin
828 if UR_Is_Zero (X) then
829 Exp := Emin;
830 end if;
832 -- Set exponent such that the radix point will be directly following the
833 -- mantissa after scaling.
835 if Denorm_On_Target or Exp /= Emin then
836 Exp := Exp - Mantissa;
837 else
838 Exp := Exp - 1;
839 end if;
841 Frac := Scaling (RT, X, -Exp);
842 New_Frac := Ceiling (RT, Frac);
844 if New_Frac = Frac then
845 if New_Frac = Scaling (RT, -Ureal_1, Mantissa - 1) then
846 New_Frac := New_Frac + Scaling (RT, Ureal_1, Uint_Minus_1);
847 else
848 New_Frac := New_Frac + Ureal_1;
849 end if;
850 end if;
852 return Scaling (RT, New_Frac, Exp);
853 end Succ;
855 ----------------
856 -- Truncation --
857 ----------------
859 function Truncation (RT : R; X : T) return T is
860 pragma Warnings (Off, RT);
861 begin
862 return UR_From_Uint (UR_Trunc (X));
863 end Truncation;
865 -----------------------
866 -- Unbiased_Rounding --
867 -----------------------
869 function Unbiased_Rounding (RT : R; X : T) return T is
870 Abs_X : constant T := abs X;
871 Result : T;
872 Tail : T;
874 begin
875 Result := Truncation (RT, Abs_X);
876 Tail := Abs_X - Result;
878 if Tail > Ureal_Half then
879 Result := Result + Ureal_1;
881 elsif Tail = Ureal_Half then
882 Result := Ureal_2 *
883 Truncation (RT, (Result / Ureal_2) + Ureal_Half);
884 end if;
886 if UR_Is_Negative (X) then
887 return -Result;
888 elsif UR_Is_Positive (X) then
889 return Result;
891 -- For zero case, make sure sign of zero is preserved
893 else
894 return X;
895 end if;
896 end Unbiased_Rounding;
898 end Eval_Fat;