1/* Native implementation of soft float functions. Only a single status
2   context is supported */
3#include "softfloat.h"
4#include <math.h>
5#if defined(CONFIG_SOLARIS)
6#include <fenv.h>
7#endif
8
9void set_float_rounding_mode(int val STATUS_PARAM)
10{
11    STATUS(float_rounding_mode) = val;
12#if (defined(CONFIG_BSD) && !defined(__APPLE__) && !defined(__GLIBC__)) || \
13    (defined(CONFIG_SOLARIS) && CONFIG_SOLARIS_VERSION < 10)
14    fpsetround(val);
15#else
16    fesetround(val);
17#endif
18}
19
20#ifdef FLOATX80
21void set_floatx80_rounding_precision(int val STATUS_PARAM)
22{
23    STATUS(floatx80_rounding_precision) = val;
24}
25#endif
26
27#if defined(CONFIG_BSD) || \
28    (defined(CONFIG_SOLARIS) && CONFIG_SOLARIS_VERSION < 10)
29#define lrint(d)		((int32_t)rint(d))
30#define llrint(d)		((int64_t)rint(d))
31#define lrintf(f)		((int32_t)rint(f))
32#define llrintf(f)		((int64_t)rint(f))
33#define sqrtf(f)		((float)sqrt(f))
34#define remainderf(fa, fb)	((float)remainder(fa, fb))
35#define rintf(f)		((float)rint(f))
36#if !defined(__sparc__) && \
37    (defined(CONFIG_SOLARIS) && CONFIG_SOLARIS_VERSION < 10)
38extern long double rintl(long double);
39extern long double scalbnl(long double, int);
40
41long long
42llrintl(long double x) {
43	return ((long long) rintl(x));
44}
45
46long
47lrintl(long double x) {
48	return ((long) rintl(x));
49}
50
51long double
52ldexpl(long double x, int n) {
53	return (scalbnl(x, n));
54}
55#endif
56#endif
57
58#if defined(_ARCH_PPC)
59
60/* correct (but slow) PowerPC rint() (glibc version is incorrect) */
61static double qemu_rint(double x)
62{
63    double y = 4503599627370496.0;
64    if (fabs(x) >= y)
65        return x;
66    if (x < 0)
67        y = -y;
68    y = (x + y) - y;
69    if (y == 0.0)
70        y = copysign(y, x);
71    return y;
72}
73
74#define rint qemu_rint
75#endif
76
77/*----------------------------------------------------------------------------
78| Software IEC/IEEE integer-to-floating-point conversion routines.
79*----------------------------------------------------------------------------*/
80float32 int32_to_float32(int v STATUS_PARAM)
81{
82    return (float32)v;
83}
84
85float32 uint32_to_float32(unsigned int v STATUS_PARAM)
86{
87    return (float32)v;
88}
89
90float64 int32_to_float64(int v STATUS_PARAM)
91{
92    return (float64)v;
93}
94
95float64 uint32_to_float64(unsigned int v STATUS_PARAM)
96{
97    return (float64)v;
98}
99
100#ifdef FLOATX80
101floatx80 int32_to_floatx80(int v STATUS_PARAM)
102{
103    return (floatx80)v;
104}
105#endif
106float32 int64_to_float32( int64_t v STATUS_PARAM)
107{
108    return (float32)v;
109}
110float32 uint64_to_float32( uint64_t v STATUS_PARAM)
111{
112    return (float32)v;
113}
114float64 int64_to_float64( int64_t v STATUS_PARAM)
115{
116    return (float64)v;
117}
118float64 uint64_to_float64( uint64_t v STATUS_PARAM)
119{
120    return (float64)v;
121}
122#ifdef FLOATX80
123floatx80 int64_to_floatx80( int64_t v STATUS_PARAM)
124{
125    return (floatx80)v;
126}
127#endif
128
129/* XXX: this code implements the x86 behaviour, not the IEEE one.  */
130#if HOST_LONG_BITS == 32
131static inline int long_to_int32(long a)
132{
133    return a;
134}
135#else
136static inline int long_to_int32(long a)
137{
138    if (a != (int32_t)a)
139        a = 0x80000000;
140    return a;
141}
142#endif
143
144/*----------------------------------------------------------------------------
145| Software IEC/IEEE single-precision conversion routines.
146*----------------------------------------------------------------------------*/
147int float32_to_int32( float32 a STATUS_PARAM)
148{
149    return long_to_int32(lrintf(a));
150}
151int float32_to_int32_round_to_zero( float32 a STATUS_PARAM)
152{
153    return (int)a;
154}
155int64_t float32_to_int64( float32 a STATUS_PARAM)
156{
157    return llrintf(a);
158}
159
160int64_t float32_to_int64_round_to_zero( float32 a STATUS_PARAM)
161{
162    return (int64_t)a;
163}
164
165float64 float32_to_float64( float32 a STATUS_PARAM)
166{
167    return a;
168}
169#ifdef FLOATX80
170floatx80 float32_to_floatx80( float32 a STATUS_PARAM)
171{
172    return a;
173}
174#endif
175
176unsigned int float32_to_uint32( float32 a STATUS_PARAM)
177{
178    int64_t v;
179    unsigned int res;
180
181    v = llrintf(a);
182    if (v < 0) {
183        res = 0;
184    } else if (v > 0xffffffff) {
185        res = 0xffffffff;
186    } else {
187        res = v;
188    }
189    return res;
190}
191unsigned int float32_to_uint32_round_to_zero( float32 a STATUS_PARAM)
192{
193    int64_t v;
194    unsigned int res;
195
196    v = (int64_t)a;
197    if (v < 0) {
198        res = 0;
199    } else if (v > 0xffffffff) {
200        res = 0xffffffff;
201    } else {
202        res = v;
203    }
204    return res;
205}
206
207/*----------------------------------------------------------------------------
208| Software IEC/IEEE single-precision operations.
209*----------------------------------------------------------------------------*/
210float32 float32_round_to_int( float32 a STATUS_PARAM)
211{
212    return rintf(a);
213}
214
215float32 float32_rem( float32 a, float32 b STATUS_PARAM)
216{
217    return remainderf(a, b);
218}
219
220float32 float32_sqrt( float32 a STATUS_PARAM)
221{
222    return sqrtf(a);
223}
224int float32_compare( float32 a, float32 b STATUS_PARAM )
225{
226    if (a < b) {
227        return float_relation_less;
228    } else if (a == b) {
229        return float_relation_equal;
230    } else if (a > b) {
231        return float_relation_greater;
232    } else {
233        return float_relation_unordered;
234    }
235}
236int float32_compare_quiet( float32 a, float32 b STATUS_PARAM )
237{
238    if (isless(a, b)) {
239        return float_relation_less;
240    } else if (a == b) {
241        return float_relation_equal;
242    } else if (isgreater(a, b)) {
243        return float_relation_greater;
244    } else {
245        return float_relation_unordered;
246    }
247}
248int float32_is_signaling_nan( float32 a1)
249{
250    float32u u;
251    uint32_t a;
252    u.f = a1;
253    a = u.i;
254    return ( ( ( a>>22 ) & 0x1FF ) == 0x1FE ) && ( a & 0x003FFFFF );
255}
256
257int float32_is_quiet_nan( float32 a1 )
258{
259    float32u u;
260    uint64_t a;
261    u.f = a1;
262    a = u.i;
263    return ( 0xFF800000 < ( a<<1 ) );
264}
265
266int float32_is_any_nan( float32 a1 )
267{
268    float32u u;
269    uint32_t a;
270    u.f = a1;
271    a = u.i;
272    return (a & ~(1 << 31)) > 0x7f800000U;
273}
274
275/*----------------------------------------------------------------------------
276| Software IEC/IEEE double-precision conversion routines.
277*----------------------------------------------------------------------------*/
278int float64_to_int32( float64 a STATUS_PARAM)
279{
280    return long_to_int32(lrint(a));
281}
282int float64_to_int32_round_to_zero( float64 a STATUS_PARAM)
283{
284    return (int)a;
285}
286int64_t float64_to_int64( float64 a STATUS_PARAM)
287{
288    return llrint(a);
289}
290int64_t float64_to_int64_round_to_zero( float64 a STATUS_PARAM)
291{
292    return (int64_t)a;
293}
294float32 float64_to_float32( float64 a STATUS_PARAM)
295{
296    return a;
297}
298#ifdef FLOATX80
299floatx80 float64_to_floatx80( float64 a STATUS_PARAM)
300{
301    return a;
302}
303#endif
304#ifdef FLOAT128
305float128 float64_to_float128( float64 a STATUS_PARAM)
306{
307    return a;
308}
309#endif
310
311unsigned int float64_to_uint32( float64 a STATUS_PARAM)
312{
313    int64_t v;
314    unsigned int res;
315
316    v = llrint(a);
317    if (v < 0) {
318        res = 0;
319    } else if (v > 0xffffffff) {
320        res = 0xffffffff;
321    } else {
322        res = v;
323    }
324    return res;
325}
326unsigned int float64_to_uint32_round_to_zero( float64 a STATUS_PARAM)
327{
328    int64_t v;
329    unsigned int res;
330
331    v = (int64_t)a;
332    if (v < 0) {
333        res = 0;
334    } else if (v > 0xffffffff) {
335        res = 0xffffffff;
336    } else {
337        res = v;
338    }
339    return res;
340}
341uint64_t float64_to_uint64 (float64 a STATUS_PARAM)
342{
343    int64_t v;
344
345    v = llrint(a + (float64)INT64_MIN);
346
347    return v - INT64_MIN;
348}
349uint64_t float64_to_uint64_round_to_zero (float64 a STATUS_PARAM)
350{
351    int64_t v;
352
353    v = (int64_t)(a + (float64)INT64_MIN);
354
355    return v - INT64_MIN;
356}
357
358/*----------------------------------------------------------------------------
359| Software IEC/IEEE double-precision operations.
360*----------------------------------------------------------------------------*/
361#if defined(__sun__) && \
362    (defined(CONFIG_SOLARIS) && CONFIG_SOLARIS_VERSION < 10)
363static inline float64 trunc(float64 x)
364{
365    return x < 0 ? -floor(-x) : floor(x);
366}
367#endif
368float64 float64_trunc_to_int( float64 a STATUS_PARAM )
369{
370    return trunc(a);
371}
372
373float64 float64_round_to_int( float64 a STATUS_PARAM )
374{
375    return rint(a);
376}
377
378float64 float64_rem( float64 a, float64 b STATUS_PARAM)
379{
380    return remainder(a, b);
381}
382
383float64 float64_sqrt( float64 a STATUS_PARAM)
384{
385    return sqrt(a);
386}
387int float64_compare( float64 a, float64 b STATUS_PARAM )
388{
389    if (a < b) {
390        return float_relation_less;
391    } else if (a == b) {
392        return float_relation_equal;
393    } else if (a > b) {
394        return float_relation_greater;
395    } else {
396        return float_relation_unordered;
397    }
398}
399int float64_compare_quiet( float64 a, float64 b STATUS_PARAM )
400{
401    if (isless(a, b)) {
402        return float_relation_less;
403    } else if (a == b) {
404        return float_relation_equal;
405    } else if (isgreater(a, b)) {
406        return float_relation_greater;
407    } else {
408        return float_relation_unordered;
409    }
410}
411int float64_is_signaling_nan( float64 a1)
412{
413    float64u u;
414    uint64_t a;
415    u.f = a1;
416    a = u.i;
417    return
418           ( ( ( a>>51 ) & 0xFFF ) == 0xFFE )
419        && ( a & LIT64( 0x0007FFFFFFFFFFFF ) );
420
421}
422
423int float64_is_quiet_nan( float64 a1 )
424{
425    float64u u;
426    uint64_t a;
427    u.f = a1;
428    a = u.i;
429
430    return ( LIT64( 0xFFF0000000000000 ) < (uint64_t) ( a<<1 ) );
431
432}
433
434int float64_is_any_nan( float64 a1 )
435{
436    float64u u;
437    uint64_t a;
438    u.f = a1;
439    a = u.i;
440
441    return (a & ~(1ULL << 63)) > LIT64 (0x7FF0000000000000 );
442}
443
444#ifdef FLOATX80
445
446/*----------------------------------------------------------------------------
447| Software IEC/IEEE extended double-precision conversion routines.
448*----------------------------------------------------------------------------*/
449int floatx80_to_int32( floatx80 a STATUS_PARAM)
450{
451    return long_to_int32(lrintl(a));
452}
453int floatx80_to_int32_round_to_zero( floatx80 a STATUS_PARAM)
454{
455    return (int)a;
456}
457int64_t floatx80_to_int64( floatx80 a STATUS_PARAM)
458{
459    return llrintl(a);
460}
461int64_t floatx80_to_int64_round_to_zero( floatx80 a STATUS_PARAM)
462{
463    return (int64_t)a;
464}
465float32 floatx80_to_float32( floatx80 a STATUS_PARAM)
466{
467    return a;
468}
469float64 floatx80_to_float64( floatx80 a STATUS_PARAM)
470{
471    return a;
472}
473
474/*----------------------------------------------------------------------------
475| Software IEC/IEEE extended double-precision operations.
476*----------------------------------------------------------------------------*/
477floatx80 floatx80_round_to_int( floatx80 a STATUS_PARAM)
478{
479    return rintl(a);
480}
481floatx80 floatx80_rem( floatx80 a, floatx80 b STATUS_PARAM)
482{
483    return remainderl(a, b);
484}
485floatx80 floatx80_sqrt( floatx80 a STATUS_PARAM)
486{
487    return sqrtl(a);
488}
489int floatx80_compare( floatx80 a, floatx80 b STATUS_PARAM )
490{
491    if (a < b) {
492        return float_relation_less;
493    } else if (a == b) {
494        return float_relation_equal;
495    } else if (a > b) {
496        return float_relation_greater;
497    } else {
498        return float_relation_unordered;
499    }
500}
501int floatx80_compare_quiet( floatx80 a, floatx80 b STATUS_PARAM )
502{
503    if (isless(a, b)) {
504        return float_relation_less;
505    } else if (a == b) {
506        return float_relation_equal;
507    } else if (isgreater(a, b)) {
508        return float_relation_greater;
509    } else {
510        return float_relation_unordered;
511    }
512}
513int floatx80_is_signaling_nan( floatx80 a1)
514{
515    floatx80u u;
516    uint64_t aLow;
517    u.f = a1;
518
519    aLow = u.i.low & ~ LIT64( 0x4000000000000000 );
520    return
521           ( ( u.i.high & 0x7FFF ) == 0x7FFF )
522        && (uint64_t) ( aLow<<1 )
523        && ( u.i.low == aLow );
524}
525
526int floatx80_is_quiet_nan( floatx80 a1 )
527{
528    floatx80u u;
529    u.f = a1;
530    return ( ( u.i.high & 0x7FFF ) == 0x7FFF ) && (uint64_t) ( u.i.low<<1 );
531}
532
533int floatx80_is_any_nan( floatx80 a1 )
534{
535    floatx80u u;
536    u.f = a1;
537    return ((u.i.high & 0x7FFF) == 0x7FFF) && ( u.i.low<<1 );
538}
539
540#endif
541