1 
2 #include <float.h>
3 #include <math.h>
4 #include <stdint.h>
5 
6 
7 #define REAL             float
8 #define REAL_ABS         fabsf
9 
10 #define REAL_MIN_NORMAL  FLT_MIN
11 #define REAL_EPSILON     FLT_EPSILON
12 #define REAL_MAX         FLT_MAX
13 #define REAL_MANT_DIG    FLT_MANT_DIG
14 
15 #define FEQREL           feqrelf
16 #include "feqrel_source.c"
17 
18 union float_t {
19     float f;
20     uint32_t w;
21 };
22 
23 int
identicalf(float x,float y)24 identicalf (float x, float y)
25 {
26     union float_t ux = { x };
27     union float_t uy = { y };
28     return ux.w == uy.w;
29 }
30 
31 float
copysignf(float x,float y)32 copysignf (float x, float y)
33 {
34     union float_t ux = { x };
35     union float_t uy = { y };
36     union float_t uz;
37     uint32_t val  = ux.w & 0x7FFFFFFF;
38     uint32_t sign = uy.w & 0x80000000;
39     uz.w = sign | val;
40     return uz.f;
41 }
42 
43 /* ported from tango/math/IEEE.d nextupf */
44 float
ieeesuccf(float x)45 ieeesuccf (float x)
46 {
47     union float_t ps = { x };
48 
49     if ((ps.w & 0x7F800000) == 0x7F800000) {
50         /* First, deal with NANs and infinity */
51         if (x == -INFINITY) return -REAL_MAX;
52         return x; /* +INF and NAN are unchanged. */
53     }
54     if (ps.w & 0x80000000)  { /* Negative number */
55         if (ps.w == 0x80000000) { /* it was negative zero */
56             ps.w = 0x00000001; /* change to smallest subnormal */
57             return ps.f;
58         }
59         --ps.w;
60     } else { /* Positive number */
61         ++ps.w;
62     }
63     return ps.f;
64 }
65 
66 /* ported from tango/math/IEEE.d nextdownf */
67 float
ieeepredf(float x)68 ieeepredf (float x)
69 {
70     return -ieeesuccf(-x);
71 }
72 
73 /* ported from tango/math/IEEE.d */
74 float
ieeemeanf(float x,float y)75 ieeemeanf (float x, float y)
76 {
77     if (!((x>=0 && y>=0) || (x<=0 && y<=0))) return NAN;
78 
79     union float_t ul;
80     union float_t xl = { x };
81     union float_t yl = { y };
82     uint32_t m = ((xl.w & 0x7FFFFFFF) + (yl.w & 0x7FFFFFFF)) >> 1;
83     m |= (xl.w & 0x80000000);
84     ul.w = m;
85 
86     return ul.f;
87 }
88 
89 float
mknanf(uint32_t payload)90 mknanf (uint32_t payload)
91 {
92     union float_t ux = { NAN };
93 
94     /* get sign, exponent, and quiet bit from NAN */
95     ux.w &= 0xFFC00000;
96 
97     /* ignore sign, exponent, and quiet bit in payload */
98     payload &= 0x003FFFFF;
99     ux.w |= payload;
100 
101     return ux.f;
102 }
103 
104 uint32_t
getnanf(float x)105 getnanf (float x)
106 {
107     union float_t payload = { x };
108 
109     /* clear sign, exponent, and quiet bit */
110     payload.w &= 0x003FFFFF;
111     return payload.w;
112 }
113