]> www.ginac.de Git - cln.git/blob - src/float/sfloat/elem/cl_SF_from_I.cc
* Removed internal gmp/ directory and other traces of it like $GMP_INCLUDES.
[cln.git] / src / float / sfloat / elem / cl_SF_from_I.cc
1 // cl_I_to_SF().
2
3 // General includes.
4 #include "cl_sysdep.h"
5
6 // Specification.
7 #include "cl_SF.h"
8
9
10 // Implementation.
11
12 #include "cl_integer.h"
13 #include "cl_I.h"
14 #include "cl_DS.h"
15
16 const cl_SF cl_I_to_SF (const cl_I& x)
17 {
18 // Methode:
19 // x=0 -> Ergebnis 0.0
20 // Merke Vorzeichen von x.
21 // x:=(abs x)
22 // Exponent:=(integer-length x)
23 //   Greife die 18 höchstwertigen Bits heraus (angeführt von einer 1).
24 //   Runde das letzte Bit weg:
25 //     Bit 0 = 0 -> abrunden,
26 //     Bit 0 = 1 und Rest =0 -> round-to-even,
27 //     Bit 0 = 1 und Rest >0 -> aufrunden.
28 //   Dabei um ein Bit nach rechts schieben.
29 //   Bei Aufrundung auf 2^17 (rounding overflow) Mantisse um 1 Bit nach rechts
30 //     schieben und Exponent incrementieren.
31       if (eq(x,0)) { return SF_0; }
32       var cl_signean sign = -(cl_signean)minusp(x); // Vorzeichen
33       var cl_I abs_x = (sign==0 ? x : -x);
34       var uintL exp = integer_length(abs_x); // (integer-length x)
35       // NDS zu |x|>0 bilden:
36       var const uintD* MSDptr;
37       var uintC len;
38       I_to_NDS_nocopy(abs_x, MSDptr=,len=,,cl_false,);
39       // MSDptr/len/LSDptr ist die NDS zu x, len>0.
40       // Führende Digits holen: Brauche SF_mant_len+1 Bits, dazu intDsize
41       // Bits (die NDS kann mit bis zu intDsize Nullbits anfangen).
42       // Dann werden diese Bits um (exp mod intDsize) nach rechts geschoben.
43       var uintD msd = msprefnext(MSDptr); // erstes Digit
44       #if (intDsize==64)
45       var uintD msdd = 0; // weiteres Digit
46       if (--len == 0) goto ok;
47       msdd = msprefnext(MSDptr);
48       #else // (intDsize<=32)
49       var uint32 msdd = 0; // weitere min(len-1,32/intDsize) Digits
50       #define NEXT_DIGIT(i)  \
51         { if (--len == 0) goto ok;                                   \
52           msdd |= (uint32)msprefnext(MSDptr) << (32-(i+1)*intDsize); \
53         }
54       DOCONSTTIMES(32/intDsize,NEXT_DIGIT);
55       #undef NEXT_DIGIT
56       #endif
57       --len; ok:
58       #if (intDsize==64)
59       // Die NDS besteht aus msd, msdd, und len weiteren Digits.
60       // Das höchste in 2^intDsize*msd+msdd gesetzte Bit ist Bit Nummer
61       // intDsize-1 + (exp mod intDsize).
62       var uintL shiftcount = exp % intDsize;
63       var uint64 mant = // führende 64 Bits
64         (shiftcount==0
65          ? msdd
66          : ((msd << (64-shiftcount)) | (msdd >> shiftcount))
67         );
68       // Das höchste in mant gesetzte Bit ist Bit Nummer 63.
69       if ( ((mant & bit(62-SF_mant_len)) ==0) // Bit 46 =0 -> abrunden
70            || ( ((mant & (bit(62-SF_mant_len)-1)) ==0) // Bit 46 =1 und Bits 45..0 =0
71                 && ((msdd & (bit(shiftcount)-1)) ==0) // und weitere Bits aus msdd =0
72                 && (!test_loop_msp(MSDptr,len)) // und alle weiteren Digits =0
73                 // round-to-even, je nach Bit 47 :
74                 && ((mant & bit(63-SF_mant_len)) ==0)
75          )    )
76         // abrunden
77         { mant = mant >> (63-SF_mant_len); }
78         else
79         // aufrunden
80         { mant = mant >> (63-SF_mant_len);
81           mant += 1;
82           if (mant >= bit(SF_mant_len+1)) // rounding overflow?
83             { mant = mant>>1; exp = exp+1; }
84         }
85       #else
86       // Die NDS besteht aus msd, msdd, und len weiteren Digits.
87       // Das höchste in 2^32*msd+msdd gesetzte Bit ist Bit Nummer
88       // 31 + (exp mod intDsize).
89       var uintL shiftcount = exp % intDsize;
90       var uint32 mant = // führende 32 Bits
91         (shiftcount==0
92          ? msdd
93          : (((uint32)msd << (32-shiftcount)) | (msdd >> shiftcount))
94         );
95       // Das höchste in mant gesetzte Bit ist Bit Nummer 31.
96       if ( ((mant & bit(30-SF_mant_len)) ==0) // Bit 14 =0 -> abrunden
97            || ( ((mant & (bit(30-SF_mant_len)-1)) ==0) // Bit 14 =1 und Bits 13..0 =0
98                 && ((msdd & (bit(shiftcount)-1)) ==0) // und weitere Bits aus msdd =0
99                 && (!test_loop_msp(MSDptr,len)) // und alle weiteren Digits =0
100                 // round-to-even, je nach Bit 15 :
101                 && ((mant & bit(31-SF_mant_len)) ==0)
102          )    )
103         // abrunden
104         { mant = mant >> (31-SF_mant_len); }
105         else
106         // aufrunden
107         { mant = mant >> (31-SF_mant_len);
108           mant += 1;
109           if (mant >= bit(SF_mant_len+1)) // rounding overflow?
110             { mant = mant>>1; exp = exp+1; }
111         }
112       #endif
113       return encode_SF(sign,(sintL)exp,mant);
114 }