diff --git a/bake-test.cmd b/bake-test.cmd deleted file mode 100644 index 448f301..0000000 --- a/bake-test.cmd +++ /dev/null @@ -1 +0,0 @@ -clang test\test_%test%.c utf8.obj -Iinc -g -lciabatta.lib \ No newline at end of file diff --git a/inc/math.h b/inc/math.h index 4b91d55..a8dbb86 100644 --- a/inc/math.h +++ b/inc/math.h @@ -23,11 +23,11 @@ typedef double double_t; #define math_errhandling (MATH_ERRNO | MATH_ERREXCEPT) // Classification -#define FP_INFINITE 0 -#define FP_NAN 1 +#define FP_ZERO 0 +#define FP_SUBNORMAL 1 #define FP_NORMAL 2 -#define FP_SUBNORMAL 3 -#define FP_ZERO 4 +#define FP_INFINITE 4 +#define FP_NAN 5 int _fpclassify(double); int _fpclassifyf(float); int _fpclassifyl(long double); diff --git a/inc/stdio.h b/inc/stdio.h index fe7925e..7929dc3 100644 --- a/inc/stdio.h +++ b/inc/stdio.h @@ -27,7 +27,7 @@ typedef struct FILE FILE; typedef struct { - int64_t offset; + unsigned long long offset; mbstate_t mbstate; } fpos_t; diff --git a/src/ciabatta.c b/src/ciabatta.c index bc9834d..416d19f 100644 --- a/src/ciabatta.c +++ b/src/ciabatta.c @@ -46,7 +46,7 @@ #include "math/basic.c" #include "math/division.c" #include "math/gen_math.c" -#include "math/ieee754.c" +#include "math/bits.c" #include "math/round.c" #include "conv/digits.c" #include "conv/strpfx.c" diff --git a/src/math/bits.c b/src/math/bits.c new file mode 100644 index 0000000..baeee7d --- /dev/null +++ b/src/math/bits.c @@ -0,0 +1,106 @@ + +// This thing doesn't compile branchlessly on MSVC lmao +// This is a branchless version of fpclassify +// The way it works is by assuming the following structure of float class +// bit 2 | bit 1 | bit 0 +// e n m +// e is 1 for infinities and nan, 0 otherwise +// n is 1 for normal numbers, 0 otherwise +// m is 1 for nan and subnormal numbers, 0 otherwise +// This leaves the following values for classes: +// 0 - FP_ZERO +// 1 - FP_SUBNORMAL +// 2 - FP_NORMAL +// 4 - FP_INFINITE +// 5 - FP_NAN +int _fpclassify(f64 x) { + // First, extract bits + u64 bits = F64_BITS(x); + i64 exp = F64_BEXP(bits); + i64 mant = F64_MANT(bits); + // Get the number that's only zero when exp = 0x7ff + i64 ee = F64_BEXP_MAX - exp; + // Invert it to get number that's 1 iff x is infinity or nan + i64 e = !(ee); + // The value (ee*exp) is zero in two cases: + // exp = 0x7ff, (infinities, nan) + // exp = 0 (zero, subnormal) + // We negate this to make it so that this value is 1 only if number is not + // normal + i64 nn = !(ee*exp); + // Negate the previous thing. Now n is 1 iff number is normal + i64 n = !nn; + // m is 1 if mantissa is nonzero and the number is not normal + i64 m = !!mant & nn; + // Construct the float class + return (e<<2) | (n<<1) | m; +} + +// Same recipe as above, different constants +int _fpclassifyf(float x) { + u64 bits = F32_BITS(x); + i64 exp = F32_BEXP(bits); + i64 mant = F32_MANT(bits); + i64 ee = F32_BEXP_MAX - exp; + i64 e = !(ee); + i64 nn = !(ee*exp); + i64 n = !nn; + i64 m = !!mant & nn; + return (e<<2) | (n<<1) | m; +} + +int _fpclassifyl(fl64 x) { + return _fpclassify((f64)x); +} + + +int _signbit(f64 x) { + union { + f64 d; + uint64_t i; + } y = { x }; + return y.i>>63; +} + +int _signbitf(float x) { + union { + float f; + uint32_t i; + } y = { x }; + return y.i>>31; +} + +int _signbitl(fl64 x) { + return _signbit(x); +} + +f64 copysign(f64 x, f64 y) { + union {f64 f; uint64_t i;} ux={x}, uy={y}; + ux.i &= ~(1ULL<<63); + ux.i |= uy.i & (1ULL<<63); + return ux.f; +} + +float copysignf(float x, float y) { + union {float f; uint32_t i;} ux={x}, uy={y}; + ux.i &= 0x7fffffff; + ux.i |= uy.i & 0x80000000; + return ux.f; +} + +fl64 copysignl(fl64 x, fl64 y) { + return copysign(x, y); +} + +f64 nan(const char *s) { + return NAN; +} + +float nanf(const char *s) { + return NAN; +} + +fl64 nanl(const char *s) { + return NAN; +} + diff --git a/src/math/ieee754.c b/src/math/ieee754.c deleted file mode 100644 index bc93923..0000000 --- a/src/math/ieee754.c +++ /dev/null @@ -1,77 +0,0 @@ - -#include -#include -#include -#include - -int _fpclassify(double x) { - union {double f; uint64_t i;} u = {x}; - int e = u.i>>52 & 0x7ff; - if (!e) return u.i<<1 ? FP_SUBNORMAL : FP_ZERO; - if (e==0x7ff) return u.i<<12 ? FP_NAN : FP_INFINITE; - return FP_NORMAL; -} - -int _fpclassifyf(float x) { - union {float f; uint32_t i;} u = {x}; - int e = u.i>>23 & 0xff; - if (!e) return u.i<<1 ? FP_SUBNORMAL : FP_ZERO; - if (e==0xff) return u.i<<9 ? FP_NAN : FP_INFINITE; - return FP_NORMAL; -} - -int _fpclassifyl(long double x) { - return _fpclassify(x); -} - - -int _signbit(double x) { - union { - double d; - uint64_t i; - } y = { x }; - return y.i>>63; -} - -int _signbitf(float x) { - union { - float f; - uint32_t i; - } y = { x }; - return y.i>>31; -} - -int _signbitl(long double x) { - return _signbit(x); -} - -double copysign(double x, double y) { - union {double f; uint64_t i;} ux={x}, uy={y}; - ux.i &= ~(1ULL<<63); - ux.i |= uy.i & (1ULL<<63); - return ux.f; -} - -float copysignf(float x, float y) { - union {float f; uint32_t i;} ux={x}, uy={y}; - ux.i &= 0x7fffffff; - ux.i |= uy.i & 0x80000000; - return ux.f; -} - -long double copysignl(long double x, long double y) { - return copysign(x, y); -} - -double nan(const char *s) { - return NAN; -} - -float nanf(const char *s) { - return NAN; -} - -long double nanl(const char *s) { - return NAN; -} - diff --git a/src/util.c b/src/util.c index abcb215..b517e6e 100644 --- a/src/util.c +++ b/src/util.c @@ -44,8 +44,30 @@ typedef wchar_t wchar; #define STR_(a) #a #define STR(a) STR_(a) -#define DOUBLE_BITS(x) ((union {double f; u64 i;}){x}).i -#define DOUBLE_CONS(x) ((union {double f; u64 i;}){x}).f -#define FLOAT_BITS(x) ((union {float f; u32 i;}){x}).i -#define FLOAT_CONS(x) ((union {float f; u32 i;}){x}).f +#define F64_BITS(x) ((union {f64 f; u64 i;}){x}).i +#define F64_CONS(x) ((union {f64 f; u64 i;}){x}).f +#define F32_BITS(x) ((union {f32 f; u32 i;}){x}).i +#define F32_CONS(x) ((union {f32 f; u32 i;}){x}).f +#define F64_MANT_MASK UINT64_C(0xfffffffffffff) +#define F64_MANT_MAX UINT64_C(0xfffffffffffff) +#define F64_MANT_BITS 52 +#define F64_BEXP_BITS 11 +#define F64_BEXP_MASK 0x7ff +#define F64_BEXP_MAX 0x7ff + +#define F32_MANT_MASK 0x7fffff +#define F32_MANT_MAX 0x7fffff +#define F32_MANT_BITS 23 +#define F32_BEXP_BITS 8 +#define F32_BEXP_MASK 0xff +#define F32_BEXP_MAX 0xff + + +#define F64_SIGN(bits) (bits >> (F64_MANT_BITS + F64_BEXP_BITS)) +#define F64_BEXP(bits) ((bits >> F64_MANT_BITS) & F64_BEXP_MASK) +#define F64_MANT(bits) ((bits) & F64_MANT_MASK) + +#define F32_SIGN(bits) (bits >> (F32_MANT_BITS + F32_BEXP_BITS)) +#define F32_BEXP(bits) ((bits >> F32_MANT_BITS) & F32_BEXP_MASK) +#define F32_MANT(bits) ((bits) & F32_MANT_MASK) diff --git a/test.cmd b/test.cmd new file mode 100644 index 0000000..09d65ed --- /dev/null +++ b/test.cmd @@ -0,0 +1,2 @@ + +clang test\%1 utf8.obj -Iinc -g -lciabatta.lib -nostdlib \ No newline at end of file diff --git a/test/test_math.c b/test/math.c similarity index 99% rename from test/test_math.c rename to test/math.c index ab73bb0..a1ee11a 100644 --- a/test/test_math.c +++ b/test/math.c @@ -5,6 +5,7 @@ #include #include #include +#include const char *show_classification(double x) { switch(fpclassify(x)) {