Add a D implementation for std.math.ilogb, templated for different floating point types.

This commit is contained in:
Johan Engelen 2015-04-12 23:43:50 +02:00
parent b34ff324f0
commit 98638abfc4

View file

@ -2595,6 +2595,22 @@ unittest
}
}
static import core.bitop;
static if (size_t.sizeof == 4)
{
private int bsr_ulong(ulong x) @trusted pure nothrow @nogc
{
size_t msb = x >> 32;
size_t lsb = cast(size_t) x;
if (msb)
return core.bitop.bsr(msb) + 32;
else
return core.bitop.bsr(lsb);
}
}
else
private alias bsr_ulong = core.bitop.bsr;
/******************************************
* Extracts the exponent of x as a signed integral value.
*
@ -2608,74 +2624,154 @@ unittest
* $(TR $(TD $(NAN)) $(TD FP_ILOGBNAN) $(TD no))
* )
*/
int ilogb(real x) @trusted nothrow @nogc
int ilogb(T)(const T x) @trusted nothrow @nogc
if(isFloatingPoint!T)
{
version (Win64_DMD_InlineAsm)
alias F = floatTraits!T;
union floatBits
{
asm pure nothrow @nogc
T rv;
ushort[T.sizeof/2] vu;
uint[T.sizeof/4] vui;
static if(T.sizeof >= 8)
long[T.sizeof/8] vl;
}
floatBits y = void;
y.rv = x;
int ex = y.vu[F.EXPPOS_SHORT] & F.EXPMASK;
static if (F.realFormat == RealFormat.ieeeExtended)
{
if (ex)
{
naked ;
fld real ptr [RCX] ;
fxam ;
fstsw AX ;
and AH,0x45 ;
cmp AH,0x40 ;
jz Lzeronan ;
cmp AH,5 ;
jz Linfinity ;
cmp AH,1 ;
jz Lzeronan ;
fxtract ;
fstp ST(0) ;
fistp dword ptr 8[RSP] ;
mov EAX,8[RSP] ;
ret ;
Lzeronan:
mov EAX,0x80000000 ;
fstp ST(0) ;
ret ;
Linfinity:
mov EAX,0x7FFFFFFF ;
fstp ST(0) ;
ret ;
// If exponent is non-zero
if (ex == F.EXPMASK) // infinity or NaN
{
if (y.vl[0] & 0x7FFF_FFFF_FFFF_FFFF) // NaN
return FP_ILOGBNAN;
else // +-infinity
return int.max;
}
else
{
return ex - F.EXPBIAS - 1;
}
}
else if (!y.vl[0])
{
// vf is +-0.0
return FP_ILOGB0;
}
else
{
// subnormal
uint msb = y.vui[MANTISSA_MSB];
uint lsb = y.vui[MANTISSA_LSB];
if (msb)
return ex - F.EXPBIAS - T.mant_dig + 1 + 32 + core.bitop.bsr(msb);
else
return ex - F.EXPBIAS - T.mant_dig + 1 + core.bitop.bsr(lsb);
}
}
else version (CRuntime_Microsoft)
else static if (F.realFormat == RealFormat.ieeeQuadruple)
{
int res;
asm pure nothrow @nogc
if (ex) // If exponent is non-zero
{
fld real ptr [x] ;
fxam ;
fstsw AX ;
and AH,0x45 ;
cmp AH,0x40 ;
jz Lzeronan ;
cmp AH,5 ;
jz Linfinity ;
cmp AH,1 ;
jz Lzeronan ;
fxtract ;
fstp ST(0) ;
fistp res ;
mov EAX,res ;
jmp Ldone ;
Lzeronan:
mov EAX,0x80000000 ;
fstp ST(0) ;
jmp Ldone ;
Linfinity:
mov EAX,0x7FFFFFFF ;
fstp ST(0) ;
Ldone: ;
if (ex == F.EXPMASK)
{
// infinity or NaN
if (y.vl[MANTISSA_LSB] | ( y.vl[MANTISSA_MSB] & 0x0000_FFFF_FFFF_FFFF)) // NaN
return FP_ILOGBNAN;
else // +- infinity
return int.max;
}
else
{
return ex - F.EXPBIAS - 1;
}
}
else if ((y.vl[MANTISSA_LSB] | (y.vl[MANTISSA_MSB] & 0x0000_FFFF_FFFF_FFFF)) == 0)
{
// vf is +-0.0
return FP_ILOGB0;
}
else
{
// subnormal
ulong msb = y.vl[MANTISSA_MSB] & 0x0000_FFFF_FFFF_FFFF;
ulong lsb = y.vl[MANTISSA_LSB];
if (msb)
return ex - F.EXPBIAS - T.mant_dig + 1 + bsr_ulong(msb) + 64;
else
return ex - F.EXPBIAS - T.mant_dig + 1 + bsr_ulong(lsb);
}
}
else
return core.stdc.math.ilogbl(x);
else static if (F.realFormat == RealFormat.ieeeDouble)
{
if (ex) // If exponent is non-zero
{
if (ex == F.EXPMASK) // infinity or NaN
{
if ((y.vl[0] & 0x7FFF_FFFF_FFFF_FFFF) == 0x7FF0_0000_0000_0000) // +- infinity
return int.max;
else // NaN
return FP_ILOGBNAN;
}
else
{
return ((ex - F.EXPBIAS) >> 4) - 1;
}
}
else if (!(y.vl[0] & 0x7FFF_FFFF_FFFF_FFFF))
{
// vf is +-0.0
return FP_ILOGB0;
}
else
{
// subnormal
uint msb = y.vui[MANTISSA_MSB] & F.MANTISSAMASK_INT;
uint lsb = y.vui[MANTISSA_LSB];
if (msb)
return ((ex - F.EXPBIAS) >> 4) - T.mant_dig + 1 + core.bitop.bsr(msb) + 32;
else
return ((ex - F.EXPBIAS) >> 4) - T.mant_dig + 1 + core.bitop.bsr(lsb);
}
}
else static if (F.realFormat == RealFormat.ieeeSingle)
{
if (ex) // If exponent is non-zero
{
if (ex == F.EXPMASK) // infinity or NaN
{
if ((y.vui[0] & 0x7FFF_FFFF) == 0x7F80_0000) // +- infinity
return int.max;
else // NaN
return FP_ILOGBNAN;
}
else
{
return ((ex - F.EXPBIAS) >> 7) - 1;
}
}
else if (!(y.vui[0] & 0x7FFF_FFFF))
{
// vf is +-0.0
return FP_ILOGB0;
}
else
{
// subnormal
uint mantissa = y.vui[0] & F.MANTISSAMASK_INT;
return ((ex - F.EXPBIAS) >> 7) - T.mant_dig + 1 + core.bitop.bsr(mantissa);
}
}
else // static if (F.realFormat == RealFormat.ibmExtended)
{
core.stdc.math.ilogbl(x);
}
}
alias FP_ILOGB0 = core.stdc.math.FP_ILOGB0;
@ -2683,34 +2779,51 @@ alias FP_ILOGBNAN = core.stdc.math.FP_ILOGBNAN;
@trusted nothrow @nogc unittest
{
assert(ilogb(real.nan) == FP_ILOGBNAN);
assert(ilogb(-real.nan) == FP_ILOGBNAN);
assert(ilogb(0.0) == FP_ILOGB0);
assert(ilogb(-0.0) == FP_ILOGB0);
assert(ilogb(real.infinity) == int.max);
assert(ilogb(-real.infinity) == int.max);
assert(ilogb(2.0) == 1);
assert(ilogb(2.0001) == 1);
assert(ilogb(1.9999) == 0);
assert(ilogb(0.5) == -1);
assert(ilogb(123.123) == 6);
assert(ilogb(-123.123) == 6);
assert(ilogb(0.123) == -4);
assert(ilogb(-double.min_normal) == -1022);
import std.typetuple, std.typecons;
foreach (F; TypeTuple!(float, double, real))
{
alias T = Tuple!(F, int);
T[13] vals = // x, ilogb(x)
[
T( F.nan , FP_ILOGBNAN ),
T( -F.nan , FP_ILOGBNAN ),
T( F.infinity, int.max ),
T( -F.infinity, int.max ),
T( 0.0 , FP_ILOGB0 ),
T( -0.0 , FP_ILOGB0 ),
T( 2.0 , 1 ),
T( 2.0001 , 1 ),
T( 1.9999 , 0 ),
T( 0.5 , -1 ),
T( 123.123 , 6 ),
T( -123.123 , 6 ),
T( 0.123 , -4 ),
];
foreach(elem; vals)
{
assert(ilogb(elem[0]) == elem[1]);
}
}
// min_normal and subnormals
assert(ilogb(-float.min_normal) == -126);
// subnormals
assert(ilogb(nextUp(-double.min_normal)) == -1023);
assert(ilogb(nextUp(-0.0)) == -1074);
assert(ilogb(nextUp(-float.min_normal)) == -127);
assert(ilogb(nextUp(-0.0F)) == -149);
static if (floatTraits!(real).realFormat == RealFormat.ieeeExtended) {
assert(ilogb(nextUp(-float(0.0))) == -149);
assert(ilogb(-double.min_normal) == -1022);
assert(ilogb(nextUp(-double.min_normal)) == -1023);
assert(ilogb(nextUp(-double(0.0))) == -1074);
static if (floatTraits!(real).realFormat == RealFormat.ieeeExtended)
{
assert(ilogb(-real.min_normal) == -16382);
assert(ilogb(nextUp(-real.min_normal)) == -16383);
assert(ilogb(nextUp(-0.0L)) == -16445);
} else static if (floatTraits!(real).realFormat == RealFormat.ieeeDouble) {
assert(ilogb(nextUp(-real(0.0))) == -16445);
}
else static if (floatTraits!(real).realFormat == RealFormat.ieeeDouble)
{
assert(ilogb(-real.min_normal) == -1022);
assert(ilogb(nextUp(-real.min_normal)) == -1023);
assert(ilogb(nextUp(-0.0L)) == -1074);
assert(ilogb(nextUp(-real(0.0))) == -1074);
}
}