Implement `decodeDouble_Int64#` primop
[ghc.git] / rts / StgPrimFloat.c
1 /* -----------------------------------------------------------------------------
2 *
3 * (c) Lennart Augustsson
4 * (c) The GHC Team, 1998-2000
5 *
6 * Miscellaneous support for floating-point primitives
7 *
8 * ---------------------------------------------------------------------------*/
9
10 #include "PosixSource.h"
11 #include "Rts.h"
12
13 #include "StgPrimFloat.h"
14
15 #include <math.h>
16 #include <float.h>
17
18 #define IEEE_FLOATING_POINT 1
19
20 #if FLT_RADIX != 2
21 # error FLT_RADIX != 2 not supported
22 #endif
23
24 /*
25 * Encoding and decoding Doubles. Code based on the HBC code
26 * (lib/fltcode.c).
27 */
28
29 #if IEEE_FLOATING_POINT
30 #define MY_DMINEXP ((DBL_MIN_EXP) - (DBL_MANT_DIG) - 1)
31 /* DMINEXP is defined in values.h on Linux (for example) */
32 #define DHIGHBIT 0x00100000
33 #define DMSBIT 0x80000000
34
35 #define MY_FMINEXP ((FLT_MIN_EXP) - (FLT_MANT_DIG) - 1)
36 #define FHIGHBIT 0x00800000
37 #define FMSBIT 0x80000000
38 #endif
39
40 #if defined(WORDS_BIGENDIAN) || defined(FLOAT_WORDS_BIGENDIAN)
41 #define L 1
42 #define H 0
43 #else
44 #define L 0
45 #define H 1
46 #endif
47
48 #define __abs(a) (( (a) >= 0 ) ? (a) : (-(a)))
49
50 /* Special version for words */
51 StgDouble
52 __word_encodeDouble (W_ j, I_ e)
53 {
54 StgDouble r;
55
56 r = (StgDouble)j;
57
58 /* Now raise to the exponent */
59 if ( r != 0.0 ) /* Lennart suggests this avoids a bug in MIPS's ldexp */
60 r = ldexp(r, e);
61
62 return r;
63 }
64
65 /* Special version for small Integers */
66 StgDouble
67 __int_encodeDouble (I_ j, I_ e)
68 {
69 StgDouble r;
70
71 r = (StgDouble)__abs(j);
72
73 /* Now raise to the exponent */
74 if ( r != 0.0 ) /* Lennart suggests this avoids a bug in MIPS's ldexp */
75 r = ldexp(r, e);
76
77 /* sign is encoded in the size */
78 if (j < 0)
79 r = -r;
80
81 return r;
82 }
83
84 /* Special version for small Integers */
85 StgFloat
86 __int_encodeFloat (I_ j, I_ e)
87 {
88 StgFloat r;
89
90 r = (StgFloat)__abs(j);
91
92 /* Now raise to the exponent */
93 if ( r != 0.0 ) /* Lennart suggests this avoids a bug in MIPS's ldexp */
94 r = ldexp(r, e);
95
96 /* sign is encoded in the size */
97 if (j < 0)
98 r = -r;
99
100 return r;
101 }
102
103 /* Special version for small positive Integers */
104 StgFloat
105 __word_encodeFloat (W_ j, I_ e)
106 {
107 StgFloat r;
108
109 r = (StgFloat)j;
110
111 /* Now raise to the exponent */
112 if ( r != 0.0 ) /* Lennart suggests this avoids a bug in MIPS's ldexp */
113 r = ldexp(r, e);
114
115 return r;
116 }
117
118 /* This only supports IEEE floating point */
119
120 void
121 __decodeDouble_2Int (I_ *man_sign, W_ *man_high, W_ *man_low, I_ *exp, StgDouble dbl)
122 {
123 /* Do some bit fiddling on IEEE */
124 unsigned int low, high; /* assuming 32 bit ints */
125 int sign, iexp;
126 union { double d; unsigned int i[2]; } u; /* assuming 32 bit ints, 64 bit double */
127
128 ASSERT(sizeof(unsigned int ) == 4 );
129 ASSERT(sizeof(dbl ) == 8 );
130 ASSERT(sizeof(dbl ) == SIZEOF_DOUBLE);
131
132 u.d = dbl; /* grab chunks of the double */
133 low = u.i[L];
134 high = u.i[H];
135
136 if (low == 0 && (high & ~DMSBIT) == 0) {
137 *man_low = 0;
138 *man_high = 0;
139 *exp = 0L;
140 } else {
141 iexp = ((high >> 20) & 0x7ff) + MY_DMINEXP;
142 sign = high;
143
144 high &= DHIGHBIT-1;
145 if (iexp != MY_DMINEXP) /* don't add hidden bit to denorms */
146 high |= DHIGHBIT;
147 else {
148 iexp++;
149 /* A denorm, normalize the mantissa */
150 while (! (high & DHIGHBIT)) {
151 high <<= 1;
152 if (low & DMSBIT)
153 high++;
154 low <<= 1;
155 iexp--;
156 }
157 }
158 *exp = (I_) iexp;
159 *man_low = low;
160 *man_high = high;
161 *man_sign = (sign < 0) ? -1 : 1;
162 }
163 }
164
165 /* This is expected to replace uses of __decodeDouble_2Int() in the long run */
166 StgInt
167 __decodeDouble_Int64 (StgInt64 *const mantissa, const StgDouble dbl)
168 {
169 if (dbl) {
170 int exp = 0;
171 *mantissa = (StgInt64)scalbn(frexp(dbl, &exp), DBL_MANT_DIG);
172 return exp-DBL_MANT_DIG;
173 } else {
174 *mantissa = 0;
175 return 0;
176 }
177 }
178
179 /* Convenient union types for checking the layout of IEEE 754 types -
180 based on defs in GNU libc <ieee754.h>
181 */
182
183 void
184 __decodeFloat_Int (I_ *man, I_ *exp, StgFloat flt)
185 {
186 /* Do some bit fiddling on IEEE */
187 int high, sign; /* assuming 32 bit ints */
188 union { float f; int i; } u; /* assuming 32 bit float and int */
189
190 ASSERT(sizeof(int ) == 4 );
191 ASSERT(sizeof(flt ) == 4 );
192 ASSERT(sizeof(flt ) == SIZEOF_FLOAT );
193
194 u.f = flt; /* grab the float */
195 high = u.i;
196
197 if ((high & ~FMSBIT) == 0) {
198 *man = 0;
199 *exp = 0;
200 } else {
201 *exp = ((high >> 23) & 0xff) + MY_FMINEXP;
202 sign = high;
203
204 high &= FHIGHBIT-1;
205 if (*exp != MY_FMINEXP) /* don't add hidden bit to denorms */
206 high |= FHIGHBIT;
207 else {
208 (*exp)++;
209 /* A denorm, normalize the mantissa */
210 while (! (high & FHIGHBIT)) {
211 high <<= 1;
212 (*exp)--;
213 }
214 }
215 *man = high;
216 if (sign < 0)
217 *man = - *man;
218 }
219 }
220
221
222 // Local Variables:
223 // mode: C
224 // fill-column: 80
225 // indent-tabs-mode: nil
226 // c-basic-offset: 4
227 // buffer-file-coding-system: utf-8-unix
228 // End: