]> www.ginac.de Git - cln.git/blob - src/float/ffloat/conv/cl_RA_to_float.cc
877f82489772a4da61dd809652667095fbf02de2
[cln.git] / src / float / ffloat / conv / cl_RA_to_float.cc
1 // cl_float_approx().
2
3 // General includes.
4 #include "cl_sysdep.h"
5
6 // Specification.
7 #include "cl_rational.h"
8
9
10 // Implementation.
11
12 #include "cl_FF.h"
13 #include "cl_RA.h"
14 #include "cl_integer.h"
15 #include "cl_I.h"
16 #include "cl_F.h"
17
18 float cl_float_approx (const cl_RA& x)
19 {
20 // Method: same as cl_RA_to_FF().
21       if (integerp(x)) {
22         DeclareType(cl_I,x);
23         return cl_float_approx(x);
24       }
25  {    // x Ratio
26       DeclareType(cl_RT,x);
27       union { ffloat eksplicit; float machine_float; } u;
28       var cl_I a = numerator(x); // +/- a
29       var const cl_I& b = denominator(x); // b
30       var cl_signean sign = -(cl_signean)minusp(a); // Vorzeichen
31       if (!(sign==0)) { a = -a; } // Betrag nehmen, liefert a
32       var sintL lendiff = (sintL)integer_length(a) // (integer-length a)
33                           - (sintL)integer_length(b); // (integer-length b)
34       if (lendiff > FF_exp_high-FF_exp_mid) // Exponent >= n-m > Obergrenze ?
35         { u.eksplicit = make_FF_word(sign,bit(FF_exp_len)-1,0); // Infinity
36           return u.machine_float;
37         }
38       if (lendiff < FF_exp_low-FF_exp_mid-2) // Exponent <= n-m+2 < Untergrenze ?
39         { u.eksplicit = make_FF_word(sign,0,0); // 0.0
40           return u.machine_float;
41         }
42       var cl_I zaehler;
43       var cl_I nenner;
44       if (lendiff >= FF_mant_len+2)
45         // n-m-25>=0
46         { nenner = ash(b,lendiff - (FF_mant_len+2)); // (ash b n-m-25)
47           zaehler = a; // a
48         }
49         else
50         { zaehler = ash(a,(FF_mant_len+2) - lendiff); // (ash a -n+m+25)
51           nenner = b; // b
52         }
53       // Division zaehler/nenner durchführen:
54       var cl_I_div_t q_r = cl_divide(zaehler,nenner);
55       var cl_I& q = q_r.quotient;
56       var cl_I& r = q_r.remainder;
57       // 2^24 <= q < 2^26, also ist q Fixnum oder Bignum mit bn_minlength Digits.
58       var uint32 mant = ((FF_mant_len+3 < cl_value_len)
59                           ? FN_to_UL(q)
60                           : cl_I_to_UL(q)
61                         );
62       if (mant >= bit(FF_mant_len+2))
63         // 2^25 <= q < 2^26, schiebe um 2 Bits nach rechts
64         { var uintL rounding_bits = mant & (bit(2)-1);
65           lendiff = lendiff+1; // Exponent := n-m+1
66           mant = mant >> 2;
67           if ( (rounding_bits < bit(1)) // 00,01 werden abgerundet
68                || ( (rounding_bits == bit(1)) // 10
69                     && (eq(r,0)) // und genau halbzahlig (r=0)
70                     && ((mant & bit(0)) ==0) // -> round-to-even
71              )    )
72             // abrunden
73             goto ab;
74             else
75             // aufrunden
76             goto auf;
77         }
78         else
79         { var uintL rounding_bit = mant & bit(0);
80           mant = mant >> 1;
81           if ( (rounding_bit == 0) // 0 wird abgerundet
82                || ( (eq(r,0)) // genau halbzahlig (r=0)
83                     && ((mant & bit(0)) ==0) // -> round-to-even
84              )    )
85             // abrunden
86             goto ab;
87             else
88             // aufrunden
89             goto auf;
90         }
91       auf:
92       mant += 1;
93       if (mant >= bit(FF_mant_len+1)) // rounding overflow?
94         { mant = mant>>1; lendiff = lendiff+1; }
95       ab:
96       // Fertig.
97       if (lendiff < (sintL)(FF_exp_low-FF_exp_mid))
98         { u.eksplicit = make_FF_word(sign,0,0); }
99       else if (lendiff > (sintL)(FF_exp_high-FF_exp_mid))
100         { u.eksplicit = make_FF_word(sign,bit(FF_exp_len)-1,0); } // Infinity
101       else
102         { u.eksplicit = make_FF_word(sign,lendiff+FF_exp_mid,mant); }
103       return u.machine_float;
104 }}