1const std = @import("../std.zig");
2const math = std.math;
3const expect = std.testing.expect;
4const expectEqual = std.testing.expectEqual;
5const expectApproxEqAbs = std.testing.expectApproxEqAbs;
6
7pub fn Frexp(comptime T: type) type {
8 return struct {
9 significand: T,
10 exponent: i32,
11 };
12}
13
14/// Breaks x into a normalized fraction and an integral power of two.
15/// f == frac * 2^exp, with |frac| in the interval [0.5, 1).
16///
17/// Special Cases:
18/// - frexp(+-0) = +-0, 0
19/// - frexp(+-inf) = +-inf, 0
20/// - frexp(nan) = nan, undefined
21pub fn frexp(x: anytype) Frexp(@TypeOf(x)) {
22 const T: type = @TypeOf(x);
23
24 const bits: comptime_int = @typeInfo(T).float.bits;
25 const Int: type = @Int(.unsigned, bits);
26
27 const exp_bits: comptime_int = math.floatExponentBits(T);
28 const mant_bits: comptime_int = math.floatMantissaBits(T);
29 const frac_bits: comptime_int = math.floatFractionalBits(T);
30 const exp_min: comptime_int = math.floatExponentMin(T);
31
32 const ExpInt: type = @Int(.unsigned, exp_bits);
33 const MantInt: type = @Int(.unsigned, mant_bits);
34 const FracInt: type = @Int(.unsigned, frac_bits);
35
36 const unreal_exponent: comptime_int = (1 << exp_bits) - 1;
37 const bias: comptime_int = (1 << (exp_bits - 1)) - 2;
38 const exp_mask: comptime_int = unreal_exponent << mant_bits;
39 const zero_exponent: comptime_int = bias << mant_bits;
40 const sign_mask: comptime_int = 1 << (bits - 1);
41 const not_exp: comptime_int = ~@as(Int, exp_mask);
42 const ones_place: comptime_int = mant_bits - frac_bits;
43 const extra_denorm_shift: comptime_int = 1 - ones_place;
44
45 var result: Frexp(T) = undefined;
46 var v: Int = @bitCast(x);
47
48 const m: MantInt = @truncate(v);
49 const e: ExpInt = @truncate(v >> mant_bits);
50
51 switch (e) {
52 0 => {
53 if (m != 0) {
54 // subnormal
55 const offset = @clz(m);
56 const shift = offset + extra_denorm_shift;
57
58 v &= sign_mask;
59 v |= zero_exponent;
60 v |= math.shl(MantInt, m, shift);
61
62 result.exponent = exp_min - @as(i32, offset) + ones_place;
63 } else {
64 // +-0 = (+-0, 0)
65 result.exponent = 0;
66 }
67 },
68 unreal_exponent => {
69 // +-nan -> {+-nan, undefined}
70 result.exponent = undefined;
71
72 // +-inf -> {+-inf, 0}
73 if (@as(FracInt, @truncate(v)) == 0)
74 result.exponent = 0;
75 },
76 else => {
77 // normal
78 v &= not_exp;
79 v |= zero_exponent;
80 result.exponent = @as(i32, e) - bias;
81 },
82 }
83
84 result.significand = @bitCast(v);
85 return result;
86}
87
88/// Generate a namespace of tests for frexp on values of the given type
89fn FrexpTests(comptime Float: type) type {
90 return struct {
91 const T = Float;
92 test "normal" {
93 const epsilon = 1e-6;
94 var r: Frexp(T) = undefined;
95
96 r = frexp(@as(T, 1.3));
97 try expectApproxEqAbs(0.65, r.significand, epsilon);
98 try expectEqual(1, r.exponent);
99
100 r = frexp(@as(T, 78.0234));
101 try expectApproxEqAbs(0.609558, r.significand, epsilon);
102 try expectEqual(7, r.exponent);
103
104 r = frexp(@as(T, -1234.5678));
105 try expectEqual(11, r.exponent);
106 try expectApproxEqAbs(-0.602816, r.significand, epsilon);
107 }
108 test "max" {
109 const exponent = math.floatExponentMax(T) + 1;
110 const significand = 1.0 - math.floatEps(T) / 2;
111 const r: Frexp(T) = frexp(math.floatMax(T));
112 try expectEqual(exponent, r.exponent);
113 try expectEqual(significand, r.significand);
114 }
115 test "min" {
116 const exponent = math.floatExponentMin(T) + 1;
117 const r: Frexp(T) = frexp(math.floatMin(T));
118 try expectEqual(exponent, r.exponent);
119 try expectEqual(0.5, r.significand);
120 }
121 test "subnormal" {
122 const normal_min_exponent = math.floatExponentMin(T) + 1;
123 const exponent = normal_min_exponent - math.floatFractionalBits(T);
124 const r: Frexp(T) = frexp(math.floatTrueMin(T));
125 try expectEqual(exponent, r.exponent);
126 try expectEqual(0.5, r.significand);
127 }
128 test "zero" {
129 var r: Frexp(T) = undefined;
130
131 r = frexp(@as(T, 0.0));
132 try expectEqual(0, r.exponent);
133 try expect(math.isPositiveZero(r.significand));
134
135 r = frexp(@as(T, -0.0));
136 try expectEqual(0, r.exponent);
137 try expect(math.isNegativeZero(r.significand));
138 }
139 test "inf" {
140 var r: Frexp(T) = undefined;
141
142 r = frexp(math.inf(T));
143 try expectEqual(0, r.exponent);
144 try expect(math.isPositiveInf(r.significand));
145
146 r = frexp(-math.inf(T));
147 try expectEqual(0, r.exponent);
148 try expect(math.isNegativeInf(r.significand));
149 }
150 test "nan" {
151 const r: Frexp(T) = frexp(math.nan(T));
152 try expect(math.isNan(r.significand));
153 }
154 };
155}
156
157// Generate tests for each floating point type
158comptime {
159 for ([_]type{ f16, f32, f64, f80, f128 }) |T| {
160 _ = FrexpTests(T);
161 }
162}
163
164test frexp {
165 inline for ([_]type{ f16, f32, f64, f80, f128 }) |T| {
166 const max_exponent = math.floatExponentMax(T) + 1;
167 const min_exponent = math.floatExponentMin(T) + 1;
168 const truemin_exponent = min_exponent - math.floatFractionalBits(T);
169
170 var result: Frexp(T) = undefined;
171 comptime var x: T = undefined;
172
173 // basic usage
174 // value -> {significand, exponent},
175 // value == significand * (2 ^ exponent)
176 x = 1234.5678;
177 result = frexp(x);
178 try expectEqual(11, result.exponent);
179 try expectApproxEqAbs(0.602816, result.significand, 1e-6);
180 try expectEqual(x, math.ldexp(result.significand, result.exponent));
181
182 // float maximum
183 x = math.floatMax(T);
184 result = frexp(x);
185 try expectEqual(max_exponent, result.exponent);
186 try expectEqual(1.0 - math.floatEps(T) / 2, result.significand);
187 try expectEqual(x, math.ldexp(result.significand, result.exponent));
188
189 // float minimum
190 x = math.floatMin(T);
191 result = frexp(x);
192 try expectEqual(min_exponent, result.exponent);
193 try expectEqual(0.5, result.significand);
194 try expectEqual(x, math.ldexp(result.significand, result.exponent));
195
196 // float true minimum
197 // subnormal -> {normal, exponent}
198 x = math.floatTrueMin(T);
199 result = frexp(x);
200 try expectEqual(truemin_exponent, result.exponent);
201 try expectEqual(0.5, result.significand);
202 try expectEqual(x, math.ldexp(result.significand, result.exponent));
203
204 // infinity -> {infinity, zero} (+)
205 result = frexp(math.inf(T));
206 try expectEqual(0, result.exponent);
207 try expect(math.isPositiveInf(result.significand));
208
209 // infinity -> {infinity, zero} (-)
210 result = frexp(-math.inf(T));
211 try expectEqual(0, result.exponent);
212 try expect(math.isNegativeInf(result.significand));
213
214 // zero -> {zero, zero} (+)
215 result = frexp(@as(T, 0.0));
216 try expectEqual(0, result.exponent);
217 try expect(math.isPositiveZero(result.significand));
218
219 // zero -> {zero, zero} (-)
220 result = frexp(@as(T, -0.0));
221 try expectEqual(0, result.exponent);
222 try expect(math.isNegativeZero(result.significand));
223
224 // nan -> {nan, undefined}
225 result = frexp(math.nan(T));
226 try expect(math.isNan(result.significand));
227 }
228}